There is a collection of animated Lorenz Attractors with interactive code created with 140 (or less!) characters of JavaScript:
https://www.dwitter.net/h/lorenz
Disclaimer: I created some of them :-)
glChAoSP /wglChAoSP is RealTime 3D strange attractors GPU explorer and hypercomplex fractals - over 200 Million particles.
It's freeware, Open Source, and both: native Multi Platform and WebGL... and display not only Lorenz Attractor, but over 100 object types between Attractors, Hypercomplex Fractals, and DLA/DLA3D (Diffusion Limited Aggregation)
Just one thing: on a touchscreen interface with Firefox, zooming out works well but zooming in is sometimes jittery or not registering at all. (Panning and rotating work fine.)
There is no advantage to using a better ODE solver on this example because the sensitivity to errors is so high. You can measure this using the uncertainty quantification methods in the Julia ODE solver suite:
where divergence on the Lorenz attactor tends to occur by t=80 or so even with accuracy of 1e-16.
But the funny thing about chaotic problems is that the shadowing theorem holds, which is:
>Although a numerically computed chaotic trajectory diverges exponentially from the true trajectory with the same initial coordinates, there exists an errorless trajectory with a slightly different initial condition that stays near ("shadows") the numerically computed one.
So you might as well just use Euler's method with high error, because it's backwards stable for this calculation, i.e. it gives a trajectory on the attractor, just the wrong trajectory. But since every method gives an O(1) error wrong trajectory after a short finite time, you might as well use the cheapest most error prone but convergent method.
While we are both here, if you are not aware of this paper it might give you a giggle (3500th order simulations on a supercomputer): https://arxiv.org/abs/1305.4222
They simulated Lorenz reliably up to 10000 time units; I managed 1400 units using MPFR on my NUC using my own code (500th order, took about 13 hours!).
Haha great! This is 100% a "might as well put it on Arxiv since no one reviewer would ever see this as significant" lol. It's at least very fun. I'm going to have to save this one.
That paper needs an update. Taylor methods utilize AD effectively to get something that seems efficient. But there's a lot of theory to explain why they aren't efficient in practice, mostly because of truncation error coefficient sizes. Also, since the time of the article, there have been much better high order explicit RK methods that have been derived, and improved schemes for adaptive order explicit extrapolation. When we benchmark these against Taylor methods which use an optimal form of AD (provided by TaylorIntegration.jl), we still see about a 10x performance difference when at the edge of Float64 range:
As you get outside of Float64 range, the parallelized explicit extrapolation methods seem to do better. This is something we're going to be getting a lot more precise benchmarks in JuliaDiffEq now though, since we have optimized generic implementations of each of the methods and can thus use Float128 and arbitrary precision on all of the methods simultaneously to really figure out the cutoff.
But thinking about efficiency and accuracy misses the point of Taylor methods. The Taylor methods have a way of getting rigorous error bounds on the integration (like interval arithmetic, but for the continuous ODE), which is the real reason why they can be useful, with the extra cost of course.
Erm, sure I know that name from somewhere! I won't dispute your claims about alternative methods, but I've found it hard to find any open documentation about this stuff, so I am unaware of it myself. Then again, I hereby announce that I am an _amateur_ at this, so I'm not taking criticism personally!
In fact one of the reasons I wrote my stuff (FWIW in c and Python) is because I couldn't follow the code in TaylorIntegration.jl ;) My stuff is _really_ simple, and I know _exactly_ how it works.
High order means much less dependence on step size. Automatic (exact) differentiation means no roundoff from finite differences with small step size.
They use a very clever and involved code generation approach, I have just hardcoded the recurrences in my implementation (much simpler). https://github.com/m4r35n357/ODE-Playground
It is not a black box approach like RK4, it depends on the model using a (hardly at all) restricted set of functions. Piecewise definitions are OK though.