From: spellucci@mathematik.tu-darmstadt.de (Peter Spellucci)
Subject: Re: looking for efficient method to solve an nonlinear ODE system
Date: 20 Jan 2000 10:39:21 GMT
Newsgroups: sci.math.num-analysis
Summary: [missing]

In article <864t37$sgs$1@bird.wu-wien.ac.at>,
 hh@hal. (Herbert Haas) writes:
|> 
|> Hello
|> 
|> I have been asking this question recently but I think there have been 
|> something gone wrong with my posting...
|> 
|> I need a more efficient method than Runke Kutta 4th order 
|> (which I use currently)
|> to solve an autonomous system like
|> 
|> x1(t)' = f1 (x1(t), x2(t), x3(t), ...)
|> x2(t)' = f2 (x1(t), x2(t), x3(t), ...)
|> x3(t)' = f3 (x1(t), x2(t), x3(t), ...) 
|> 
|> where I know the initial conditions at t=0: x1(0), x2(0), x3(0),....
|> and from which I want to calculate future data
|> 
the answer depends ........... very very much on the properties of
your system. If you observe "inefficiency", what does this really mean?
rk4 obtains order 4 with 4 values of F=(f1,f2,...) and very little 
overhead. but rk4 as such (i am afraid you indeed mean that first rk4 which 
is "the" rk in old books) has no built in mechanism for error control, 
using halving (or doubling) the spepsize increases effort by 50%.
"better" are dormand prince 4/5 (with 6 f) or dormand/prince 7/8 

for codes see the book of Hairer, Norsett and Wanner,
solving ordinary differential equaitions I,II (springer)
see the URL below 

which allow easily for a stepsize control via the local error. 
but, e.g. dormand/prince 4/5 (n.b. there exists an even better pair now
by shampine) need 6 F's per step, so again this effort is "high".
You can obtain methods of arbitrarily high order with only 1 F per step:
these are the Adams-Bashforth formulae. But these are never used in practice,
since their domain of absolute stability is so small.

here we come to a second point, which in practice is of even more importance
than pure local precision:
practically important ODE systems are mostly stable or even dissipative.
And, of course, one wishes that the discretized scheme is so too. 
roughly that means for rk4 that  
    time step size * Norm(Jacobian F) <= 2.8
so, if you have what is named a "stiff system", your time step size must be very
small. But there are methods which avoid such restrictive conditions.
these are all implicit and require the solution of equations to obtain the
next x-value. if the ODE is nonlinear, then these equations are nonlinear too
and a new problem arises: solver (ususally simplified Newton) plus 
initial guess (usually x_previous) plus stepsize restriction to make
Newton convergent. (But in case of stiffness, this is a secondary problem).
If F is weakly nonlinear, say
   F(x) = Ax + g(x)   with g(x) smooth, norm(Jacobian(g)) small,
then Rosenbrock-Wanner codes are very good, but otherwise, for a strongly
nonlinear stiff
and dissipative system, implicit Runge Kutta based on RadauII integration
is the choice . For codes see
 http://www.unige.ch/math/folks/hairer/

Even harder becomes the job if the solution X contains highly oscillatory
components without damping. then little choice remains ...
implicit midpoint rule, trapezoidal rule ,    ...

special problems need special treatment ...
hamiltonian systems -> symplectic integrators
periodic solutions -> P-stable methods   

for software see the URL above and of course
http://www.netlib.org -> ode, toms, 

hope that helps
peter