  Following the theory presented in Part I, where we derived the
Euler-Lagrange, nonlinear, second-order kinematic ray tracing equation for
smooth heterogeneous general anisotropic media, this part is devoted for its
numerical finite-element solution. For a given initial-guess trajectory between
two fixed endpoints, discretized with a set of nodes, we update the location
and direction of the ray trajectories at the nodes, to obtain the nearest
stationary ray path. Starting with the Euler-Lagrange equation derived in Part
I, we apply the weak formulation and the Galerkin method to reduce this
second-order, ordinary differential equation into a nonlinear, local,
first-order, weighted residual algebraic equation set. The solution is based on
a finite element approach with the Hermite polynomial interpolation, for
computing the ray characteristics between the nodes. The Hermite interpolation
is a natural choice, since, in addition to the nodal locations, also the ray
directions at the nodes are independent degrees of freedom of the finite
element method. This is in particular important in anisotropic media. We then
construct the target function to be minimized, where we distinguish between two
functions triggered by the two possible types of stationary rays: minimum
traveltime and saddle-point solutions, mainly in cases where the stationary
paths include caustics. The target functions include two penalty terms related
to two essential constraints. The first is related to the distribution of the
nodes along the ray and the second enforces the normalization of the ray
direction to a unit length. The minimization process involves the computation
of the traveltime gradient vector and the Hessian matrix. The latter is used
also to analyze the type of the stationary ray. Finally, we demonstrate the
efficiency and accuracy of the proposed method along three canonical examples.
