When is the Hessian contracted with a vector field a closed form? Suppose M is a Riemannian manifold and $g, f : M \rightarrow \mathbb{R}$ are smooth functions.
When is the $1$-form Hess$^f(\nabla g, -)$ closed?
I'm looking for simple conditions involving f,g and possibly the geometry of M.
Thank you.
 A: Let us change the notation. Using abstract indices is just a taste preference: you can convert into your choice afterwards.
Actually, I would also change the acting characters, as the letter $g$ in Riemannian geometry (and General Relativity!) usually means the Riemannian metric on the manifold $M$ in the consideration. So, let us talk about functions $\phi, \psi \colon M \to \mathbb{R}$.
The Hessian $\mathrm{Hess}^{\phi}$ of a function $f$ will be written simply as $(\nabla \nabla \phi)_{ab}$ or, sligthly abusively, as $\nabla_a \nabla_b \phi$. The condition $\nabla_a \nabla_b \phi = \nabla_b \nabla_b \phi$ (for any function $\phi$) is precisely the torsion freeness of the Levi-Civita connection $\nabla$ that is tacitly used by default in Riemannian geometry. This way we see that the Hessian is a symmetric $2$-tensor.
The gradient $\nabla \psi$ of a function $\psi$ in the abstract index notation is written as $\nabla^a \psi$ meaning $\nabla^a \psi := g^{a b} \nabla_b \psi$ because $\nabla_a \psi$ is just the differential of $\psi$.
Now your quantity $\mathrm{Hess}^{\phi}(\nabla \psi,-)$ is represented by a $1$-form $\omega_b := (\nabla^a \psi) \nabla_a \nabla_b \phi$, and you are asking for the condition $\mathrm{d} \omega = 0$ in terms of the functions $\phi,\psi$ and the geometry of the manifold, that is the metric $g$.
The exterior covariant derivative $\mathrm{d} \omega$ is a $2$-form that can be expressed in terms of any torsion free connection $\nabla$ as 
$$
(\mathrm{d} \omega)_{a b} = \nabla_a \omega_b - \nabla_b \omega_a
$$ (the Levi-Civita connection $\nabla$ of $g$ will do).
So, we can calculate
$$
\begin{align}
(\mathrm{d} \omega)_{a b} & = \nabla_a \omega_b - \nabla_b \omega_a = \nabla_a \Big((\nabla^c \psi) \nabla_c \nabla_b \phi\Big) - \nabla_b \Big((\nabla^c \psi) \nabla_c \nabla_a \phi\Big) \\
& =  \Big((\nabla_a\nabla^c \psi) \nabla_c \nabla_b \phi + (\nabla^c \psi) \nabla_a \nabla_c \nabla_b \phi\Big) - (a \leftrightarrow b)
\end{align}
$$
where $(a \leftrightarrow b)$ means the first half of the expression with the indices $a$ and $b$ interchanged.
The term $(\nabla_a\nabla^c \psi) \nabla_c \nabla_b \phi$ is a symmetric $2$-tensor and will vanish after the antisymmetrization "$-(a \leftrightarrow b)$", so we only need to look at the term $(\nabla^c \psi) \nabla_a \nabla_c \nabla_b \phi$ in which we would like to commute the covariant derivatives to obtain $(\nabla^c \psi) \nabla_c \nabla_a \nabla_b \phi$ but this will  produce a curvature term:
$$
\nabla_{a} \nabla_{c} \nabla_{b} \phi = \nabla_{c} \nabla_{a} \nabla_{b} \phi - \mathrm{R}_{a c}{}^d{}_b \nabla_{d} \phi
$$
At this point, using the torsion-freeness of the Levi-Civita connection ($\nabla_{[a} \nabla_{b]} \phi = 0$) we see that
$$
(\mathrm{d} \omega)_{a b} = - 2\, (\nabla^c \psi) \mathrm{R}_{[a| c}{}^d{}_{|b]} \nabla_{d} \phi = 2 \, \mathrm{R}_{c[ab]d}(\nabla^c \phi)(\nabla^d \psi)
$$
For instance, the condition $(\mathrm{d} \omega)_{a b} = 0$ holds when $\mathrm{R}=0$ (the Riemannian structure is flat). 
Now I change the notation back to present the answer in a consumable form.

Proposition. If $\phi$ and $\psi$ are smooth functions on a Riemannian manifold $(M,g)$, the $1$-form $\mathrm{Hess}^{\phi}(\nabla \psi,-)$ is closed if and only if 
  $$
\mathrm{R}(\nabla \phi, X, Y, \nabla \psi) - \mathrm{R}(\nabla \phi, Y, X, \nabla \psi) = 0 \tag{*}
$$
  for all tangent vectors $X,Y$ at each point of the manifold $M$. 
  Here $\mathrm{R}$ is the Riemannian curvature of metric $g$ and $\nabla \phi$ is the gradient of function $\phi$.

This analysis can be continued, of course. One may notice that the condition $(*)$ is purely local, so the gradients can be replaced by arbitrary vectors $Z,W$. I leave this to the interested reader :-) 
