How to take the "gradient" of a matrix? I want to find $(D^2 F)$ where $\vec{F}(\vec{x}) = \frac{\vec{x}}{\|\vec{x}\|}$ ($F$ is a row vector and $Du$ is a "column vector", where $u\in \mathbb{R}$ ).
I know that 
$$(DF)_{ij} = \frac{\delta_{ij}\|x\|^2 - 2 x_ix_j}{\|x\|^4}\qquad\rightarrow\qquad (DF) = \frac{I}{\|x\|^2} - \frac{2xx^T}{\|x\|^4}.$$
How exactly do I find $(D^2F)$? I tried reading the section on wikipedia and it hasn't been productive so far.
 A: Maybe it helps when you consider derivatives as linear operators. This means if you have
$$
F\colon \mathbb R^n \to \mathbb R^n
$$
you consider
$$
 DF\colon \mathbb R^n \to L(\mathbb R^n,\mathbb R^n),
$$
where $L(A,B)$ is the set of all linear maps from $A$ to $B$. Usually, $L(\mathbb R^n,\mathbb R^n)$ is identified with the set of matrices $\mathbb R^{n\times n}$.
Now consider $D^2F = D(DF)$ as
$$
D^2F = D(DF)\colon \mathbb R^n \to L\bigl(\mathbb R^n,L(\mathbb R^n,\mathbb R^n)\bigr)
$$
where $L\bigl(\mathbb R^n,L(\mathbb R^n,\mathbb R^n)\bigr)$ is the set of linear maps from $\mathbb R^n$ into the set of linear mappings from $\mathbb R^n$ into $\mathbb R^n$. You could identify this as $\mathbb R^{n\times n\times n}$.
This means $D^2F$ can be considered some kind of tensor.
The easiest way to tackle your question would be to calculate the components $(D^2F)_{ijk}$ of $D^2F$. Since you already have $DF$ in its matrix form you can differentiate it wrt. to $x_k$. You have $\|x\|^2 = x^Tx$
\begin{align*}
\frac{\partial}{\partial x_k} x^Tx &= 2x_k \\
\frac{\partial}{\partial x_k} \frac{I}{x^Tx} &= \frac{-2Ix_k}{(x^Tx)^2}
= \frac{-2Ix_k}{\|x\|^4}\\
\frac{\partial}{\partial x_k} (x^Tx)^2 &= 4x_k (x^Tx) = 4x_k \|x\|^2 \\
\frac{\partial}{\partial x_k} xx^T &= e_kx^T + x e_k^T \\
\frac{\partial}{\partial x_k} \frac{2xx^T}{(x^Tx)^2}
&= 2\frac{(e_kx^T + x e_k^T)(x^Tx)^2 - xx^T 4x_k\|x\|^2}{(x^Tx)^4} \\
\end{align*}
If you put it together you get
$$
\frac{\partial}{\partial x_k} DF(x)
= \frac{-2Ix_k}{\|x\|^4} - 2\frac{(e_kx^T + x e_k^T)\|x\|^2 - 4x_kxx^T }{\|x\|^6}
$$
You now have the matrix $(D^2F)_{\bullet,\bullet,k} = \frac{\partial}{\partial x_k} DF(x)$. You could write $D^2F = \left(\frac{\partial}{\partial x_k} DF(x)\right)_{k=1}^n$.
Hope this gives some insight.
(Please point out any calculation errors, if present)
A: (I assume that $F$ should be $\frac{x}{\|x\|^2}$ and not $\frac{x}{\|x\|}$.)
When you take the second vector derivative of a vector, you will get a tensor of rank 3. We use simplified Einstein notation (where we don't care about distinguishing upper and lower indices):
$$ F_k = \frac{x_k}{\|x\|^2}. $$
Applying $\partial_j$ gives
$$
\begin{align}
\partial_j F_k
&= \frac{\|x\|^2 \partial_j x_k - x_k \partial_j \|x\|^2}{\|x\|^4} \\
&= \frac{\|x\|^2 \delta_{jk} - x_k \partial_j (x_m x_m)}{\|x\|^4} \\
&= \frac{\delta_{jk}}{\|x\|^2} - \frac{2x_j x_k}{\|x\|^4}. \\
\end{align}
$$
where we have used $\partial_j(x_m x_m) = 2 x_m \delta_{jm} = 2x_j$.
Next we apply $\partial_i$, doing derivatives in terms of $\|x\|^2 = x_m x_m$ where necessary:
$$
\begin{align}
\partial_i \partial_j F_k
&= \partial_i \left( \frac{\delta_{jk}}{\|x\|^2} - \frac{2x_j x_k}{\|x\|^4} \right) \\
&= -\frac{\delta_{jk}}{\|x\|^4} \partial_i \|x\|^2
   - \frac{2}{\|x\|^8}
     \left(
       \|x\|^4 \partial_i (x_j x_k) - x_j x_k \partial_i (\|x\|^2)^2
     \right) \\
&= -\frac{\delta_{jk}}{\|x\|^4} \partial_i (x_m x_m)
   - \frac{2}{\|x\|^8}
     \left(
       \|x\|^4 (x_j \delta_{ik} + x_k \delta_{ij})
      - x_j x_k \cdot 2 \|x\|^2 \partial_i \|x\|^2
     \right) \\
&= -\frac{\delta_{jk}}{\|x\|^4} \cdot 2x_i
   - \frac{2}{\|x\|^8}
     \left(
       \|x\|^4 (x_j \delta_{ik} + x_k \delta_{ij})
      - x_j x_k \cdot 2 \|x\|^2 \cdot 2 x_i
     \right) \\
&= -2 \cdot \frac{x_i \delta_{jk} + x_j \delta_{ki} + x_k \delta_{ij}}{\|x\|^4}
+ 8 \cdot \frac{x_i x_j x_k}{\|x\|^6}
\end{align}
$$
Updated. Now with a visual representation of the result:

A: It is a good question. You would need some way to way to express it, so in some sense it would become more of a question of how to express it. We will maybe need some way to express "transpose to third dimension". This would be able to express $(D^2F)_{..k}$ in Wauzl's answer as a kind of generalized "outer product" between the nabla operator and a matrix. I think tensor algebra provides such notation, but I am not very confident explaining it.
