/* file LE.c 

 It computes the Lyapunov exponents 
 of a dynamical system 
  dx/dt = f(x)
 by simultaneous numerical integration of original and variational equations.

 This demo in written for the 3-D Lorenz system and
 it was intended to be easily modifiable.

*/

# include "easynum.c"    /* order is significant */
# include "model.c"

/*
You can retrieve the files easynum.c and model.c by ftp from
ftp://tiger.vscht.cz/pub/xpplot/easynum.tar.Z

For plotting data you can try xpplot.tar.Z
*/

main ()
{
 real 
  s  = 10,     /* parameters */
  r  = 28,
  b  = 8.0/3.0,

  TT = 100,    /* transient time */
  T  = 100;    /* time to integrate */
 vector
  x,           /* vector of the state variables */
  le,          /* vector of resulting Lyapunov exponents */
  pars;        /* vector of parameters */
 int 
  output_wanted = 0, /* try 1 to see the results of integration  */
  steps = 10,  /* number of integration steps between re-orthonormalizations */
  nle = 3,     /* number of L.E. wanted */
  d = 3;       /* dimension of the original system 
                  dimension of the full system will be dim = d * (nle+1)   */

/* store parameters to the vector to pass them for the integrator */
pars[0] = s;
pars[1] = r;
pars[2] = b;

/* set initial conditions */
 x[0] = 0.12;
 x[1] = 0.34;
 x[2] = 0.56;

/* compute the Lyapunov exponents */
LE (x,le,rhs,pars,T,TT,output_wanted,steps,nle,d);

/* print them */
writevector ("LEs are: ",le,nle);
}
/*********************************************/