A Taylor integrator compatible with SymPy

[Francesco Biscani]

Hello,

we just released the latest version of our Taylor integrator heyoka.py:

https://github.com/bluescarni/heyoka.py

heyoka.py is an implementation of Taylor's method for the numerical integration of systems of ODEs based on automatic differentiation and just-in-time compilation via LLVM.

Current features include:

- support for both double-precision and extended-precision floating-point types,
- the ability to maintain machine precision accuracy over tens of billions of timesteps,

- high-precision zero-cost dense output,
- accurate and reliable event detection,

- excellent performance,

- batch mode integration to harness the power of modern SIMD instruction sets.

heyoka.py needs to represent the ODEs symbolically in order to apply the automatic differentiation rules necessary for an efficient implementation of Taylor's method. For this purpose, heyoka.py uses its own expression system, but in recent versions we added the ability to convert heyoka.py's symbolic expressions to/from SymPy. Here's a simple example of interoperability between heyoka.py and SymPy:

https://bluescarni.github.io/heyoka.py/notebooks/sympy_interop.html

Here instead is a non-trivial example where the equations of motion are formulated via SymPy's classical mechanics module and then integrated via heyoka.py:

https://bluescarni.github.io/heyoka.py/notebooks/tides_spokes.html

This second example also shows how the common subexpression elimination capabilities of heyoka.py were able to drastically simplify highly-complex Lagrangian equations.

As a long-time observer/user of SymPy, I thought that other SymPy users might find this project interesting. I am also looking for feedback on our SymPy conversions facilities, as this is my first time digging into the SymPy expression system internals.

Thanks and kind regards,

  Francesco

[Chris Smith]

I had a brief look at the tides and spokes. My impression was that it was well presented.

I am curious on the large expressions for which you used subexpression elimination and wonder if you have tried `lambdify(..., cse=True)` on them with any favorable results.

I look forward to spending some more time looking at what you have done.

Best regards,

/c

[Aaron Meurer]

On Thu, Sep 30, 2021 at 11:59 AM Chris Smith <smi...@gmail.com> wrote:
>
> I had a brief look at the tides and spokes. My impression was that it was well presented.
>
> I am curious on the large expressions for which you used subexpression elimination and wonder if you have tried `lambdify(..., cse=True)` on them with any favorable results.

The cse argument to lambdify is new in SymPy 1.9, which is still in
pre-release (as an aside, it looks like the release note entry for
https://github.com/sympy/sympy/pull/21546 wasn't very accurate).

Aaron Meurer

>
> I look forward to spending some more time looking at what you have done.
>
> Best regards,
> /c
>
> On Wednesday, September 29, 2021 at 4:23:17 PM UTC-5 blues...@gmail.com wrote:
>>
>> Hello,
>>
>> we just released the latest version of our Taylor integrator heyoka.py:
>>
>> https://github.com/bluescarni/heyoka.py
>>
>> heyoka.py is an implementation of Taylor's method for the numerical integration of systems of ODEs based on automatic differentiation and just-in-time compilation via LLVM.
>>
>> Current features include:
>>
>> - support for both double-precision and extended-precision floating-point types,
>> - the ability to maintain machine precision accuracy over tens of billions of timesteps,
>> - high-precision zero-cost dense output,
>> - accurate and reliable event detection,
>> - excellent performance,
>> - batch mode integration to harness the power of modern SIMD instruction sets.
>>
>> heyoka.py needs to represent the ODEs symbolically in order to apply the automatic differentiation rules necessary for an efficient implementation of Taylor's method. For this purpose, heyoka.py uses its own expression system, but in recent versions we added the ability to convert heyoka.py's symbolic expressions to/from SymPy. Here's a simple example of interoperability between heyoka.py and SymPy:
>>
>> https://bluescarni.github.io/heyoka.py/notebooks/sympy_interop.html
>>
>> Here instead is a non-trivial example where the equations of motion are formulated via SymPy's classical mechanics module and then integrated via heyoka.py:
>>
>> https://bluescarni.github.io/heyoka.py/notebooks/tides_spokes.html
>>
>> This second example also shows how the common subexpression elimination capabilities of heyoka.py were able to drastically simplify highly-complex Lagrangian equations.
>>
>> As a long-time observer/user of SymPy, I thought that other SymPy users might find this project interesting. I am also looking for feedback on our SymPy conversions facilities, as this is my first time digging into the SymPy expression system internals.
>>
>> Thanks and kind regards,
>>
>> Francesco
>

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/d607d3af-4b4d-4e10-a291-e3790df8bc05n%40googlegroups.com.

[Francesco Biscani]

Hi Chris,

thanks for the kind words.

I haven't tried the cse=True approach in SymPy, I'll check it out when 1.9 comes out. I would imagine that the outcome will be similar to the CSE in our expression system, but it will be interesting to know if there's something new we can learn.

Kind regards,

  Francesco.

--

[Jason Moore]

Hi,

Thanks for sharing. I'll try it out with some other mechanics problems. Looks nice!

Jason

moorepants.info
+01 530-601-9791

--

You received this message because you are subscribed to the Google Groups "sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHExjCv%3DwxH6ZRpkkMqJDMmXg-h1Q43Q5dZq4ahYcGE8MGfFHA%40mail.gmail.com.

[Francesco Biscani]

Hi Jason,

cheers! The interoperability with SymPy and its mechanics module is a work in progress, so please let me know what you think we could improve on.

Kind regards,

  Francesco.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAP7f1Aj2cFtgQBS_5oqPoEq9DCyF%2BJfZRr-W3UaLi_WfBaPNag%40mail.gmail.com.

[Jason Moore]

Francesco,

I'm prepping a new course I teach at TU Delft in multibody dynamics and one week is about numerical integration. I'll try things out when I'm working on those lessons.

Jason

moorepants.info
+01 530-601-9791

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHExjCuG0Z_nP-7Q3Tba%2B_3KfB2K%2B1TJQZOZoSNB3mtxdw498A%40mail.gmail.com.

[Francesco Biscani]

Cheers, let me know if you have any issue/question.

Kind regards,

  Francesco.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAP7f1AjAukR1k1jGucmULgbzDjvcZFFuvqHxbHxCMHnp3Z7boA%40mail.gmail.com.