[proofplan] The strategy is to show that locally, on each small subinterval, $\gamma$ coincides with a smooth radial geodesic from a normal-coordinate chart centred at a point on $\gamma$. Without loss of generality, assume $\gamma$ has unit speed (rescale time). Fix any $t_0 \in [0, 1]$ and let $p_0 := \gamma(t_0)$. Choose $\varepsilon_0 > 0$ such that $\exp_{p_0}\big|_{B(0, \varepsilon_0)}$ is a diffeomorphism onto an open neighbourhood $U_0$ of $p_0$, by [Exponential Map as a Local Diffeomorphism](/theorems/2712). Choose $\delta > 0$ small enough that $\gamma([t_0 - \delta, t_0 + \delta] \cap [0, 1]) \subseteq U_0$ and that the length of $\gamma$ on this subinterval is less than $\varepsilon_0$. The minimality and unit speed of $\gamma$ guarantee $\operatorname{length}(\gamma|_{[t_0, t_0 + \delta]}) = \delta = d(\gamma(t_0), \gamma(t_0 + \delta))$, so $\gamma|_{[t_0, t_0 + \delta]}$ is a length minimizer to a point in $U_0$ at distance $\delta < \varepsilon_0$. Apply [Geodesics Minimize Length Locally](/theorems/2720) part (3): the unique constant-speed length minimizer of length $< \varepsilon_0$ from $p_0$ to its image is the radial geodesic $t \mapsto \exp_{p_0}((t - t_0)\cdot v_0)$ for some $v_0 \in T_{p_0}M$ with $|v_0| = 1$. Hence $\gamma$ coincides with this smooth radial geodesic on $[t_0, t_0 + \delta]$. Since $t_0$ was arbitrary, $\gamma$ is locally a smooth geodesic, hence smooth and a geodesic on $[0, 1]$. [/proofplan]

[step:Reduce to the case of unit-speed parametrisation] By hypothesis, $\gamma$ is a piecewise $C^1$ minimal geodesic traversed at constant speed $c \ge 0$. (Here "minimal geodesic" means: $\operatorname{length}(\gamma) = d(\gamma(0), \gamma(1))$, and $\gamma$ minimises length among piecewise $C^1$ curves from $\gamma(0)$ to $\gamma(1)$.) If $c = 0$, then $\gamma$ is the constant curve at some point, which is $C^\infty$ and a (degenerate) geodesic. Assume henceforth $c > 0$. Reparametrize: define \begin{align*} \bar\gamma: [0, c] &\to M \\ s &\mapsto \gamma(s/c). \end{align*} Then $|\dot{\bar\gamma}(s)|_{g_{\bar\gamma(s)}} = (1/c)\, |\dot\gamma(s/c)|_{g_{\gamma(s/c)}} = (1/c) \cdot c = 1$, so $\bar\gamma$ has unit speed. The length and minimality are preserved (length is invariant under reparametrization; the image and minimality are unchanged). It suffices to prove $\bar\gamma$ is $C^\infty$ and a geodesic, since then $\gamma(t) = \bar\gamma(ct)$ is the affine reparametrization, which is also $C^\infty$ and a geodesic (the geodesic equation $\nabla_{\dot\gamma} \dot\gamma = 0$ is preserved under affine reparametrization). Henceforth, assume $\gamma: [0, L] \to M$ has unit speed, where $L = \operatorname{length}(\gamma) = d(\gamma(0), \gamma(L))$. The interval $[0, 1]$ in the theorem statement is replaced by $[0, L]$ to ease notation; the conclusion is unaffected. [/step]

[step:For each $t_0 \in [0, L]$, choose a normal-coordinate chart and a small subinterval] Fix $t_0 \in [0, L)$. We will show $\gamma$ is smooth on a closed interval $[t_0, t_0 + \delta]$ for some $\delta > 0$. Set $p_0 := \gamma(t_0) \in M$. By [Exponential Map as a Local Diffeomorphism](/theorems/2712), there exists $\varepsilon_0 > 0$ such that \begin{align*} \exp_{p_0}\big|_{B(0, \varepsilon_0)}: B(0, \varepsilon_0) \subseteq T_{p_0}M \to U_0 \subseteq M \end{align*} is a diffeomorphism onto an open neighbourhood $U_0$ of $p_0$. By continuity of $\gamma$ and openness of $U_0$, there exists $\delta_1 > 0$ such that $\gamma([t_0, t_0 + \delta_1] \cap [0, L]) \subseteq U_0$. Set \begin{align*} \delta := \min\bigl(\delta_1, \, \tfrac{1}{2}\varepsilon_0, \, L - t_0\bigr) > 0. \end{align*} Then $[t_0, t_0 + \delta] \subseteq [0, L]$, $\gamma([t_0, t_0 + \delta]) \subseteq U_0$, and the length of $\gamma|_{[t_0, t_0 + \delta]}$ — which equals $\delta$ since $\gamma$ has unit speed — satisfies $\delta \le \tfrac{1}{2}\varepsilon_0 < \varepsilon_0$. [/step]

[step:Show that $\gamma|_{[t_0, t_0 + \delta]}$ is a length minimizer of length $< \varepsilon_0$ from $p_0$ to $\gamma(t_0 + \delta)$] Set $q_0 := \gamma(t_0 + \delta) \in U_0$. We show $\operatorname{length}(\gamma|_{[t_0, t_0 + \delta]}) = d(p_0, q_0) = \delta$. Since $\gamma$ has unit speed, \begin{align*} \operatorname{length}(\gamma|_{[t_0, t_0 + \delta]}) = \int_{t_0}^{t_0 + \delta} |\dot\gamma(t)|\, d\mathcal{L}^1(t) = \int_{t_0}^{t_0 + \delta} 1\, d\mathcal{L}^1(t) = \delta. \end{align*} Now we show $d(p_0, q_0) = \delta$. The inequality $d(p_0, q_0) \le \delta$ holds because $\gamma|_{[t_0, t_0 + \delta]}$ is a piecewise $C^1$ curve from $p_0$ to $q_0$ of length $\delta$, and $d(p_0, q_0)$ is by definition the infimum of lengths of such curves. For the reverse inequality, suppose for contradiction that $d(p_0, q_0) < \delta$. Then there exists $\tilde\gamma_0 \in \Omega(p_0, q_0)$ (piecewise $C^1$) with $\operatorname{length}(\tilde\gamma_0) < \delta$. Construct a competitor curve from $\gamma(0)$ to $\gamma(L)$ by concatenating: \begin{align*} \tilde\gamma := \gamma|_{[0, t_0]} \cdot \tilde\gamma_0 \cdot \gamma|_{[t_0 + \delta, L]}, \end{align*} where $\cdot$ denotes the concatenation of piecewise $C^1$ curves (with parameter rescaling so the result is defined on $[0, L]$). This is piecewise $C^1$ and goes from $\gamma(0)$ to $\gamma(L)$, with total length \begin{align*} \operatorname{length}(\tilde\gamma) = t_0 + \operatorname{length}(\tilde\gamma_0) + (L - t_0 - \delta) < t_0 + \delta + (L - t_0 - \delta) = L. \end{align*} But this contradicts $L = d(\gamma(0), \gamma(L))$ (the minimality hypothesis on $\gamma$). Hence $d(p_0, q_0) \ge \delta$, and combined with $\le$, $d(p_0, q_0) = \delta$. We have established \begin{align*} \operatorname{length}(\gamma|_{[t_0, t_0 + \delta]}) = \delta = d(p_0, q_0) < \varepsilon_0. \end{align*} So $\gamma|_{[t_0, t_0 + \delta]}$ is a length-minimizing piecewise $C^1$ curve from $p_0$ to $q_0$, with length less than $\varepsilon_0$. [/step]

[step:Apply the local minimization theorem to identify $\gamma|_{[t_0, t_0 + \delta]}$ as a smooth radial geodesic] By [Geodesics Minimize Length Locally](/theorems/2720) part (1), applied with $p = p_0$ and $\varepsilon = \varepsilon_0$: there is a unique geodesic $\eta \in \Omega(p_0, q_0)$ with $\operatorname{length}(\eta) < \varepsilon_0$, and $\eta$ is the unique constant-speed length-minimizing curve from $p_0$ to $q_0$. We must verify that $\gamma|_{[t_0, t_0 + \delta]}$ has constant speed. By construction, $\gamma|_{[t_0, t_0 + \delta]}$ is a unit-speed reparametrisation of itself; since unit speed is constant, it has constant speed. By the uniqueness clause of part (1) of [Geodesics Minimize Length Locally](/theorems/2720), the constant-speed length-minimising curve from $p_0$ to $q_0$ on $[t_0, t_0 + \delta]$ (after affine reparametrisation to $[0, 1]$) is uniquely determined as the radial geodesic \begin{align*} s \mapsto \exp_{p_0}(s \cdot a_0), \quad s \in [0, 1], \end{align*} where $a_0 := (\exp_{p_0}\big|_{B(0, \varepsilon_0)})^{-1}(q_0) \in B(0, \varepsilon_0)$. Translating back to the original parametrisation on $[t_0, t_0 + \delta]$, \begin{align*} \gamma(t) = \exp_{p_0}\bigl((t - t_0) \cdot v_0\bigr) \quad \text{for } t \in [t_0, t_0 + \delta], \end{align*} where $v_0 := a_0/\delta \in T_{p_0}M$ (so that the affine rescaling $s = (t - t_0)/\delta$ converts $s \cdot a_0$ to $(t - t_0)\, v_0$). Note $|v_0|_{g_{p_0}} = |a_0|_{g_{p_0}}/\delta = \delta/\delta = 1$, so $v_0$ is a unit vector. The map $t \mapsto \exp_{p_0}((t - t_0)v_0)$ is $C^\infty$ on its domain — the exponential map is smooth on its domain of definition by the standard theory of smooth dependence of geodesics on initial data. Hence $\gamma$ is $C^\infty$ on $[t_0, t_0 + \delta]$, and being expressed as a radial geodesic, it satisfies $\nabla_{\dot\gamma}\dot\gamma = 0$ on $(t_0, t_0 + \delta)$. So $\gamma$ is a smooth geodesic on $[t_0, t_0 + \delta]$.[/step]

[guided]The application of [Geodesics Minimize Length Locally](/theorems/2720) part (1) is the heart of the proof. We have a curve $\gamma|_{[t_0, t_0 + \delta]}$ that achieves $\operatorname{length} = d(p_0, q_0) = \delta < \varepsilon_0$. The theorem then says: such a length-minimising curve, when reparametrised at constant speed, must be the unique radial geodesic from $p_0$ to $q_0$ in the normal chart. Why does $\gamma|_{[t_0, t_0 + \delta]}$ have constant speed? Because $\gamma$ has unit speed by the reduction in Step 1, and unit-speed curves have constant speed. Why does $\gamma|_{[t_0, t_0 + \delta]}$ minimise length among piecewise $C^1$ curves from $p_0$ to $q_0$, not just have a length equal to $d(p_0, q_0)$? These are the same statement: $d(p_0, q_0)$ is, by definition, the infimum of lengths of piecewise $C^1$ curves from $p_0$ to $q_0$. A curve realising this infimum is a length minimiser. The conclusion that $\gamma$ is smooth on $[t_0, t_0 + \delta]$ uses the fact that $\exp_{p_0}$ is $C^\infty$, established when constructing the exponential map from the geodesic flow. The composition of a $C^\infty$ map with the linear (hence $C^\infty$) parametrisation $t \mapsto (t - t_0) v_0$ is $C^\infty$, so $\gamma$ is $C^\infty$ on the chosen subinterval. The last step — verifying the geodesic equation $\nabla_{\dot\gamma} \dot\gamma = 0$ — follows because radial geodesics from a point are by definition geodesics, satisfying the geodesic equation in the normal chart.[/guided]

[step:Patch the local conclusions to obtain $C^\infty$ and the geodesic equation on all of $[0, L]$] For each $t_0 \in [0, L)$, the previous step gives an open interval $(t_0 - \delta, t_0 + \delta) \cap [0, L]$ (using the same construction with $\delta$ extended to a two-sided interval; the same argument applies for $t_0 \in (0, L]$ on $[t_0 - \delta, t_0]$) on which $\gamma$ is $C^\infty$ and satisfies $\nabla_{\dot\gamma}\dot\gamma = 0$. The collection of these open intervals covers $[0, L]$. Hence $\gamma$ is $C^\infty$ on $[0, L]$ and satisfies the geodesic equation everywhere. To make the patching explicit: $\gamma$ being $C^\infty$ and a geodesic on each open subinterval of a cover means $\gamma$ is $C^\infty$ on the union $[0, L]$ (smoothness is local) and satisfies the geodesic equation pointwise everywhere on $[0, L]$ (each point lies in some open subinterval where the equation holds, and the equation is pointwise). The endpoint behaviour: at $t = 0$ and $t = L$, smoothness in the one-sided sense follows from smoothness on a one-sided open neighbourhood by the same argument applied to the one-sided radial geodesic chart. The geodesic equation is satisfied at the endpoints by continuity of $\nabla_{\dot\gamma}\dot\gamma$ (which is continuous because $\gamma$ is $C^\infty$, so $\dot\gamma$ is continuous, and $\nabla$ acts smoothly). We have proved: the unit-speed reparametrisation $\bar\gamma$ of the original $\gamma$ is $C^\infty$ on its full parameter interval and is a geodesic. Reparametrising back via the affine map $t \mapsto t/c$ (with $c$ the original constant speed), we conclude that the original $\gamma: [0, 1] \to M$ is $C^\infty$ and a geodesic.[/step]

[guided]The patching argument is local-to-global. We have shown that for every parameter $t_0$, $\gamma$ is smooth on a neighbourhood of $t_0$ and satisfies the geodesic equation there. Two facts then complete the proof: 1. **Smoothness is a local property.** A function is $C^\infty$ on a set if and only if it is $C^\infty$ on a neighbourhood of each point in the set. So smoothness on each neighbourhood implies smoothness on the union. 2. **The geodesic equation is pointwise.** At each $t \in [0, L]$, the equation $\nabla_{\dot\gamma}\dot\gamma|_t = 0$ refers to the value of a vector field at a single point. We have established this value is zero on a neighbourhood of $t$, hence at $t$. Both arguments use that the open neighbourhoods constructed in Step 4 form an open cover of the compact $[0, L]$, but in fact the local-to-global passage does not need compactness — pointwise validity on a neighbourhood of every point is exactly equivalent to pointwise validity everywhere. A subtle issue is whether the smoothness improvement is automatic given the geodesic equation. The answer is yes: any solution of the geodesic equation, even a priori only $C^2$, is automatically $C^\infty$ by the smoothness of the connection coefficients (the Christoffel symbols are smooth, so a $C^2$ solution of $\ddot x_i + \Gamma_{jk}^i(\gamma)\dot x_j \dot x_k = 0$ inherits regularity from the right-hand side). This is a standard result on smooth dependence in ODE theory.[/guided]