All the problem, in my view, arises form the fact that the EL equations are introduced in a too sloppy way. (The variational approach makes even more obscure an obscure setup.)
Actually,
The coordinates $\dot{q}^k$ and $q^k$ are independent and they become dependent only "on-shell", i.e., on the curve that solvese the equations of motion. Only then $\dot{q}$ becomes the time derivative of $q$ and this should be expliticly stated!
The EL equations are 2n not n: in coordinatess one looks for a curve $$\gamma : I \ni t \mapsto (t, q(t),\dot{q}(t)) \in \mathbb{R}\times \mathbb{R}^n\times \mathbb{R}^n$$
such that (notice the order of the various operations, the partial derivatives are computed before evaluating the Lagrangian on the curve)
$$\left \{ \begin{array}{rl}
\left.\frac{d}{dt}\right|_{\gamma(t)} \frac{\partial L(t,q, \dot{q})}{\partial \dot{q}^k} -\left.\frac{\partial L(t,q, \dot{q})}{\partial q^k}\right|_{\gamma(t)}&=0\\
\left.\frac{d}{dt}\right|_{\gamma(t)} q^k &= \left.\dot{q}^k\right|_{\gamma(t)}
\end{array}
\right.$$
I stress that these equations are 1st order ordinary differential equations.
In particular the Lagrangian is a function of $2n+1$ independent variables:
$$t, q^1,\ldots, q^n, \dot{q}^1, \ldots, \dot{q}^n$$
A natural way to globally fix all this approach is the use of a jet bundle viewing the spacetime of kinetic states, locally determined by natural coordinates $(t, q^1,\ldots, q^n, \dot{q}^1, \ldots, \dot{q}^n)$ a as a fiber bundle over the real line.
In natural charts over the spacetime of kinetic states the transformation maps are
$$\overline{t} = t+c $$
$$\overline{q}^k = \overline{q}^k(t,q)$$
$$ \dot{\overline{q}}^k =\frac{\partial \overline{q}^k}{\partial t}+ \sum_{j=1}^n \frac{\partial \overline{q}^k }{\partial q^j} \dot{q}^j$$
These charts form a preferred atlas on the spacetime of kinetic states and define it. Notice that, in the transformation laws above, there is no given curve to use just because $\dot{q}^k$ is not the derivative, along a curve, of $q^k$. It is an independent variable, the said relation holds on the solutions of the motion equations however.
The EL equations written as above, when assuming that the Lagrangian is a scalar field on the spacetime of kinetic states, are invariant under the transformation laws written above.
A less elaborate way is to use the tengent bundle of the configuration space, but this approach is limitative, because also the configuration space does not exist individually, but is a canonical fiber of a fiber bundle: the configuration spacetime.
This approach, as it removes all notational ambiguities, makes very clear the passage to the Hamiltonian formalism.
The Hamiltonian formalism assumes tha there is a space, the spacetime of phases whose local coordinates are, in fact, $(t,q,p)$. These natural charts form a preferred atlas on the spacetime of phases and define it.
These local coordinates are associated to the Lagrangian local coordinates
$(t,q,\dot{q})$ by means of the Legendre diffeomorphism when a Lagrangian $L(t,q,\dot{q}))$ is given
$$Leg: (t,q,\dot{q}) \to (t(t,q,\dot{q}), q((t,q,\dot{q})), p((t,q,\dot{q}))$$
defined as
$$t=t$$
$$ q^k= q^k$$
$$p_k = \frac{\partial L(t,q,\dot{q})}{\partial \dot{q}^k}\:.$$
The basic fact is now that $$\gamma : I \ni t \mapsto (t, q(t),\dot{q}(t)) \in \mathbb{R}\times \mathbb{R}^n\times \mathbb{R}^n$$ satisfies the EL equations written above if and only if
$$\hat{\gamma}(t) := Leg(\gamma(t))$$
satisfies the Hamiltonian equations $$\left.\frac{d}{dt}\right|_{\hat{\gamma}(t)}q^k = \left.\frac{\partial H(t,q,p)}{\partial p_k}\right|_{\hat{\gamma}(t)}$$
$$\left.\frac{d}{dt}\right|_{\hat{\gamma}(t)} p_k= -\left.\frac{\partial H(t,q,p)}{\partial q^k}\right|_{\hat{\gamma}(t)}$$
Notice that I accurately avoided to indicate the temporal derivatives with the dot, since the dot denotes a Lagrangian coordinate! For instance $\dot{p}$ simply does not exist in this formalism: in particular it is not a Hamiltonian coordinate!
The Hamiltonian function $H(t,q,p)$ is constructed out of $L(t,q, \dot{q})$ by inverting $Leg$ in particular so that $\dot{q}^k=\dot{q}^k(t,q,p)$ is a known function:
$$H(t,q,p) := \sum_{k=1}^np_k \dot{q}^k(t,q,p)- L(t,q, \dot{q}(t,q,p))$$
The spacetime of phases has a nature that is independent of any Lagrangian one uses. Natural local coordinate systems are connected by the transformation laws
$$\overline{t} = t+c $$
$$\overline{q}^k = \overline{q}^k(t,q)$$
$$ \overline{p}_k =\sum_{j=1}^n \frac{\partial q^j }{\partial \overline{q}^k} p_j\:.$$
These relations guarantee that the Legendre transformation is actually global and not only local.
