[proofplan]
The regularity hypothesis says that the equation $p_n=\ell$ can be solved locally for the cyclic-coordinate velocity $w=\dot q_n$ by the implicit function theorem. The Routhian is then the partial Legendre transform $L-\ell w$ restricted to this local momentum level. Differentiating $R_\ell$ shows that all derivatives of the branch $\psi_\ell$ cancel because $\partial_wL=\ell$ on the level set. Therefore the Euler-Lagrange equations for $R_\ell$ are exactly the non-cyclic Euler-Lagrange equations for $L$ restricted to $p_n=\ell$, while the cyclic equation is precisely conservation of $p_n$.
[/proofplan]
[step:Solve the fixed momentum equation for the cyclic velocity]
Define the map $F:O\to\mathbb{R}$ by
\begin{align*}
F(x,r,v,w)=p_n(x,r,v,w)-\ell.
\end{align*}
Since $L$ is $C^2$ and $p_n=\partial_wL$, the momentum component $p_n:O\to\mathbb{R}$ is $C^1$, so $F$ is $C^1$. Then $F(x_0,r_0,v_0,w_0)=0$ and
\begin{align*}
\partial_wF(x_0,r_0,v_0,w_0)=\partial_w p_n(x_0,r_0,v_0,w_0)\neq 0.
\end{align*}
By the [Implicit Function Theorem](/theorems/52), after shrinking to neighbourhoods $A$ of $(x_0,v_0)$, $J$ of $r_0$, and $B$ of $w_0$ so that $A\times J\times B\subset O$, there is a unique $C^1$ map $\phi_\ell:A\times J\to B$ such that
\begin{align*}
F(x,r,v,w)=0 \quad \iff \quad w=\phi_\ell(x,v,r)
\end{align*}
for all $(x,r,v,w)\in A\times J\times B$. Since $\partial_rL=0$ on the coordinate neighbourhood and $L\in C^2$, differentiating the identity $\partial_rL=0$ with respect to $w$ and using equality of mixed partial derivatives gives
\begin{align*}
\partial_rp_n=\partial_r\partial_wL=\partial_w\partial_rL=0.
\end{align*}
Because $J$ is an interval, $r\mapsto p_n(x,r,v,w)$ is constant on $J$ for each fixed $(x,v,w)\in A\times B$. Thus, for fixed $(x,v)\in A$, the same $w\in B$ solves $F(x,r,v,w)=0$ for every $r\in J$; the local uniqueness in the implicit function theorem forces $\phi_\ell(x,v,r)$ to be independent of $r$. Define $\psi_\ell:A\to B$ by $\psi_\ell(x,v)=\phi_\ell(x,v,r_0)$. Then
\begin{align*}
F(x,r,v,w)=0 \quad \iff \quad w=\psi_\ell(x,v)
\end{align*}
for all $(x,r,v,w)\in A\times J\times B$.
[guided]
We first isolate the equation that defines the momentum level. Define $F:O\to\mathbb{R}$ by
\begin{align*}
F(x,r,v,w)=p_n(x,r,v,w)-\ell.
\end{align*}
Because $L$ is $C^2$ and $p_n=\partial_wL$, the momentum component $p_n:O\to\mathbb{R}$ is $C^1$. Hence $F$ is also $C^1$, as required for the implicit function theorem. The point $(x_0,r_0,v_0,w_0)$ lies on the level set $p_n=\ell$, so
\begin{align*}
F(x_0,r_0,v_0,w_0)=0.
\end{align*}
The derivative of $F$ in the unknown variable $w=\dot q_n$ is
\begin{align*}
\partial_wF(x_0,r_0,v_0,w_0)=\partial_w p_n(x_0,r_0,v_0,w_0)=\partial_w^2L(x_0,r_0,v_0,w_0).
\end{align*}
This number is nonzero by hypothesis. Therefore the [Implicit Function Theorem](/theorems/52) applies to the equation $F=0$ with $w$ as the solved variable. It gives neighbourhoods $A\subset\mathbb{R}^{n-1}\times\mathbb{R}^{n-1}$, $J\subset\mathbb{R}$, and $B\subset\mathbb{R}$, shrunk so that $A\times J\times B\subset O$, together with a unique $C^1$ function $\phi_\ell:A\times J\to B$, such that
\begin{align*}
p_n(x,r,v,w)=\ell \quad \iff \quad w=\phi_\ell(x,v,r)
\end{align*}
for all $(x,r,v,w)\in A\times J\times B$.
The next point is to prove that this branch actually has no dependence on $r$. Because $q_n$ is a cyclic coordinate, $\partial_rL=0$ throughout the coordinate neighbourhood. Since $L\in C^2$, the mixed second partial derivatives in the $r$ and $w$ variables commute. Differentiating the identity $\partial_rL=0$ with respect to $w$ therefore gives
\begin{align*}
\partial_rp_n=\partial_r\partial_wL=\partial_w\partial_rL=0.
\end{align*}
Since $J$ is an interval, the identity $\partial_rp_n=0$ implies that the function $r\mapsto p_n(x,r,v,w)$ is constant on $J$ for every fixed $(x,v,w)\in A\times B$. Thus, for fixed $(x,v)\in A$, if a value $w\in B$ solves $p_n(x,r_1,v,w)=\ell$ at one $r_1\in J$, then the same $w$ solves $p_n(x,r_2,v,w)=\ell$ for every $r_2\in J$. The implicit-function theorem also gives local uniqueness of the solution in $B$, so $\phi_\ell(x,v,r_1)=\phi_\ell(x,v,r_2)$ for all $r_1,r_2\in J$. Therefore the branch is independent of $r$. Define $\psi_\ell:A\to B$ by $\psi_\ell(x,v)=\phi_\ell(x,v,r_0)$. Then
\begin{align*}
p_n(x,r,v,w)=\ell \quad \iff \quad w=\psi_\ell(x,v)
\end{align*}
for all $(x,r,v,w)\in A\times J\times B$.
[/guided]
[/step]
[step:Compute the derivatives of the Routhian and cancel the branch terms]
If $n=1$, there are no non-cyclic coordinate derivatives to compute in this step. Otherwise fix $i\in\{1,\dots,n-1\}$. Define the local Routhian $R_\ell:A\to\mathbb{R}$ by
\begin{align*}
R_\ell(x,v)=L(x,r,v,\psi_\ell(x,v))-\ell\,\psi_\ell(x,v).
\end{align*}
The expression is independent of the choice of $r\in J$ because $\partial_rL=0$. Since $p_n(x,r,v,\psi_\ell(x,v))=\ell$, the chain rule gives
\begin{align*}
\partial_{x_i}R_\ell(x,v)=\partial_{x_i}L(x,r,v,\psi_\ell(x,v)).
\end{align*}
Indeed, the terms involving $\partial_{x_i}\psi_\ell$ combine to
\begin{align*}
\partial_wL(x,r,v,\psi_\ell(x,v))\,\partial_{x_i}\psi_\ell(x,v)-\ell\,\partial_{x_i}\psi_\ell(x,v)=0.
\end{align*}
The same computation in the velocity variable gives
\begin{align*}
\partial_{v_i}R_\ell(x,v)=\partial_{v_i}L(x,r,v,\psi_\ell(x,v)).
\end{align*}
[guided]
The purpose of this computation is to show that differentiating after restricting to the momentum level introduces no extra force terms. If $n=1$, there are no non-cyclic coordinate derivatives to compute. Otherwise fix $i\in\{1,\dots,n-1\}$ and define the local Routhian $R_\ell:A\to\mathbb{R}$ by
\begin{align*}
R_\ell(x,v)=L(x,r,v,\psi_\ell(x,v))-\ell\,\psi_\ell(x,v).
\end{align*}
This definition is independent of the selected $r\in J$ because $q_n$ is cyclic, so $\partial_rL=0$, and the branch $\psi_\ell$ is independent of $r$ by the previous step.
We now differentiate with respect to $x_i$. The chain rule applies because $L$ is $C^2$ and $\psi_\ell$ is $C^1$. It gives
\begin{align*}
\partial_{x_i}R_\ell(x,v)=\partial_{x_i}L(x,r,v,\psi_\ell(x,v))+\partial_wL(x,r,v,\psi_\ell(x,v))\,\partial_{x_i}\psi_\ell(x,v)-\ell\,\partial_{x_i}\psi_\ell(x,v).
\end{align*}
On the chosen momentum level, $p_n(x,r,v,\psi_\ell(x,v))=\ell$, and by definition $p_n=\partial_wL$. Hence the two terms involving $\partial_{x_i}\psi_\ell$ combine to
\begin{align*}
\partial_wL(x,r,v,\psi_\ell(x,v))\,\partial_{x_i}\psi_\ell(x,v)-\ell\,\partial_{x_i}\psi_\ell(x,v)=0.
\end{align*}
Therefore
\begin{align*}
\partial_{x_i}R_\ell(x,v)=\partial_{x_i}L(x,r,v,\psi_\ell(x,v)).
\end{align*}
The same chain-rule computation in the velocity variable gives
\begin{align*}
\partial_{v_i}R_\ell(x,v)=\partial_{v_i}L(x,r,v,\psi_\ell(x,v))+\partial_wL(x,r,v,\psi_\ell(x,v))\,\partial_{v_i}\psi_\ell(x,v)-\ell\,\partial_{v_i}\psi_\ell(x,v).
\end{align*}
Again $\partial_wL=\ell$ on the branch, so the terms involving $\partial_{v_i}\psi_\ell$ cancel and
\begin{align*}
\partial_{v_i}R_\ell(x,v)=\partial_{v_i}L(x,r,v,\psi_\ell(x,v)).
\end{align*}
[/guided]
[/step]
[step:Project full solutions to the reduced Euler-Lagrange equations]
Let $q:I\to\mathbb{R}^n$, $q(t)=(x(t),r(t))$, be a $C^2$ solution of the Euler-Lagrange equations for $L$ with
\begin{align*}
p_n(x(t),r(t),\dot x(t),\dot r(t))=\ell
\end{align*}
for all $t\in I$. By the defining property of $\psi_\ell$,
\begin{align*}
\dot r(t)=\psi_\ell(x(t),\dot x(t)).
\end{align*}
Using the derivative identities from the previous step along this curve, for each $i\in\{1,\dots,n-1\}$ we have
\begin{align*}
\partial_{x_i}R_\ell(x(t),\dot x(t))=\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t))
\end{align*}
and
\begin{align*}
\partial_{v_i}R_\ell(x(t),\dot x(t))=\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
The Euler-Lagrange equation for $L$ in the $q_i$ coordinate gives
\begin{align*}
\frac{d}{dt}\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
Substitution yields
\begin{align*}
\frac{d}{dt}\partial_{v_i}R_\ell(x(t),\dot x(t))=\partial_{x_i}R_\ell(x(t),\dot x(t)).
\end{align*}
Thus $x$ satisfies the Euler-Lagrange equations for $R_\ell$.
[guided]
We start with a full trajectory $q:I\to\mathbb{R}^n$ written as $q(t)=(x(t),r(t))$, where $x:I\to\mathbb{R}^{n-1}$ denotes the non-cyclic coordinates and $r:I\to\mathbb{R}$ denotes the cyclic coordinate. Assume $q$ is a $C^2$ solution of the Euler-Lagrange equations for $L$ and that it stays on the momentum level:
\begin{align*}
p_n(x(t),r(t),\dot x(t),\dot r(t))=\ell
\end{align*}
for every $t\in I$.
Because the curve lies in the local branch domain and $p_n=\ell$, the defining property of $\psi_\ell$ gives
\begin{align*}
\dot r(t)=\psi_\ell(x(t),\dot x(t)).
\end{align*}
This is the point where the full trajectory is restricted to the solved momentum level. Along this branch, the derivative identities for the Routhian become, for each $i\in\{1,\dots,n-1\}$,
\begin{align*}
\partial_{x_i}R_\ell(x(t),\dot x(t))=\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t))
\end{align*}
and
\begin{align*}
\partial_{v_i}R_\ell(x(t),\dot x(t))=\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
The full Euler-Lagrange equation in the non-cyclic coordinate $q_i$ states that
\begin{align*}
\frac{d}{dt}\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
Substituting the two Routhian derivative identities into this equation gives
\begin{align*}
\frac{d}{dt}\partial_{v_i}R_\ell(x(t),\dot x(t))=\partial_{x_i}R_\ell(x(t),\dot x(t)).
\end{align*}
Since this holds for every $i\in\{1,\dots,n-1\}$, the projected curve $x$ satisfies the Euler-Lagrange equations for $R_\ell$.
[/guided]
[/step]
[step:Reconstruct a full solution from a reduced solution]
Conversely, let $x:I\to\mathbb{R}^{n-1}$ be a $C^2$ solution of the Euler-Lagrange equations for $R_\ell$ with $(x(t),\dot x(t))\in A$. Fix $t_0\in I$ and $r_0'\in J$. Since $x$ is $C^2$ and $\psi_\ell$ is $C^1$, the right-hand side $t\mapsto\psi_\ell(x(t),\dot x(t))$ is $C^1$, so the scalar reconstruction equation has a local $C^2$ solution. On any subinterval of $I$ containing $t_0$ on which this solution exists and remains in $J$, let $r$ denote the resulting map into $J$ and write
\begin{align*}
\dot r(t)=\psi_\ell(x(t),\dot x(t)), \qquad r(t_0)=r_0'.
\end{align*}
Define $q:I\to\mathbb{R}^n$ by $q(t)=(x(t),r(t))$. Since $x$ is $C^2$ and $\psi_\ell$ is $C^1$, the map $t\mapsto\psi_\ell(x(t),\dot x(t))$ is $C^1$; hence the solution $r$ is $C^2$ and $q$ is a $C^2$ curve. Then
\begin{align*}
p_n(x(t),r(t),\dot x(t),\dot r(t))=\ell
\end{align*}
because $\dot r(t)=\psi_\ell(x(t),\dot x(t))$.
For $i\in\{1,\dots,n-1\}$, the derivative identities for $R_\ell$ give
\begin{align*}
\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{v_i}R_\ell(x(t),\dot x(t))
\end{align*}
and
\begin{align*}
\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{x_i}R_\ell(x(t),\dot x(t)).
\end{align*}
Since $x$ satisfies the Euler-Lagrange equations for $R_\ell$, it follows that
\begin{align*}
\frac{d}{dt}\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
Thus the non-cyclic Euler-Lagrange equations for $L$ hold.
For the cyclic coordinate, the Euler-Lagrange equation is
\begin{align*}
\frac{d}{dt}\partial_wL(x(t),r(t),\dot x(t),\dot r(t))=\partial_rL(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
The left-hand side is $\frac{d}{dt}\ell=0$, and the right-hand side is $0$ because $q_n$ is cyclic. Hence the cyclic Euler-Lagrange equation also holds. Therefore $q=(x,r)$ solves the full Euler-Lagrange equations for $L$ and remains on the momentum level $p_n=\ell$ for as long as the chosen coordinate chart and branch are defined.
[guided]
Now we reverse the construction. Let $x:I\to\mathbb{R}^{n-1}$ be a $C^2$ solution of the Euler-Lagrange equations for $R_\ell$ with $(x(t),\dot x(t))\in A$ for all $t\in I$. Choose an initial time $t_0\in I$ and an initial cyclic coordinate $r_0'\in J$. Because $x$ is $C^2$ and $\psi_\ell$ is $C^1$, the function $t\mapsto\psi_\ell(x(t),\dot x(t))$ is $C^1$. Hence the scalar ordinary differential equation
\begin{align*}
\dot r(t)=\psi_\ell(x(t),\dot x(t)), \qquad r(t_0)=r_0'
\end{align*}
has a local $C^2$ solution. On any subinterval of $I$ containing $t_0$ on which this solution exists and remains in $J$, write this solution as $r$ mapping that subinterval into $J$. The curve $q$ defined by $q(t)=(x(t),r(t))$ is then $C^2$ on that subinterval.
The reconstruction equation gives $\dot r(t)=\psi_\ell(x(t),\dot x(t))$. By the defining property of the branch $\psi_\ell$, this implies
\begin{align*}
p_n(x(t),r(t),\dot x(t),\dot r(t))=\ell
\end{align*}
for every $t$ for which the branch is defined.
We next verify the non-cyclic Euler-Lagrange equations for $L$. For each $i\in\{1,\dots,n-1\}$, the derivative identities for $R_\ell$ give
\begin{align*}
\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{v_i}R_\ell(x(t),\dot x(t))
\end{align*}
and
\begin{align*}
\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{x_i}R_\ell(x(t),\dot x(t)).
\end{align*}
Since $x$ solves the Euler-Lagrange equations for $R_\ell$,
\begin{align*}
\frac{d}{dt}\partial_{v_i}R_\ell(x(t),\dot x(t))=\partial_{x_i}R_\ell(x(t),\dot x(t)).
\end{align*}
Substituting the derivative identities into this equality yields
\begin{align*}
\frac{d}{dt}\partial_{\dot q_i}L(x(t),r(t),\dot x(t),\dot r(t))=\partial_{x_i}L(x(t),r(t),\dot x(t),\dot r(t)).
\end{align*}
Thus all non-cyclic Euler-Lagrange equations for $L$ hold.
It remains to check the cyclic equation. Since $p_n=\partial_wL$ and the reconstructed curve satisfies $p_n=\ell$, the left-hand side of the cyclic Euler-Lagrange equation is
\begin{align*}
\frac{d}{dt}\partial_wL(x(t),r(t),\dot x(t),\dot r(t))=\frac{d}{dt}\ell=0.
\end{align*}
The right-hand side is
\begin{align*}
\partial_rL(x(t),r(t),\dot x(t),\dot r(t))=0
\end{align*}
because $q_n$ is cyclic. Therefore the cyclic Euler-Lagrange equation also holds, and $q=(x,r)$ is a full local solution on the fixed momentum level for as long as the selected chart and branch remain defined.
[/guided]
[/step]