numode.spad line 1 [edit on github]
This package is a suite of functions for the numerical integration of an ordinary differential equation of n variables: dy/dx = f(y, x)y is an n-vector All the routines are based on a 4-th order Runge-Kutta kernel. These routines generally have as arguments: n, the number of dependent variables; x1, the initial point; h, the step size; y, a vector of initial conditions of length n which upon exit contains the solution at x1 + h; derivs, a function which computes the right hand side of the ordinary differential equation: derivs(dydx, y, x) computes dydx, a vector which contains the derivative information. In order of increasing complexity: rk4(y, n, x1, h, derivs) advances the solution vector to x1 + h and return the values in y. rk4(y, n, x1, h, derivs, t1, t2, t3, t4) is the same as rk4(y, n, x1, h, derivs) except that you must provide 4 scratch arrays t1-t4 of size n. Starting with y at x1, rk4f(y, n, x1, x2, ns, derivs) uses ns fixed steps of a 4-th order Runge-Kutta integrator to advance the solution vector to x2 and return the values in y. Argument x2, is the final point, and ns, the number of steps to take. rk4qc(y, n, x1, step, eps, yscal, derivs) takes a 5-th order Runge-Kutta step with monitoring of local truncation to ensure accuracy and adjust stepsize. The function takes two half steps and one full step and scales the difference in solutions at the final point. If the error is within eps, the step is taken and the result is returned. If the error is not within eps, the stepsize if decreased and the procedure is tried again until the desired accuracy is reached. Upon input, an trial step size must be given and upon return, an estimate of the next step size to use is returned as well as the step size which produced the desired accuracy. The scaled error is computed as error = MAX(ABS((y2steps(i) - y1step(i))/yscal(i))) and this is compared against eps. If this is greater than eps, the step size is reduced accordingly to hnew = 0.9 * hdid * (error/eps)^(-1/4) If the error criterion is satisfied, then we check if the step size was too fine and return a more efficient one. If error > eps * (6.0E-04) then the next step size should be hnext = 0.9 * hdid * (error/\spadeps)^(-1/5) Otherwise hnext = 4.0 * hdid is returned. A more detailed discussion of this and related topics can be found in the book "Numerical Recipes" by W.Press, B.P. Flannery, S.A. Teukolsky, W.T. Vetterling published by Cambridge University Press. Argument step is a record of 3 floating point numbers (to_try , did , next), eps is the required accuracy, yscal is the scaling vector for the difference in solutions. On input, step.to_try should be the guess at a step size to achieve the accuracy. On output, step.did contains the step size which achieved the accuracy and step.next is the next step size to use. rk4qc(y, n, x1, step, eps, yscal, derivs, t1, t2, t3, t4, t5, t6, t7) is the same as rk4qc(y, n, x1, step, eps, yscal, derivs) except that the user must provide the 7 scratch arrays t1-t7 of size n. rk4a(y, n, x1, x2, eps, h, ns, derivs) is a driver program which uses rk4qc to integrate n ordinary differential equations starting at x1 to x2, keeping the local truncation error to within eps by changing the local step size. The scaling vector is defined as yscal(i) = abs(y(i)) + abs(h*dydx(i)) + tiny where y(i) is the solution at location x, dydx is the ordinary differential equation's right hand side, h is the current step size and tiny is 10 times the smallest positive number representable. The user must supply an estimate for a trial step size and the maximum number of calls to rk4qc to use. Argument x2 is the final point, eps is local truncation, ns is the maximum number of call to rk4qc to use.
rk4(y, n, x1, h, derivs) uses a 4-th order Runge-Kutta method to numerically integrate the ordinary differential equation dy/dx = f(y, x) of n variables, where y is an n-vector. Argument y is a vector of initial conditions of length n which upon exit contains the solution at x1 + h, n is the number of dependent variables, x1 is the initial point, h is the step size, and derivs is a function which computes the right hand side of the ordinary differential equation. For details, see NumericalOrdinaryDifferentialEquations.
rk4(y, n, x1, h, derivs, t1, t2, t3, t4) is the same as rk4(y, n, x1, h, derivs) except that you must provide 4 scratch arrays t1-t4 of size n. For details, see NumericalOrdinaryDifferentialEquations.
rk4a(y, n, x1, x2, eps, h, ns, derivs) is a driver function for the numerical integration of an ordinary differential equation dy/dx = f(y, x) of n variables, where y is an n-vector using a 4-th order Runge-Kutta method. For details, see NumericalOrdinaryDifferentialEquations.
rk4f(y, n, x1, x2, ns, derivs) uses a 4-th order Runge-Kutta method to numerically integrate the ordinary differential equation dy/dx = f(y, x) of n variables, where y is an n-vector. Starting with y at x1, this function uses ns fixed steps of a 4-th order Runge-Kutta integrator to advance the solution vector to x2 and return the values in y. For details, see NumericalOrdinaryDifferentialEquations.
rk4qc(y, n, x1, step, eps, yscal, derivs) is a subfunction for the numerical integration of an ordinary differential equation dy/dx = f(y, x) of n variables, where y is an n-vector using a 4-th order Runge-Kutta method. This function takes a 5-th order Runge-Kutta step with monitoring of local truncation to ensure accuracy and adjust stepsize. For details, see NumericalOrdinaryDifferentialEquations.
rk4qc(y, n, x1, step, eps, yscal, derivs, t1, t2, t3, t4, t5, t6, t7) is a subfunction for the numerical integration of an ordinary differential equation dy/dx = f(y, x) of n variables, where y is an n-vector using a 4-th order Runge-Kutta method. This function takes a 5-th order Runge-Kutta step with monitoring of local truncation to ensure accuracy and adjust stepsize. For details, see NumericalOrdinaryDifferentialEquations.