Just another site

An improved PARI/gp script for solving the Lorenz ODE system

leave a comment »

Today, I’m posting a variant of the PARI/gp script
to solve the Lorenz system of differential equations.

The use of factorials and binomial coefficients has been
eliminated, so I’d expect it to run somewhat faster than
the script from an earlier post.


The PARI/gp script:


Lorenz(X0, Y0, Z0) =
order = 300;
dt = 1.0/250 ;
for(J=1,order, fa[J] = factorial(J));
x0 = X0; y0 = Y0; z0 = Z0;
yc = y0;
zc = z0;
sig = 10.0;
rho = 28.0;
beta = 8.0/3 ;
for(Q = 1,1000,
for(L = 1, 250,
yp[1] = yc;
zp[1] = zc;
for(J=2, order+1,
xp[J] = (sig*(yp[J-1]-xp[J-1]))/(J-1) ;
yp[J] = (rho*xp[J-1] – sum(K=0,J-2, xp[K+1]*zp[J-1-K]) -yp[J-1])/(J-1);
zp[J] = (sum(K=0,J-2, xp[K+1]*yp[J-1-K]) – beta*zp[J-1])/(J-1); );
xn = xc + sum(J=1, order, xp[J+1]*dt^J);
yn = yc + sum(J=1, order, yp[J+1]*dt^J);
zn = zc + sum(J=1, order, zp[J+1]*dt^J) ;
xc = xn;
yc = yn;
zc = zn; );
print(Q, ” “, precision(xc, 20)))


Results for x(100), x(200), x(300) and x(400), which agree with
the Kehlett and Logg reference solution for those values:


? Lorenz(1.0, 0.0, 0.0)

100 12.2860066770189848088523959400428665893
200 -0.428878362621396885476845988868367367420
300 -7.9465362714960512478195581280017205595
400 5.9432455860803396010549930121329814024


Advancing the numerical solution by one time unit,
for example from t = 100.0 to t = 101.0, takes approximately
50 seconds. In the inner loop “for L = 1,250”,
one can see that there is no factorial function and
there are no binomial coefficients. That’s why I consider
this script an improvement over the earlier one.

I believe there is a trade-off that can be made, for the
same precision, between a smaller time-step and a smaller
value of the order of the Taylor series.

So, to minimize the total time, one would look at
pairs of time-step and Taylor series order giving the
same precision , e.g. at t = 1000, or even t = 250,
to find the pair which runs in the least amount of time.

Personal note:
the file is in the directory :
on my PC.


Written by meditationatae

August 11, 2017 at 3:53 pm

Posted in History

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: