  The anisotropic diffusion equation is imperative in understanding cosmic ray
diffusion across the Galaxy, the heliosphere, and its interplay with the
ambient magnetic field. This diffusion term contributes to the highly stiff
nature of the CR transport equation. In order to conduct numerical simulations
of time-dependent cosmic ray transport, implicit integrators have been
traditionally favoured over the CFL-bound explicit integrators in order to be
able to take large step sizes. We propose exponential methods that directly
compute the exponential of the matrix to solve the linear anisotropic diffusion
equation. These methods allow us to take even larger step sizes; in certain
cases, we are able to choose a step size as large as the simulation time, i.e.,
only one time step. This can substantially speed-up the simulations whilst
generating highly accurate solutions (l2 error $\leq 10^{-10}$). Additionally,
we test an approach based on extracting a constant diffusion coefficient from
the anisotropic diffusion equation, where the constant coefficient term is
solved implicitly or exponentially and the remainder is treated using some
explicit method. We find that this approach, for homogeneous linear problems,
is unable to improve on the exponential-based methods that directly evaluate
the matrix exponential.
