heyoka

heyoka#

The heyókȟa […] is a kind of sacred clown in the culture of the Sioux (Lakota and Dakota people) of the Great Plains of North America. The heyoka is a contrarian, jester, and satirist, who speaks, moves and reacts in an opposite fashion to the people around them.

heyoka is a C++ library for the integration of ordinary differential equations (ODEs) via Taylor’s method, based on automatic differentiation techniques and aggressive just-in-time compilation via LLVM. Notable features include:

  • support for single-precision, double-precision, extended-precision (80-bit and 128-bit), and arbitrary-precision floating-point types,

  • high-precision zero-cost dense output,

  • accurate and reliable event detection,

  • builtin support for analytical mechanics - bring your own Lagrangians/Hamiltonians and let heyoka formulate and solve the equations of motion,

  • builtin support for high-order variational equations - compute not only the solution, but also its partial derivatives,

  • builtin support for machine learning applications via neural network models,

  • the ability to maintain machine precision accuracy over tens of billions of timesteps,

  • batch mode integration to harness the power of modern SIMD instruction sets (including AVX/AVX2/AVX-512/Neon/VSX),

  • ensemble simulations and automatic parallelisation.

If you prefer using Python rather than C++, heyoka can be used from Python via heyoka.py, its Python bindings.

If you are using heyoka as part of your research, teaching, or other activities, we would be grateful if you could star the repository and/or cite our work. For citation purposes, you can use the following BibTex entry, which refers to the heyoka paper (arXiv preprint):

@article{10.1093/mnras/stab1032,
    author = {Biscani, Francesco and Izzo, Dario},
    title = "{Revisiting high-order Taylor methods for astrodynamics and celestial mechanics}",
    journal = {Monthly Notices of the Royal Astronomical Society},
    volume = {504},
    number = {2},
    pages = {2614-2628},
    year = {2021},
    month = {04},
    issn = {0035-8711},
    doi = {10.1093/mnras/stab1032},
    url = {https://doi.org/10.1093/mnras/stab1032},
    eprint = {https://academic.oup.com/mnras/article-pdf/504/2/2614/37750349/stab1032.pdf}
}

heyoka’s novel event detection system is described in the following paper (arXiv preprint):

@article{10.1093/mnras/stac1092,
    author = {Biscani, Francesco and Izzo, Dario},
    title = "{Reliable event detection for Taylor methods in astrodynamics}",
    journal = {Monthly Notices of the Royal Astronomical Society},
    volume = {513},
    number = {4},
    pages = {4833-4844},
    year = {2022},
    month = {04},
    issn = {0035-8711},
    doi = {10.1093/mnras/stac1092},
    url = {https://doi.org/10.1093/mnras/stac1092},
    eprint = {https://academic.oup.com/mnras/article-pdf/513/4/4833/43796551/stac1092.pdf}
}

As a simple example, consider the ODE system corresponding to the pendulum,

\[\begin{split}\begin{cases} x^\prime = v \\ v^\prime = -9.8 \sin x \end{cases}\end{split}\]

with initial conditions

\[\begin{split}\begin{cases} x\left( 0 \right) = 0.05 \\ v\left( 0 \right) = 0.025 \end{cases}\end{split}\]

Here’s how the ODE system is defined and numerically integrated in heyoka:

#include <iostream>

#include <heyoka/heyoka.hpp>

using namespace heyoka;

int main()
{
    // Create the symbolic variables x and v.
    auto [x, v] = make_vars("x", "v");

    // Create the integrator object
    // in double precision.
    auto ta = taylor_adaptive<double>{// Definition of the ODE system:
                                      // x' = v
                                      // v' = -9.8 * sin(x)
                                      {prime(x) = v, prime(v) = -9.8 * sin(x)},
                                      // Initial conditions
                                      // for x and v.
                                      {0.05, 0.025}};

    // Integrate for 10 time units.
    ta.propagate_for(10.);

    // Print the state vector.
    std::cout << "x(10) = " << ta.get_state()[0] << '\n';
    std::cout << "v(10) = " << ta.get_state()[1] << '\n';
}

Output:

x(10) = 0.0487397
y(10) = 0.0429423

heyoka is released under the MPL-2.0 license. The authors are Francesco Biscani and Dario Izzo (European Space Agency).