[proofplan]
We expand the total derivative in the Euler-Lagrange equation using the chain rule in the local coordinates $(q,v,t)$. The coefficients of the acceleration components $\ddot q_j$ are exactly the entries of the fiber Hessian with respect to the velocity variables. All remaining terms depend only on $(q,\dot q,t)$ and derivatives of $\ell$, so they form the vector $F$. If the fiber Hessian is invertible, the resulting linear system has the unique solution $a=H^{-1}F$, and smoothness follows from the adjugate formula for the inverse matrix.
[/proofplan]
[step:Expand the Euler-Lagrange derivative in local coordinates]
Let $\ell:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}$ denote the local representative of $L$ in the chart $(U,\varphi)$. Fix a smooth curve $q:I\to \varphi(U)$, where $I\subset \mathbb{R}$ is an interval. Define its velocity map $v:I\to \mathbb{R}^n$ by $v(t)=\dot q(t)$ for $t\in I$. Define its acceleration map $a:I\to \mathbb{R}^n$ by $a(t)=\ddot q(t)$ for $t\in I$.
For each $i\in\{1,\dots,n\}$, define the smooth function
\begin{align*}
G_i:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}
\end{align*}
by setting $G_i(q,v,t)=\frac{\partial \ell}{\partial v_i}(q,v,t)$. Applying the coordinate form of the [Chain Rule for Maps Between Euclidean Spaces](/theorems/323) to the composite map $t\mapsto G_i(q(t),v(t),t)$ gives
\begin{align*}
\frac{d}{dt}\frac{\partial \ell}{\partial v_i}(q(t),\dot q(t),t)=\sum_{k=1}^n \frac{\partial^2 \ell}{\partial q_k\partial v_i}(q(t),\dot q(t),t)\dot q_k(t)+\sum_{j=1}^n \frac{\partial^2 \ell}{\partial v_j\partial v_i}(q(t),\dot q(t),t)\ddot q_j(t)+\frac{\partial^2 \ell}{\partial t\partial v_i}(q(t),\dot q(t),t).
\end{align*}
[guided]
Fix a smooth curve $q:I\to \varphi(U)$. Its velocity is the map $v:I\to\mathbb{R}^n$ defined by $v(t)=\dot q(t)$ for $t\in I$, and its acceleration is the map $a:I\to\mathbb{R}^n$ defined by $a(t)=\ddot q(t)$ for $t\in I$. In coordinates, the Euler-Lagrange expression differentiates the function $\partial \ell/\partial v_i$ along the curve $t\mapsto (q(t),v(t),t)$.
For each $i\in\{1,\dots,n\}$, define
\begin{align*}
G_i:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}
\end{align*}
by setting $G_i(q,v,t)=\frac{\partial \ell}{\partial v_i}(q,v,t)$. The composite map $t\mapsto G_i(q(t),v(t),t)$ has three kinds of dependence on $t$: through $q(t)$, through $v(t)$, and through the explicit time coordinate $t$. Applying the coordinate form of the [Chain Rule for Maps Between Euclidean Spaces](/theorems/323) therefore gives
\begin{align*}
\frac{d}{dt}G_i(q(t),v(t),t)=\sum_{k=1}^n \frac{\partial G_i}{\partial q_k}(q(t),v(t),t)\dot q_k(t)+\sum_{j=1}^n \frac{\partial G_i}{\partial v_j}(q(t),v(t),t)\dot v_j(t)+\frac{\partial G_i}{\partial t}(q(t),v(t),t).
\end{align*}
Since $G_i=\partial \ell/\partial v_i$ and $\dot v_j(t)=\ddot q_j(t)$, this becomes
\begin{align*}
\frac{d}{dt}\frac{\partial \ell}{\partial v_i}(q(t),\dot q(t),t)=\sum_{k=1}^n \frac{\partial^2 \ell}{\partial q_k\partial v_i}(q(t),\dot q(t),t)\dot q_k(t)+\sum_{j=1}^n \frac{\partial^2 \ell}{\partial v_j\partial v_i}(q(t),\dot q(t),t)\ddot q_j(t)+\frac{\partial^2 \ell}{\partial t\partial v_i}(q(t),\dot q(t),t).
\end{align*}
This computation isolates exactly where the acceleration appears: only the derivatives with respect to the velocity variables contribute factors of $\ddot q_j(t)$.
[/guided]
[/step]
[step:Collect the non-acceleration terms into the force vector]
Define the fiber Hessian map
\begin{align*}
H:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}^{n\times n}
\end{align*}
by setting $H(q,v,t)_{ij}=\frac{\partial^2 \ell}{\partial v_j\partial v_i}(q,v,t)$.
Substitute the expansion into the Euler-Lagrange equation
\begin{align*}
\frac{d}{dt}\frac{\partial \ell}{\partial v_i}(q(t),\dot q(t),t)-\frac{\partial \ell}{\partial q_i}(q(t),\dot q(t),t)=0.
\end{align*}
Moving all terms not involving $\ddot q$ to the right-hand side gives
\begin{align*}
\sum_{j=1}^n \frac{\partial^2 \ell}{\partial v_j\partial v_i}(q(t),\dot q(t),t)\ddot q_j(t)=\frac{\partial \ell}{\partial q_i}(q(t),\dot q(t),t)-\sum_{k=1}^n \frac{\partial^2 \ell}{\partial q_k\partial v_i}(q(t),\dot q(t),t)\dot q_k(t)-\frac{\partial^2 \ell}{\partial t\partial v_i}(q(t),\dot q(t),t).
\end{align*}
For each $i\in\{1,\dots,n\}$, define the smooth function
\begin{align*}
F_i:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}
\end{align*}
by setting
\begin{align*}
F_i(q,v,t)=\frac{\partial \ell}{\partial q_i}(q,v,t)-\sum_{k=1}^n \frac{\partial^2 \ell}{\partial q_k\partial v_i}(q,v,t)v_k-\frac{\partial^2 \ell}{\partial t\partial v_i}(q,v,t).
\end{align*}
Thus the coefficient matrix of the acceleration vector is the formally declared fiber Hessian $H$, and the right-hand side is exactly $F_i(q(t),\dot q(t),t)$.
[/step]
[step:Solve the fiber Hessian linear system where the Hessian is invertible]
Define the smooth vector-valued map
\begin{align*}
F:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}^n
\end{align*}
by setting $F(q,v,t)=(F_1(q,v,t),\dots,F_n(q,v,t))$.
At a point $(q,v,t)$ where $H(q,v,t)$ is invertible, the Euler-Lagrange equations form the linear system
\begin{align*}
H(q,v,t)a=F(q,v,t)
\end{align*}
for the unknown acceleration vector $a\in\mathbb{R}^n$. Since $H(q,v,t)$ is invertible, this system has the unique solution
\begin{align*}
a=H(q,v,t)^{-1}F(q,v,t).
\end{align*}
Therefore the Euler-Lagrange equations determine $\ddot q$ uniquely at every point where the fiber Hessian is invertible.
[guided]
We now rewrite the collected Euler-Lagrange equations as a finite-dimensional linear system. The unknown is the acceleration vector $a=(a_1,\dots,a_n)\in\mathbb{R}^n$, where $a_j$ represents the value of $\ddot q_j$ at the fixed position, velocity, and time $(q,v,t)$. The coefficient matrix is the fiber Hessian map
\begin{align*}
H:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}^{n\times n}
\end{align*}
defined by $H(q,v,t)_{ij}=\frac{\partial^2 \ell}{\partial v_j\partial v_i}(q,v,t)$. The non-acceleration terms are collected into the smooth functions
\begin{align*}
F_i:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}
\end{align*}
defined by
\begin{align*}
F_i(q,v,t)=\frac{\partial \ell}{\partial q_i}(q,v,t)-\sum_{k=1}^n \frac{\partial^2 \ell}{\partial q_k\partial v_i}(q,v,t)v_k-\frac{\partial^2 \ell}{\partial t\partial v_i}(q,v,t).
\end{align*}
Define the vector-valued map
\begin{align*}
F:\varphi(U)\times\mathbb{R}^n\times\mathbb{R}\to\mathbb{R}^n
\end{align*}
by setting $F(q,v,t)=(F_1(q,v,t),\dots,F_n(q,v,t))$.
To see independently why these definitions give the collected system, start from the $i$th Euler-Lagrange equation
\begin{align*}
\frac{d}{dt}\frac{\partial \ell}{\partial v_i}(q(t),\dot q(t),t)-\frac{\partial \ell}{\partial q_i}(q(t),\dot q(t),t)=0.
\end{align*}
By the [Chain Rule for Maps Between Euclidean Spaces](/theorems/323), applied to the composite map $t\mapsto G_i(q(t),v(t),t)$ with $G_i=\partial \ell/\partial v_i$, we have
\begin{align*}
\frac{d}{dt}\frac{\partial \ell}{\partial v_i}(q(t),\dot q(t),t)=\sum_{k=1}^n \frac{\partial^2 \ell}{\partial q_k\partial v_i}(q(t),\dot q(t),t)\dot q_k(t)+\sum_{j=1}^n \frac{\partial^2 \ell}{\partial v_j\partial v_i}(q(t),\dot q(t),t)\ddot q_j(t)+\frac{\partial^2 \ell}{\partial t\partial v_i}(q(t),\dot q(t),t).
\end{align*}
Substituting this full expansion into the Euler-Lagrange equation and moving every term not involving $\ddot q$ to the right-hand side gives
\begin{align*}
\sum_{j=1}^n \frac{\partial^2 \ell}{\partial v_j\partial v_i}(q(t),\dot q(t),t)\ddot q_j(t)=F_i(q(t),\dot q(t),t).
\end{align*}
By the definition of $H$, the coefficient of $\ddot q_j(t)$ is $H(q(t),\dot q(t),t)_{ij}$. With these definitions, the $i$th Euler-Lagrange equation is
\begin{align*}
\sum_{j=1}^n H(q,v,t)_{ij}a_j=F_i(q,v,t).
\end{align*}
Equivalently, the full system is
\begin{align*}
H(q,v,t)a=F(q,v,t).
\end{align*}
The hypothesis says precisely that the matrix $H(q,v,t)$ is invertible. Hence multiplication by $H(q,v,t)$ is a bijective [linear map](/page/Linear%20Map) from $\mathbb{R}^n$ to $\mathbb{R}^n$, and the unique vector solving the system is
\begin{align*}
a=H(q,v,t)^{-1}F(q,v,t).
\end{align*}
This is the point of the fiber Hessian condition: it rules out both non-uniqueness of accelerations and failure of existence for the acceleration determined by a given position, velocity, and time.
[/guided]
[/step]
[step:Prove smooth dependence of the acceleration on position velocity and time]
Fix a point $(q_0,v_0,t_0)\in\varphi(U)\times\mathbb{R}^n\times\mathbb{R}$ at which $H(q_0,v_0,t_0)$ is invertible. The determinant $\det H(q,v,t)$ is a polynomial expression in the entries of $H(q,v,t)$. Since the entries of $H$ are smooth, the composite function $(q,v,t)\mapsto \det H(q,v,t)$ is smooth, hence continuous. Because $\det H(q_0,v_0,t_0)\neq 0$, there exists an open neighbourhood $\Omega\subset \varphi(U)\times\mathbb{R}^n\times\mathbb{R}$ of $(q_0,v_0,t_0)$ such that $\det H(q,v,t)\neq 0$ for every $(q,v,t)\in\Omega$.
The entries of $H$ and $F$ are smooth on $\Omega$ because $\ell$ is smooth. Define
\begin{align*}
A:\Omega\to\mathbb{R}^n
\end{align*}
by setting $A(q,v,t)=H(q,v,t)^{-1}F(q,v,t)$.
Let $\operatorname{adj}:\mathbb{R}^{n\times n}\to\mathbb{R}^{n\times n}$ denote the adjugate map, whose $(i,j)$-entry is the $(j,i)$-cofactor of the input matrix. For each $(q,v,t)\in\Omega$, the [Adjugate Identity](/theorems/397) gives the inverse matrix formula
\begin{align*}
H(q,v,t)^{-1}=\frac{\operatorname{adj}(H(q,v,t))}{\det H(q,v,t)}.
\end{align*}
The entries of $\operatorname{adj}(H)$ and $\det H$ are polynomial functions of the entries of $H$, hence are smooth functions of $(q,v,t)$. Since $\det H(q,v,t)\neq 0$ on $\Omega$, division by $\det H$ preserves smoothness. Therefore $H^{-1}$ is smooth on $\Omega$, and multiplying the smooth matrix-valued map $H^{-1}$ by the smooth vector-valued map $F$ shows that $A$ is smooth. This proves that the Euler-Lagrange equations determine the acceleration uniquely as a smooth function on the open regular locus where the fiber Hessian is invertible.
[guided]
Fix a point $(q_0,v_0,t_0)\in\varphi(U)\times\mathbb{R}^n\times\mathbb{R}$ where the fiber Hessian $H(q_0,v_0,t_0)$ is invertible. Invertibility is an open condition for matrices, and here it follows directly from the determinant. Since $H$ is smooth, each entry of $H$ is smooth. The determinant is a polynomial expression in the matrix entries, so the function $(q,v,t)\mapsto\det H(q,v,t)$ is smooth and therefore continuous. At the chosen point, invertibility means
\begin{align*}
\det H(q_0,v_0,t_0)\neq 0.
\end{align*}
By continuity of the determinant function, there is an open neighbourhood $\Omega\subset \varphi(U)\times\mathbb{R}^n\times\mathbb{R}$ of $(q_0,v_0,t_0)$ such that
\begin{align*}
\det H(q,v,t)\neq 0
\end{align*}
for every $(q,v,t)\in\Omega$. Thus $H(q,v,t)$ is invertible throughout $\Omega$, not merely at the single point.
Now define the acceleration map on this neighbourhood by
\begin{align*}
A:\Omega\to\mathbb{R}^n
\end{align*}
with $A(q,v,t)=H(q,v,t)^{-1}F(q,v,t)$.
We must justify that this map is smooth. Let $\operatorname{adj}:\mathbb{R}^{n\times n}\to\mathbb{R}^{n\times n}$ denote the adjugate map, whose $(i,j)$-entry is the $(j,i)$-cofactor of the input matrix. For every $(q,v,t)\in\Omega$, the [Adjugate Identity](/theorems/397) gives the inverse matrix formula
\begin{align*}
H(q,v,t)^{-1}=\frac{\operatorname{adj}(H(q,v,t))}{\det H(q,v,t)}.
\end{align*}
The entries of $\operatorname{adj}(H(q,v,t))$ are polynomial expressions in the entries of $H(q,v,t)$, and $\det H(q,v,t)$ is also a polynomial expression in those entries. Since the entries of $H$ are smooth functions of $(q,v,t)$, both $\operatorname{adj}(H)$ and $\det H$ vary smoothly on $\Omega$. The denominator never vanishes on $\Omega$, so division by $\det H$ preserves smoothness. Therefore $H^{-1}:\Omega\to\mathbb{R}^{n\times n}$ is smooth.
Finally, $F:\Omega\to\mathbb{R}^n$ is smooth because its components are built from first and second partial derivatives of the smooth local Lagrangian $\ell$ and the coordinate functions $v_k$. Matrix-vector multiplication is a smooth bilinear operation, so the product $H^{-1}F$ is smooth. Hence $A$ is a smooth function of $(q,v,t)$, and substituting $v=\dot q$ shows that the Euler-Lagrange equations determine $\ddot q$ uniquely as a smooth function of $(q,\dot q,t)$ near every point where the fiber Hessian is invertible.
[/guided]
[/step]