I would love to understand the famous formula $g_{ij}(x) = \delta_{ij} + \frac{1}{3}R_{kijl}x^kx^l +O(\|x\|^3)$, which is valid in Riemannian normal coordinates and possibly more general situations.
I'm aware of 2 proofs: One using Jacobi fields [cf. e.g. S.Sternberg's "Curvature in Mathematics and Physics" from which the question title and formula is stolen 🙂 or cf. S.Lang's "Differential and Riemannian Manifolds"]. The other proof involves computing that $\partial_k\partial_lg_{ij}(x)$ shares some symmetries of curvature [cf. M.Spivak's "A Comprehensive Introduction to Differential Geometry, Vol. 2" where it is a several page "hairy computation" or cf. H.Weyl's 1923 edition of Riemann's Habilitationsvortrag (reprinted in a recent German book by Jürgen Jost) which I find uncomprehensible.]
Are you aware of any other proof? Are normal coordinates necessary?
While the Jacobi fields proof is short and elegant enough, it irks me that it requires "higher technology" not involved in the endproduct. Somehow the formula should be provable by pure calculus. Indeed, it is stated as an exercise in P.Petersen's "Riemannian Geometry": From the context I guess he thinks it should follow from the expression of $\partial_lg_{ij}$ as a sum of 2 Christoffel symbols and the simplified expression for curvature at $x=0$ where the Christoffel symbols vanish. Alas my attempts at this go in circles…
I find the situation quite amazing: Not many textbooks treat this fundamental and historic formula. (Estimating from the sample on my shelf it is $3/17.$ E.g. it seems it's not even in Levi-Civita's classic.)
Update/Scholium:
In classical language: The knackpoint seems to be a "differential Bianchi formula" for the Christoffel symbols at $0$. This follows from the geodesic equation. I see no other way yet.
A more modern approach minimizing (but not eliminating) the role of geodesics is in A.Gray's Tubes book. (Noted in comments. I'm waiting for www.amazon.de to deliver this treasure.)
$\bullet$ While geodesics are very geometric and normal coordinates are very practical, methinks the formula is a tad ungeometric. What I'm hoping/asking for is a coordinate-independent formula for the second derivative of $g$ in terms of a suitable "reference connection".
Best Answer
Perhaps the simplest way to understand this formula is to think about how you would go about deriving it: Try to find the 'best' coordinates you can centered on a given point and see what doesn't change in such coordinates.
Suppose $g$ is a Riemannian metric on $M$ and $p\in M$ is fixed. Start by choosing a $p$-centered local coordinate system $x = (x^i)$ on $U\subset M$ and write $$ g = g_{ij}(x)\,\mathrm{d}x^i\mathrm{d}x^j. $$ Since $\bigl(g_{ij}(0)\bigr)$ is a positive definite matrix, you can make a linear change of coordinates in $x$ so that $g_{ij}(0) = \delta_{ij}$. Call such a $p$-centered coordinate system $0$-adapted to $g$ at $p$.
Now, ask what would be the effect of expressing $g$ in the coordinates $y=(y^i)$ that are related to the coordinate $x$ by $x^i = y^i + \tfrac12a^i_{jk} y^jy^k$ for some $a^i_{jk} = a^i_{kj}$. It is easy to see by Taylor series expansion that you can uniquely choose the $a^i_{jk}$ so that, when we write $$ g = \bar g_{ij}(y)\,\mathrm{d}y^i\mathrm{d}y^j, $$ we have, for all $i$, $j$, and $k$, $$ \frac{\partial\bar g_{ij}}{\partial y^k}(0) = 0. $$ (It's clear that this is the same number of equations as unknowns for the $a^i_{jk}$, one just has to check that the inhomogeneous system of equations has only the zero solution when the inhomogeneous part is set to zero.) Call such a system of $p$-centered coordinates $1$-adapted to $g$ at $p$. Thus, for a system of coordinates $y$ that is $1$-adapted to $g$ at $p$, one has $$ g = \left(\delta_{ij} + \tfrac12 \frac{\partial^2g_{ij}}{\partial y^k\partial y^l}(0)\, y^ky^l + R^3_{ij}(y)\right) \ \mathrm{d}y^i\mathrm{d}y^j, $$ where $R^3_{ij}(y)$ vanishes to order $3$ at $y=0$.
Finally, consider what such a metric would look like in the coordinates $z = (z^i)$ that are defined by $y^i = z^i + \tfrac16 b^i_{jkl} z^jz^kz^l$ for some constants $b^i_{jkl} = b^i_{kjl} = b^i_{jlk}$. Now, there are $n^2(n{+}1)(n{+}2)/6$ unknowns $b^i_{jkl}$, but there are $n^2(n{+}1)^2/4$ quantities in the second-order Taylor expansion of $g = {\bar g}_{ij}(z)\mathrm{d}z^i\mathrm{d}z^j$, i.e., $$ g = \left(\delta_{ij} + \tfrac12 \frac{\partial^2{\bar g}_{ij}}{\partial z^k\partial z^l}(0)\, z^kz^l + {\bar R}^3_{ij}(z)\right) \ \mathrm{d}z^i\mathrm{d}z^j. $$ Thus, the equations $\frac{\partial^2{\bar g}_{ij}}{\partial z^k\partial z^l}(0)=0$, as linear equations for the $b^i_{jkl}$, are overdetermined by $$ n^2(n{+}1)^2/4 - n^2(n{+}1)(n{+}2)/6 = n^2(n^2{-}1)/12 $$ equations.
It is not hard to see that the corresponding homogeneous equations in the $b^i_{jkl}$ have only the solution $b^i_{jkl}=0$. In fact, the $b^i_{jkl}$ are uniquely determined by requiring that, when we compute the Taylor expansion about $z=0$ we get $$ g = \left(\delta_{ij} + \tfrac12 h_{ij,kl}\, z^kz^l + R^3_{ij}(z)\right) \ \mathrm{d}z^i\mathrm{d}z^j, $$ with $h_{ij,kl}+h_{ik,lj}+h_{il,jk}=0$ (which is $n^2(n{+}1)(n{+}2)/6$ independent equations on the $b^i_{jkl}$). Say that a system of coordinates $z = (z^i)$ centered at $p$ for which $g$ has its Taylor expansion at $p$ of the above form is $2$-adapted to $g$ at $p$. Two such coordinate systems at $p$ are related in the form $z^i = a^i_j\,\bar z^j + O(|{\bar z}|^4)$, where $a = (a^i_j)$ is an orthogonal matrix.
Thus, the $2$-adapted condition forces the $h_{ij,kl}$ to lie in a vector space of dimension $n^2(n^2{-}1)/12$, as explained above.
It's now a matter of linear algebra to show, as Riemann did, that these conditions imply that the $h_{ij,kl}$ can be written uniquely in the form $$ h_{ij,kl} = \tfrac13(R_{kijl}+R_{lijk}) $$ where $R_{ijkl}=-R_{jikl}=-R_{ijlk}$ and $R_{ijkl}+R_{iklj}+R_{iljk}=0$.