Second order derivatives of a function. We have to use the Hessian to calculate the second order derivatives of a function. While that is okay if the function is mapped from $\mathbb{R}^n$ to $\mathbb{R}$, how does one proceed if it is mapped from $\mathbb{R}^n \longrightarrow \mathbb{R}^m$, where $m > 1$? How does one take the derivative of f with respect to each variable when f is itself in more than one dimension?
For example, if there is a function f : $\mathbb{R}^2 \longrightarrow \mathbb{R}^3$ with $(x,y)$ mapped to $(x^2 + y, y^3, \cos(y))$, how does one calculate the Hessian?
Thank you
 A: I guess what you are looking for is the following:
Let $\vec{f}(x,y) = \begin{pmatrix} u(x,y) \\ v(x,y) \\ w(x,y)\end{pmatrix}$. Its "Hessian" is simply the "vector-valued matrix"
$$ H_{\vec{f}} = \begin{pmatrix}
\begin{pmatrix} u_{xx} \\ v_{xx} \\ w_{xx}\end{pmatrix} & 
\begin{pmatrix} u_{xy} \\ v_{xy} \\ w_{xy}\end{pmatrix} \\
\begin{pmatrix} u_{yx} \\ v_{yx} \\ w_{yx}\end{pmatrix} &
\begin{pmatrix} u_{yy} \\ v_{yy} \\ w_{yy}\end{pmatrix}
\end{pmatrix} $$
where $ u_{xy} = \frac{\partial^2}{\partial x\partial y} u$ etc. Note that each entry of the "matrix" is now a vector in its own right. 
A: The importance of the Hessian comes from the following fact:

Let $X$ and $Y$ be normed spaces, and let $L(X,Y)$ denote the space of continuous linear maps from $X$ to $Y$. There is a canonical isomorphism (and actually an isometry) $$L(X,L(X,Y)) \simeq \operatorname{Bil}(X\times X\to Y),$$ where Bil is the set of bilinear continuous maps from $X \times X$ to $Y$.

Since the second derivative belongs to $L(X,L(X,Y))$, we can identify it with a bilinear map. If $f \colon \mathbb{R}^n \to \mathbb{R}$, then the bilinear maps from $\mathbb{R}^n \times \mathbb{R}^n$ to $\mathbb{R}$ can be represented by $(n\times n)$ matrices, and you get the Hessian matrix. If the codomain of $f$ has higher dimension, matrices no longer suffice. I guess you'd need tensors. See also here.
Additional material: assume that $\mathbf{f} \colon \mathbb{R}^n \to \mathbb{R}^m$ is of class $C^2$. Then the second differential at $\mathbf{x}_0 \in \mathbb{R}^n$ acts like $$\mathrm{d}^2 \mathbf{f}(\mathbf{x}_0)\colon (\mathbf{v}_1,\mathbf{v}_2) \mapsto \sum_{i,j=1}^n \frac{\partial}{\partial x_{i}} \frac{\partial \mathbf{f}}{\partial x_{j}}(\mathbf{x}_0) v_{1,i}v_{2,j}$$ on the generic vectors $\mathbf{v}_1$, $\mathbf{v}_2 \in \mathbb{R}^n$, and $v_{1,i}$, $v_{2,j}$ are the components of $\mathbf{v}_1$ and $\mathbf{v}_2$, respectively.
