In #473, we plumbed the PhaseSpacePosition arrays of shape (ndim, n) all the way down to the C level, so they can be used in vectorization-friendly gradient functions without any copies or transposes. However, this now means that the gradient path uses a different data layout than the energy/density/hessian paths, which is a potential maintenance burden. Ideally, we would re-plumb those paths in the same way.