  Pulsar timing arrays (PTAs) detect low-frequency gravitational waves (GWs) by
looking for correlated deviations in pulse arrival times. Current Bayesian
searches use Markov Chain Monte Carlo (MCMC) methods, which struggle to sample
the large number of parameters needed to model the PTA and GW signals. As the
data span and number of pulsars increase, this problem will only worsen. An
alternative Monte Carlo sampling method, Hamiltonian Monte Carlo (HMC),
utilizes Hamiltonian dynamics to produce sample proposals informed by
first-order gradients of the model likelihood. This in turn allows it to
converge faster to high dimensional distributions. We implement HMC as an
alternative sampling method in our search for an isotropic stochastic GW
background, and show that this method produces equivalent statistical results
to similar analyses run with standard MCMC techniques, while requiring 100-200
times fewer samples. We show that the speed of HMC sample generation scales as
$\mathcal{O}(N_\mathrm{psr}^{5/4})$ where $N_\mathrm{psr}$ is the number of
pulsars, compared to $\mathcal{O}(N_\mathrm{psr}^2)$ for MCMC methods. These
factors offset the increased time required to generate a sample using HMC,
demonstrating the value of adopting HMC techniques for PTAs.
