[proofplan] We compare the two top-degree forms $(\operatorname{div} X)\,\omega_g$ and $d(i(X)\omega_g)$ by evaluating both on a frame at an arbitrary point $p \in M$. The Leibniz rule for $\nabla$ applied to the contraction $i(X)$, combined with metric compatibility $\nabla\omega_g = 0$, reduces $\nabla_Y(i(X)\omega_g)$ to $i(\nabla_Y X)\omega_g$. Choosing a normal orthonormal frame at $p$ converts the exterior derivative $d$ into an antisymmetrisation of $\nabla$ at $p$, and the trace structure of the resulting expression collapses it to $\operatorname{tr}(\nabla X)\,\omega_g = (\operatorname{div} X)\,\omega_g$. Both sides are tensorial in $X$, so the pointwise identity at every $p$ gives the global identity. [/proofplan]

[step:Reduce to a pointwise computation in a normal orthonormal frame] Let $(M, g)$ be an oriented Riemannian manifold, $\nabla$ the Levi-Civita connection, $\omega_g$ the Riemannian volume form, $X \in \mathfrak{X}(M)$, and $i(X) : \Omega^k(M) \to \Omega^{k-1}(M)$ the interior product (contraction with $X$ in the first slot). Both sides of the identity \begin{align*} (\operatorname{div} X)\,\omega_g = d(i(X)\,\omega_g) \tag{$\star$} \end{align*} are top-degree forms on $M$, and both depend pointwise on $X$ and its first covariant derivative. Hence ($\star$) is an equality of tensorial expressions, and it suffices to verify it at an arbitrary point $p \in M$. Fix $p \in M$ and choose a **normal orthonormal frame** $\{e_1, \ldots, e_n\}$ at $p$: a local orthonormal frame on a neighbourhood of $p$ with $\nabla_Y e_i|_p = 0$ for all $Y \in T_p M$ and all $i = 1, \ldots, n$. Such a frame exists; one can construct it by parallel-transporting an orthonormal basis of $T_p M$ along radial geodesics in a normal coordinate chart centred at $p$. Let $\{\omega^1, \ldots, \omega^n\}$ be the dual coframe, so $\omega^i(e_j) = \delta_{ij}$. In this frame at $p$, \begin{align*} \omega_g|_p = \omega^1 \wedge \cdots \wedge \omega^n|_p. \end{align*} [/step]

[step:Compute $\nabla_Y(i(X)\omega_g)$ via the Leibniz rule and metric compatibility] The covariant derivative $\nabla$ extends to all tensor bundles as a derivation that commutes with contractions. Applied to the contraction $i(X)\psi$ of a vector field $X$ with a $k$-form $\psi$, this gives the Leibniz identity \begin{align*} \nabla_Y(i(X)\psi) = i(\nabla_Y X)\psi + i(X)\nabla_Y\psi \qquad \text{for all } Y \in \mathfrak{X}(M). \tag{L} \end{align*} Specialise (L) to $\psi := \omega_g$. The Levi-Civita connection is metric-compatible, and the Riemannian volume form is determined by $g$ together with the orientation; consequently $\nabla\omega_g = 0$, i.e., $\nabla_Y \omega_g = 0$ for all $Y \in \mathfrak{X}(M)$. The second term in (L) drops, leaving \begin{align*} \nabla_Y(i(X)\omega_g) = i(\nabla_Y X)\,\omega_g. \tag{L$'$} \end{align*} [/step]

[step:Express $d$ as an antisymmetrisation of $\nabla$ at $p$] For any $(n-1)$-form $\eta \in \Omega^{n-1}(M)$, the exterior derivative on a torsion-free connection satisfies the global identity \begin{align*} d\eta(Y_0, Y_1, \ldots, Y_{n-1}) = \sum_{j=0}^{n-1} (-1)^j (\nabla_{Y_j}\eta)(Y_0, \ldots, \widehat{Y_j}, \ldots, Y_{n-1}), \end{align*} where $\widehat{\cdot}$ denotes omission. This is the standard formula expressing $d$ as the antisymmetrisation of $\nabla$ (which holds because $\nabla$ is torsion-free). Apply this to $\eta := i(X)\omega_g$ at $p$ with $(Y_0, \ldots, Y_{n-1}) := (e_1, \ldots, e_n)$: \begin{align*} d(i(X)\omega_g)(e_1, \ldots, e_n)|_p = \sum_{j=1}^{n} (-1)^{j-1} (\nabla_{e_j}(i(X)\omega_g))(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p. \end{align*} By (L$'$) at $p$, $\nabla_{e_j}(i(X)\omega_g)|_p = i(\nabla_{e_j}X)\omega_g|_p$, so \begin{align*} d(i(X)\omega_g)(e_1, \ldots, e_n)|_p = \sum_{j=1}^{n} (-1)^{j-1} \big(i(\nabla_{e_j}X)\omega_g\big)(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p. \tag{D} \end{align*} [/step]

[step:Evaluate the contraction sum and identify the trace]Write $\nabla_{e_j} X|_p = \sum_{k=1}^n A_{kj}\, e_k|_p$ in the frame at $p$, so that the matrix of $\nabla X|_p : T_p M \to T_p M$ in the basis $\{e_k\}$ is $(A_{kj})$. By definition of $\operatorname{div}$, \begin{align*} (\operatorname{div} X)(p) = \operatorname{tr}(\nabla X|_p) = \sum_{k=1}^n A_{kk}. \end{align*} Compute one term of the sum (D). For fixed $j$, \begin{align*} \big(i(\nabla_{e_j}X)\omega_g\big)(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p &= \omega_g(\nabla_{e_j}X, e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p \\ &= \sum_{k=1}^n A_{kj}\, \omega_g(e_k, e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p. \end{align*} The form $\omega_g|_p = \omega^1 \wedge \cdots \wedge \omega^n|_p$ evaluated on the tuple $(e_k, e_1, \ldots, \widehat{e_j}, \ldots, e_n)$ vanishes unless $\{k, 1, \ldots, \widehat{j}, \ldots, n\} = \{1, \ldots, n\}$, i.e., unless $k = j$. When $k = j$, the tuple $(e_j, e_1, \ldots, \widehat{e_j}, \ldots, e_n)$ is the cyclic shift bringing $e_j$ to the front of $(e_1, \ldots, e_n)$, achieved by $j-1$ transpositions, so \begin{align*} \omega_g(e_j, e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p = (-1)^{j-1} \omega_g(e_1, \ldots, e_n)|_p = (-1)^{j-1}. \end{align*} Therefore \begin{align*} \big(i(\nabla_{e_j}X)\omega_g\big)(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p = (-1)^{j-1} A_{jj}. \end{align*} Substituting into (D), \begin{align*} d(i(X)\omega_g)(e_1, \ldots, e_n)|_p = \sum_{j=1}^n (-1)^{j-1} \cdot (-1)^{j-1} A_{jj} = \sum_{j=1}^n A_{jj} = (\operatorname{div} X)(p). \end{align*} On the other hand, \begin{align*} \big((\operatorname{div} X)\,\omega_g\big)(e_1, \ldots, e_n)|_p = (\operatorname{div} X)(p) \cdot \omega_g(e_1, \ldots, e_n)|_p = (\operatorname{div} X)(p). \end{align*} The two top-degree forms agree on the frame $(e_1, \ldots, e_n)$ at $p$. Since both are alternating $n$-forms on the $n$-dimensional space $T_p M$, agreement on one ordered basis forces agreement as forms at $p$: \begin{align*} \big(d(i(X)\omega_g)\big)|_p = \big((\operatorname{div} X)\,\omega_g\big)|_p. \end{align*} The point $p \in M$ was arbitrary, so this is the asserted global identity ($\star$).[/step]

[guided]We have arrived at formula (D), \begin{align*} d(i(X)\omega_g)(e_1, \ldots, e_n)|_p = \sum_{j=1}^{n} (-1)^{j-1} \big(i(\nabla_{e_j}X)\omega_g\big)(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p, \end{align*} and we want to evaluate this and recognise the right-hand side as $(\operatorname{div} X)(p)$. To do so we must extract numerical content from each $j$-th summand. The natural strategy is to expand $\nabla_{e_j} X$ in the orthonormal frame, plug into the volume form, and watch the alternating signs collide with the signs coming from the cyclic structure of $\omega_g$. **Expand $\nabla_{e_j} X$ in the frame.** The vector $\nabla_{e_j} X|_p \in T_p M$ has coordinates in the basis $\{e_k|_p\}$: \begin{align*} \nabla_{e_j} X|_p = \sum_{k=1}^n A_{kj}\, e_k|_p, \end{align*} where the matrix $(A_{kj})$ is precisely the matrix of the endomorphism $\nabla X|_p : T_p M \to T_p M$ in the basis $\{e_k\}$. Why does this help? Because the divergence is by definition the trace of this endomorphism: \begin{align*} (\operatorname{div} X)(p) = \operatorname{tr}(\nabla X|_p) = \sum_{k=1}^n A_{kk}. \end{align*} So our goal is to show that the right-hand side of (D) collapses to $\sum_k A_{kk}$ — only diagonal entries should survive. **Compute one summand of (D).** Fix $j$. Using the definition $i(V)\omega_g(\cdot) = \omega_g(V, \cdot)$ and the linearity of $\omega_g$ in its first argument, \begin{align*} \big(i(\nabla_{e_j}X)\omega_g\big)(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p &= \omega_g(\nabla_{e_j}X, e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p \\ &= \sum_{k=1}^n A_{kj}\, \omega_g(e_k, e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p. \end{align*} Now we ask: which $k$ contribute? At $p$ the form $\omega_g|_p = \omega^1 \wedge \cdots \wedge \omega^n|_p$ is alternating, so it vanishes on any tuple with a repeated entry. The tuple $(e_k, e_1, \ldots, \widehat{e_j}, \ldots, e_n)$ already contains every $e_i$ with $i \neq j$. If $k \neq j$, then $e_k$ duplicates one of these, and $\omega_g$ vanishes. The only surviving term is $k = j$ — which is exactly the diagonal entry $A_{jj}$. This is the mechanism by which off-diagonal entries are killed. **Compute the sign for $k = j$.** When $k = j$, we evaluate $\omega_g$ on $(e_j, e_1, \ldots, \widehat{e_j}, \ldots, e_n)$. To compare this with the standard ordering $(e_1, \ldots, e_n)$ on which $\omega_g|_p$ takes the value $+1$, count the transpositions needed. The tuple is obtained from $(e_1, \ldots, e_n)$ by extracting $e_j$ from its position and pulling it to the front — that is, by $j-1$ adjacent transpositions. Since $\omega_g$ is alternating, each transposition changes the sign: \begin{align*} \omega_g(e_j, e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p = (-1)^{j-1} \omega_g(e_1, \ldots, e_n)|_p = (-1)^{j-1}. \end{align*} Combining, \begin{align*} \big(i(\nabla_{e_j}X)\omega_g\big)(e_1, \ldots, \widehat{e_j}, \ldots, e_n)|_p = (-1)^{j-1} A_{jj}. \end{align*} **Watch the signs cancel.** Substitute back into (D). The summand carries two sign factors of $(-1)^{j-1}$: one from the antisymmetrisation in the $d$-formula, and one from the cyclic permutation we just computed. They square to $+1$: \begin{align*} d(i(X)\omega_g)(e_1, \ldots, e_n)|_p = \sum_{j=1}^n (-1)^{j-1} \cdot (-1)^{j-1} A_{jj} = \sum_{j=1}^n A_{jj} = (\operatorname{div} X)(p). \end{align*} This is the heart of the proof: the alternating signs from $d$ and the alternating signs from $\omega_g$ exactly compensate, leaving the bare trace. **Compare with the left-hand side.** The form $(\operatorname{div} X)\,\omega_g$ is, by definition, the function $\operatorname{div} X$ multiplied by the volume form. Evaluating on the frame at $p$, \begin{align*} \big((\operatorname{div} X)\,\omega_g\big)(e_1, \ldots, e_n)|_p = (\operatorname{div} X)(p) \cdot \omega_g(e_1, \ldots, e_n)|_p = (\operatorname{div} X)(p), \end{align*} using $\omega_g(e_1, \ldots, e_n)|_p = 1$ since $\{e_i\}$ is a positively oriented orthonormal basis at $p$. **Conclude pointwise then globally.** We have shown \begin{align*} d(i(X)\omega_g)(e_1, \ldots, e_n)|_p = (\operatorname{div} X)(p) = \big((\operatorname{div} X)\,\omega_g\big)(e_1, \ldots, e_n)|_p. \end{align*} Why does agreement on a single ordered basis suffice? Both sides are alternating $n$-forms on the $n$-dimensional vector space $T_p M$, and the space of alternating $n$-forms on an $n$-dimensional space is one-dimensional. An alternating $n$-form is determined by its value on any single ordered basis, so equality on $(e_1, \ldots, e_n)$ forces \begin{align*} \big(d(i(X)\omega_g)\big)|_p = \big((\operatorname{div} X)\,\omega_g\big)|_p. \end{align*} Finally, the point $p \in M$ was arbitrary. As discussed in the first step, both sides of ($\star$) are tensorial in $X$ and depend only on $X$ and $\nabla X$ pointwise, so a pointwise identity at every $p$ is exactly the global identity ($\star$). The proof is complete. Two retrospective remarks on the strategy. First, the choice of a normal orthonormal frame served a single purpose: it made $\nabla_{e_j} e_k|_p = 0$ so that the global Koszul-type formula for $d$ reduced at $p$ to the naive antisymmetrisation of $\nabla$. Without the normal-frame condition, (D) would have carried extra terms of the form $\eta(\ldots, \nabla_{e_j} e_k, \ldots)$ that complicate the cyclic-shift bookkeeping. Second, the entire computation is a "Cramer-style" extraction of a trace: the contraction $i(X)$ followed by $d$ produces, via the alternating structure of $\omega_g$, exactly the trace of $\nabla X$ — the same trace that defines $\operatorname{div} X$. This is why the divergence and $d \circ i(X)$ on the volume form are the same object.[/guided]