In page 16 of McDuff & Salamon's Introduction to Symplectic Topology they prove that the critical points $z = (x,y)\colon [a,b] \to \Bbb R^{2n}$ of the action integral $$\Phi_H(z) = \int_a^b \langle y, \dot{x}\rangle - H(t,x,y)\,{\rm d}t$$satisfy Hamilton's equations $\dot{x} = \partial_yH$ and $\dot{y} = -\partial_xH$. Namely, they compute the first variation of $\Phi_H(z)$ as $$\widehat{\Phi_H}(z) = \int_a^b \langle \eta,\dot{x}-\partial_yH\rangle - \langle \xi, \dot{y}+\partial_xH\rangle\,{\rm d}t,$$where $(\xi,\eta)$ is the variational vector field of a variation of $z$ with fixed endpoints. I have no problems whatsoever with the proof they present, but I am having trouble formulating a generalization of this statement in the setting of differentiable manifolds. What I have in mind is:
Let $Q$ be a differentiable manifold and $H\colon T^*Q \to \Bbb R$ a smooth Hamiltonian. Consider the action integral $$\mathscr{A}^H(x,{\sf p}) = \int_a^b \mathbb{F}H(x(t),{\sf p}(t)){\sf p}(t) - H(x(t),{\sf p}(t))\,{\rm d}t,$$where $(x,{\sf p})\colon [a,b] \to T^*Q$ is a smooth curve and $\mathbb{F}H$ is the fiber derivative of $H$, given in coordinates by $$\mathbb{F}H(x,{\sf p}) = \sum_{k=1}^n \frac{\partial H}{\partial p_k}(x,{\sf p})\frac{\partial}{\partial q^k}\bigg|_x.$$Then if $(x,{\sf p})$ is a critical point of $\mathscr{A}^H$ and we write $$(x(t),{\sf p}(t)) = (q^1(t),\ldots, q^n(t),p_1(t),\ldots, p_n(t)),$$we have Hamilton's equations $$\frac{{\rm d}q^k}{{\rm d}t}(t) = \frac{\partial H}{\partial p_k}(x(t),{\sf p}(t)) \qquad\mbox{and}\qquad \frac{{\rm d}p_k}{{\rm d}t}(t) = -\frac{\partial H}{\partial q^k}(x(t),{\sf p}(t)).$$
However, I'm not entirely sure of this, as a priori there's no relation between $x(t)$ and ${\sf p}(t)$ (except for ${\sf p}(t) \in T_{x(t)}^*Q$), in contrast with curves $(x(t),\dot{x}(t))$ in $TQ$, when analyzing the Lagrangian case. In particular, I cannot reproduce the proof given in McDuff & Salamon, since there's no apparent way to use integration by parts, as we have no $t$-derivatives in the integrand.
I don't really expect anyone to fix the statement and provide the computation (although it would be nice), but I need help figuring out what sort of relation I am missing here (with that I should be able to run the computation of the first variation myself). In particular, I am not assuming that $H$ is hyperregular, and we do not have Legendre transformations at our disposal. What to do?
Perhaps it is best to compute things intrinsically.
The action functional in the cotangent bundle is given by $$A_H(x)=\int_{[0,1]} x^* \alpha-\int_{[0,1]} H \circ x,$$ where $\alpha$ is the canonical $1$-form in the cotangent bundle. You can check that this coincides with your formula in the case of $\mathbb{R}^{2n}$.
Differentiating and using Cartan's magic formula yields that the stationary points of the action functional satisfy Hamilton's equations. (Recall that $d\alpha =\omega$.)