|
|
◆ doprin()
| void doprin |
( |
int |
n, |
|
|
void(*)(int, double, double *, double *) |
f, |
|
|
double * |
t, |
|
|
double |
y[], |
|
|
double |
yp[], |
|
|
double |
tout, |
|
|
double * |
rtol, |
|
|
double * |
atol, |
|
|
int |
itol, |
|
|
void(*)(int, double, double *, double *, int, int *) |
solout, |
|
|
int |
iout, |
|
|
double |
work[], |
|
|
int |
lwork, |
|
|
int |
iwork[], |
|
|
int |
liwork, |
|
|
int * |
info |
|
) |
| |
Initial value problem of ordinary differential equations (7(6)-th order Runge-Kutta-Nystrom method) (for second order differential equations)
NOTE - THIS ROUTINE IS DEPRECATED AND WILL BE REMOVED IN THE NEXT VERSION.
- Purpose
- This routine integrates a system of second order ordinary differential equations of the form
d2y/dt2 = f(t, y), y = y0, y' = y'0 at t = t0
where t0, y0 and y'0 are the given initial values of t, y and y', respectively. y and y' may be a vector if the above is a system of differential equations.
This routine is the code of the Runge-Kutta-Nystrom method of order 7(6) by Dormand and Prince with step size control.
See for details in the reference below.
- Parameters
-
| [in] | n | Number of differential equations. (n >= 1) |
| [in] | f | The user supplied subroutine, which calculates the derivatives of the differential equations, defined as follows. void f(int n, double t, const double y[], double ypp[])
{
ypp[i] = computed second derivative at t and y[] (for i = 0 to n-1)
}
where n is the number of equations, and ypp[] is the second derivatives at given t and y[], i.e. ypp[i] = d2yi/dt2 = fi(t, y[0], ..., y[n-1]) (i = 0 to n-1). |
| [in,out] | t | This routine integrates from t to tout. The initial point of the integration is to be given, and the last point of the final step will be returned.
[in] Initial value of the independent variable t.
[out] Last value of the independent variable t of the final step (normally equals to tout). The solution was successfully advanced to this point. It is possible to continue the integration to new point by recalling this routine with the new tout value. |
| [in,out] | y[] | Array y[ly] (ly >= n)
[in] Initial values of the dependent variables y[] at initial t.
[out] Computed solution approximation at last t (normally equals to tout). |
| [in,out] | yp[] | Array yp[lyp] (lyp >= n)
[in] Initial derivative values y'1, ..., y'n at initial t.
[out] Computed approximations of derivatives at last t (normally euqals to tout). |
| [in] | tout | Set tout to the point at which a solution is desired. Integration either forward in t (tout > t) or backward in t (tout < t) is permitted.
The routine advances the solution from t to tout using step sizes which are automatically selected so as to achieve the desired accuracy. |
| [in] | rtol | Scalar if itol = 0, or array rtol[lrtol] if itol = 1 (lrtol >= n) (rtol or all components of rtol[] >= 0)
The relative error tolerance(s). The routine keeps, roughly, the local error of y[i] below
rtol*abs(y[i]) + atol (i = 0 to n-1) (if itol = 0)
or
rtol[i]*abs(y[i]) + atol[i] (i = 0 to n-1) (if itol = 1).
Both rtol and atol, or, rtol[i] and atol[i] (i = 0 to n-1) should not be 0 at the same time. |
| [in] | atol | Scalar if itol = 0, or array atol[latol] if itol = 1 (latol >= n) (atol or all components of atol[] >= 0)
The absolute error tolerance(s). The routine keeps, roughly, the local error of y[i] below
rtol*abs(y[i]) + atol (i = 0 to n-1) (if itol = 0)
or
rtol[i]*abs(y[i]) + atol[i] (i = 0 to n-1) (if itol = 1).
Both rtol and atol, or, rtol[i] and atol[i] (i = 0 to n-1) should not be 0 at the same time. |
| [in] | itol | Specifies whether the parameters rtol and atol are scalars or arrays.
= 0: rtol and atol are scalars.
= 1: rtol and atol are arrays.
(For other values, itol = 0 is assumed.) |
| [in] | solout | The user supplied subroutine to print out the intermediate solutions, which is called after every successful step, defined as follows. void solout(int nr, double t, double y[], double yp[], int n, int *irtrn)
{
Output the y[] and yp[] values at nr-th step t. n is the order of equations.
The input value of irtrn will be 0, 1 or 2 in the first, intermediate or last call of solout, respectively. irtrn also serves to interrupt the integration. If irtrn is set to the negative value in solout, the integration will be interrupted and exit with info = 2.
}
|
| [in] | iout | Switch for calling the subroutine solout.
= 0: solout is never called
= 1: solout is used for output
(nrdens) must be specified in iwork[4].
(For other values, iout = 0 is assumed.) |
| [in,out] | work[] | Array work[lwork]
Work array.
work[0] to work[19] serve as parameters for the program. If the input parameter is set to 0, the default parameter value defined for each parameter will be loaded.
[in]
work[0]: Initial step size (default = estimated by the code)
work[1]: Maximum step size (default = tout - t)
work[3] and work[4]: Parameters for step size selection. The new step size is chosen subject to the restriction work[3] <= hnew/hold <= work[4]. (default work[3] = 0.2, work[4] = 10)
work[7]: The safety factor in step size prediction (0.0001 < work[7] < 1). If work[7] <= 0.0001 or work[7] >= 1, default value is used. (default = 0.9)
[out]
work[0]: Last step size |
| [in] | lwork | Size of array work[]. (lwork >= 8*n + 20)
If lwork < 0, abs(lwork) will be used and work[0] to work[19] will be initialized to zeros. |
| [in,out] | iwork[] | Array iwork[liwork]
Work array.
iwork[0] to iwork[19] serve as parameters for the program. If the input parameter is set to 0, the default parameter value defined for each parameter will be loaded.
[in]
iwork[1]: Maximum number of allowed steps. (default = 10000) (if iwork[1] < 0, default value will be used)
iwork[12]: Specifies when iwork[13] to iwork[17] are reset to zero.
= 0: Reset whenever this routine is called.
!= 0: Reset if info = 0.
[out]
iwork[13]: nfcn = number of function evaluations
iwork[15]: nstep = number of computed steps
iwork[16]: naccept = number of accepted steps
iwork[17]: nreject = number of rejected steps (due to error test) |
| [in] | liwork | Size of array iwork[]. (liwork >= 20)
If liwork < 0, abs(liwork) will be used and iwork[0] to iwork[19] will be initialized to zeros. |
| [in,out] | info | [in]
= 0: Initialize and start computation (Solve new problem).
= 1: Continue computation with new tout value (Resume computation of previous call).
[out]
= -1: The argument n had an illegal value. (n < 1)
= -7: The argument rtol had an illegal value (rtol < 0 or rtol[i] < 0)
= -7: The argument rtol or atol had an illegal value (rtol = 0 and atol = 0, or rtol[i] = 0 and atol[i] = 0)
= -8: The argument atol had an illegal value (atol < 0 or atol[i] < 0)
= -13: The argument lwork had an illegal value (lwork too small)
= -15: The argument liwork had an illegal value (liwork too small)
= -16: the argument info had an illegal value (info != 0 nor 1)
= 1: Successful exit
= 2: Interrupted by solout (normal return)
= 11: Convergence failure (maximum number of steps exceeded, etc.). |
- Reference
- E. Hairer, S.P. Norsett and G. Wanner, "Solving Ordinary Differential Equations I", Springer Series in Computational Mathematics, Springer-Verlag (1987)
- J. R. Dormand and P. J. Prince, "New Runge-Kutta algorithms for numerical simulation in dynamical astoronomy", Celestial Mechanics 18, 223?232 (1978)
|