Let $M = (M, g)$ be a Riemannian manifold, and let $p \in M$.

Writing $d$ for the geodesic distance in $M$, there is a function
$$
d(-, p)^2 : M \to \mathbb{R}.
$$
This function is smooth near $p$. Hence for each point $x \in M$ sufficiently close to $p$, we have the Hessian
$$
\text{Hess}_x(d(-, p)^2)
$$
(defined using the Levi-Civita connection), which is a bilinear form on $T_x M$. In particular, we can take $x$ to be equal to $p$ itself, giving a bilinear form
$$
\text{Hess}_p(d(-, p)^2)
$$
on $T_p M$. But of course, we already have another bilinear form on $T_p M$, namely, the Riemannian metric $g_p$ itself. And the fact is that up to a constant factor, these two forms are equal:
$$
g_p = \frac{1}{2} \text{Hess}_p(d(-, p)^2).
$$
**I'm looking for a reference for this fact.** For the purposes of what I'm writing, it would ideally be a reference that states this fact in the same simple direct terms as above, without involving any other differential-geometric concepts (e.g. normal coordinates).

I understand that this is a basic fact of Riemannian geometry, so I've already looked for it in various introductions to the subject, including those by do Carmo, Jost, Lee, and Petersen. But I haven't found it stated in any of those sources (which isn't to say it's not there). I have found more sophisticated stuff about $\text{Hess}_x(d(-, p)^2)$ for points $x$ *different* from $p$, but not the simple fact I'm looking for.

Requests for references often result in people giving their favourite proofs rather than a reference. While that doesn't do any harm (and can be quite interesting), I emphasize that it's a **reference** I'm looking for, not a proof.