How to find a particle's dynamics in general relativity?

As in classical mechanics, there are several approaches to determining the dynamics of a particle in general relativity. Crucially we should observe that, in general relativity, the effect of a gravitational field is entirely captured by the behaviour of the metric. Providing the particle in question experiences no other external influences, we can determine its motion using the geodesic equation as you suppose. This is one of the central axioms of general relativity:

A particle moving only under the influence of gravity follows geodesics in curved spacetime. These geodesics are timelike if the particle is massive and null if the particle is massless.

Dynamics à la Newton:

Knowing the above, we can directly write down the equation governing particle motion. In Riemannian geometry, it is natural to define a geodesic as a path that extremises the length between two fixed points. For us, this definition comes with some subtleties, due to the fact that two distinct points can be zero distance apart. For our purposes, there is a better definition:

A geodesic is a path whose tangent vector is parallel transported along it. In wordier terms, if you take the tangent vector at one point, and then carry it along the path without changing the direction it points in, the vector remains tangent to the path at the destination.

In less wordy terms, geodesics are straight lines.

Mathematically this is expressed as

$$ \nabla_U U = 0 \qquad \mathrm{where} \qquad U^\mu = \frac{\mathrm{d}x^\mu}{\mathrm{d}\lambda}\,.$$ Here $\lambda$ is a parameter along the curve and $U$ is tangent to it. If we expand out this equation using the definition of the covariant derivative, we arrive at the famous geodesic equation: $$ \frac{\mathrm{d}^2 x^\mu}{\mathrm{d}\lambda^2} + \Gamma^\mu{}_{\nu \rho}\frac{\mathrm{d}x^\nu}{\mathrm{d}\lambda}\frac{\mathrm{d}x^\rho}{\mathrm{d}\lambda} = 0 \,.$$ This equation is the analogue of Newton's first law: particles under the influence of no external force move in straight lines at constant speed. It is nothing other than

$$\vec{a} = 0\,.$$

This comparison informs us how to modify the geodesic equation when external forces are present. For definiteness, suppose that the particle is moving under the influence of an electromagnetic field, in addition to gravity. We simply take geodesic equation and replace the right hand side with force over mass. Writing the Lorentz force in the covariant form used in special relativity, and choosing $\lambda$ to be proper time $\tau$, we find that a charged particle obeys the equation: $$ \frac{\mathrm{d}^2 x^\mu}{\mathrm{d}\tau^2} + \Gamma^\mu{}_{\nu \rho}\frac{\mathrm{d}x^\nu}{\mathrm{d}\tau}\frac{\mathrm{d}x^\rho}{\mathrm{d}\tau} = \frac{q}{m} F^\mu{}_\nu \frac{\mathrm{d}x^\nu}{\mathrm{d}\tau} \,.$$ Note that if we turn off the gravitational field by setting the connection components to zero, we recover the special relativistic form of Newton's second law for a particle moving in an electromagnetic field, as we must do.

Dynamics à la Lagrange:

One of the advantages of Lagrangian over Newtonian mechanics is that it's often easier to write down some scalar 'energy' that characterises our system than it is to specify the vectorial forces acting on it. With no external influences present, we want an action that, when extremised, gives us geodesics. For massive particles, we can employ the definition of geodesic familiar from Riemannian geometry – action equals length:

$$S = m\int\mathrm{d}\tau = m\int \mathrm{d}\lambda \sqrt{-g_{\mu \nu}\frac{\mathrm{d}x^\mu}{\mathrm{d}\lambda}\frac{\mathrm{d}x^\nu}{\mathrm{d}\lambda}} \,.$$

It's a worthwhile exercise to check that this indeed yields the geodesic equation when extremised (note that one should be careful about the limits on the integral when performing the variation!).

With this free-particle action in place, it's a simple matter to add extra terms that capture the effect of interactions of the particle with external fields or other particles. If an electromagnetic field is present, the action becomes

$$ S = m\int \mathrm{d} \tau - q \int \mathrm{d} \tau A_\mu \frac{\mathrm{d}x^\mu}{\mathrm{d}\tau} \,. $$

For massless particles this picture is complicated by the fact that every null path has zero length. To talk about the Lagrangian mechanics of massless particles, we need to introduce a new action which contains a Lagrange multiplier $e(\lambda)$ defined along the worldline. The following action is valid for particles of any mass:

$$ S = \frac{1}{2} \int \mathrm{d} \lambda \left(\frac{1}{e(\lambda)} g_{\mu \nu}\frac{\mathrm{d}x^\mu}{\mathrm{d}\lambda}\frac{\mathrm{d}x^\nu}{\mathrm{d}\lambda} - m^2 e(\lambda)\right)\,.$$

This is a technical point – more information can be found in this and this stackexchange post. At the end of the day, we're defining these actions so that they return the geodesic equation when extremised. Adding extra terms to the action to account for interactions is then often easier than adding extra terms directly to the geodesic equation.

Dynamics à la Hamilton:

Perhaps some other time; this answer is already long enough. Let me just note here that it is entirely possible to formulate a Hamiltonian approach to particle dynamics in a curved spacetime, but a little tricky. There are additional subtleties associated with the fact that the conjugate momenta are not independent – they satisfy a mass-shell constraint – which compels us to include a Lagrange multiplier as above. For more information see these lecture notes by Paul Townsend.

In general relativity, in the absence of other fields, point particles move along geodesics. Note that a Newtonian potential is replaced in general relativity by a component of the metric. In the presence of other fields, the motion need no longer be geodetic.

More specifically, for a particle of constant mass $m>0$ and charge $e$, there is a simple covariant Hamiltonian description: In the presence of a gravitational field (metric tensor) $G(x)$ and an electromagnetic field $A(x)$, the Hamiltonian is given by $$H=\frac{1}{2m}(p+eA(q))^TG(q)(p+eA(q)).$$ This Hamiltonian generates the trajectory parameterized by the particle's eigentime. For an uncharged particle (or in the absence of an electromagnetic field) this reduced to $H=p^TG(q)p/2m$, which gives the geodetic motion.

Details are discussed, e.g., in Thirring's course in mathematical Physics, Part I, Chapter 5.