NUMAL Section 5.2.1.1.1.1.B
BEGIN SECTION : 5.2.1.1.1.1.B (February, 1979)
PROCEDURE : RKE.
AUTHOR: P.A. BEENTJES.
INSTITUTE: MATHEMATICAL CENTRE.
RECEIVED: 740520.
BRIEF DESCRIPTION:
RKE SOLVES AN INITIAL VALUE PROBLEM FOR A SYSTEM OF FIRST
ORDER ORDINARY DIFFERENTIAL EQUATIONS DY / DX = F(X,Y).
THE SYSTEM IS SUPPOSED TO BE NON-STIFF.
KEYWORDS:
INITIAL VALUE PROBLEM.
SYSTEM OF FIRST ORDER ORDINARY DIFFERENTIAL EQUATIONS.
CALLING SEQUENCE:
THE HEADING OF THE PROCEDURE READS:
"PROCEDURE" RKE (X, XE, N, Y, DER, DATA, FI, OUT);
"VALUE" N, FI; "INTEGER" N; "REAL" X, XE; "BOOLEAN" FI;
"ARRAY" Y, DATA;
"PROCEDURE" DER, OUT;
"CODE" 33033;
RKE INTEGRATES THE SYSTEM OF ORDINARY DIFFERENTIAL EQUATIONS
DY / DX = F(X, Y), FROM X = X0 TO X = XE WHERE Y(X0) = Y0.
THE MEANING OF THE FORMAL PARAMETERS IS:
X: <VARIABLE>;
THE INDEPENDENT VARIABLE;
ENTRY: THE INITIAL VALUE X0;
XE: <ARITHMETIC EXPRESSION>;
THE FINAL VALUE OF X;
N: <ARITHMETIC EXPRESSION>;
THE NUMBER OF EQUATIONS OF THE SYSTEM;
Y: <ARRAY IDENTIFIER>;
"ARRAY" Y[1 : N];
THE DEPENDENT VARIABLES;
ENTRY: THE INITIAL VALUES AT X = X0;
EXIT : THE VALUES OF THE SOLUTION AT X = XE;
DER: <PROCEDURE IDENTIFIER>;
THE HEADING OF THIS PROCEDURE READS:
"PROCEDURE" DER (T, V); "VALUE" T; "REAL" T; "ARRAY" V;
THIS PROCEDURE PERFORMS AN EVALUATION OF THE RIGHT HAND
SIDE OF THE SYSTEM WITH DEPENDENT VARIABLES V[1 : N] AND
INDEPENDENT VARIABLE T; UPON COMPLETION OF DER
THE RIGHT HAND SIDE SHOULD BE OVERWRITTEN ON V[1 : N];
DATA: <ARRAY IDENTIFIER>;
"ARRAY" DATA[1 : 6];
IN ARRAY DATA ONE SHOULD GIVE:
DATA[1]: THE RELATIVE TOLERANCE;
DATA[2]: THE ABSOLUTE TOLERANCE;
AFTER EACH STEP DATA[3:6] CONTAINS :
DATA[3]: THE STEPLENGTH USED FOR THE LAST STEP;
DATA[4]: THE NUMBER OF INTEGRATION STEPS PERFORMED;
DATA[5]: THE NUMBER OF INTEGRATION STEPS REJECTED;
DATA[6]: THE NUMBER OF INTEGRATION STEPS SKIPPED;
IF UPON COMPLETION OF RKE DATA[6] > 0 ,
RESULTS SHOULD BE CONSIDERED MOST CRITICALLY;
FI: <BOOLEAN EXPRESSION>;
IF FI = "TRUE" THE INTEGRATION STARTS AT X0 WITH A
TRIAL STEP XE - X0; IF FI = "FALSE" THE INTEGRATION IS
CONTINUED WITH A STEPLENGTH DATA[3] * SIGN(XE - X0);
OUT: <PROCEDURE IDENTIFIER>;
THE HEADING OF THIS PROCEDURE READS:
"PROCEDURE" OUT;
AFTER EACH INTEGRATION STEP PERFORMED OUT CAN BE USED TO
OBTAIN INFORMATION FROM THE SOLUTION PROCESS, E.G. THE VALUES
OF X, Y[1 : N] AND DATA[3 : 6]. POUT CAN ALSO BE USED TO
UPDATE DATA.
PROCEDURES USED : NONE.
REQUIRED CENTRAL MEMORY:
CIRCA 5 * N MEMORY PLACES.
METHOD AND PERFORMANCE:
THE METHOD UPON WHICH THE PROCEDUDRE IS BASED, IS A MEMBER OF A
CLASS OF FIFTH ORDER RUNGE KUTTA FORMULAS PRESENTED IN REFERENCE
[1]. AUTOMATIC STEPSIZE CONTROL IS IMPLEMENTED IN A WAY AS PROPOSED
IN REFERENCE [2]. FOR TESTRESULTS AND FURTHER INFORMATION
SEE REFERENCE [3].
REFERENCES:
[1]. R. ENGLAND.
ERROR ESTIMATES FOR RUNGE KUTTA TYPE SOLUTIONS TO SYSTEMS
OF ORDINARY DIFFERENTIAL EQUATIONS.
THE COMPUTER JOURNAL , VOLUME 12, P 166 - 169, 1969.
[2]. J.A. ZONNEVELD.
AUTOMATIC NUMERICAL INTEGRATION.
MATH. CENTRE TRACT 8(1970).
[3]. P.A. BEENTJES.
SOME SPECIAL FORMULAS OF THE ENGLAND CLASS OF FIFTH ORDER
RUNGE - KUTTA SCHEMES.
MATH. CENTRE REPORT NW 24/74.
EXAMPLE OF USE:
THE SOLUTION AT T = 1 AND T = -1 OF THE SYSTEM
DX / DT = Y - Z,
DY / DT = X * X + 2 * Y + 4 * T,
DZ / DT = X * X + 5 * X + Z * Z + 4 * T,
WITH X = Y = 0 AND Z = 2 AT T = 0,
CAN BE OBTAINED BY THE FOLLOWING PROGRAM:
"BEGIN" "REAL" T, TE; "ARRAY" Y[1 : 3], DATA[1 : 6];
"PROCEDURE" RHS(T, Y); "VALUE" T; "REAL" T; "ARRAY" Y;
"BEGIN" "REAL" XX, YY, ZZ;
XX:= Y[1]; YY:= Y[2]; ZZ:= Y[3];
Y[1]:= YY - ZZ;
Y[2]:= XX * XX + 2 * YY + 4 * T;
Y[3]:= XX * (XX + 5) + 2 * ZZ + 4 * T
"END" RHS;
"PROCEDURE" INFO;
"IF" T = TE "THEN"
"BEGIN" "REAL" ET, T2, AEX, AEY, AEZ, REX, REY, REZ;
ET:= EXP(T); T2:= 2 * T;
REX:= -ET * SIN(T2); AEX:= REX - Y[1]; REX:= ABS(AEX / REX);
REY:= ET * ET * (8 + 2 * T2 - SIN(2 * T2)) / 8 - T2 - 1;
REZ:= ET * (SIN(T2) + 2 * COS(T2)) + REY;
AEY:= REY - Y[2]; REY:= ABS(AEY / REY); AEZ:= REZ - Y[3];
REZ:= ABS(AEZ / REZ);
OUTPUT(61, "(""(" T = ")", +D, //,
"(" RELATIVE AND ABSOLUTE ERRORS IN X, Y AND Z :")", //,
"(" RE(X) RE(Y) RE(Z) AE(X) AE(Y) AE(Z)")", //,
6(B,-.2D"+D), //,
"(" NUMBER OF INTEGRATION STEPS PERFORMED :")",4ZD,/,
"(" NUMBER OF INTEGRATION STEPS SKIPPED :")",4ZD,/,
"(" NUMBER OF INTEGRATION STEPS REJECTED :")",4ZD,///")",
T, REX, REY, REZ, ABS(AEX), ABS(AEY), ABS(AEZ),
DATA[4], DATA[6], DATA[5])
"END" INFO;
TE:= 1;
LEFT:
Y[1]:= Y[2]:= 0; Y[3]:= 2; T:=0;
DATA[1]:= DATA[2]:= "-5;
RKE(T, TE, 3, Y, RHS, DATA, "TRUE", INFO);
"IF" TE = 1 "THEN" "BEGIN" TE:= -1; "GOTO" LEFT "END"
"END"
THIS PROGRAM DELIVERS:
T = +1
RELATIVE AND ABSOLUTE ERRORS IN X, Y AND Z :
RE(X) RE(Y) RE(Z) AE(X) AE(Y) AE(Z)
.37"-6 .15"-5 .13"-5 .91"-6 .13"-4 .11"-4
NUMBER OF INTEGRATION STEPS PERFORMED : 9
NUMBER OF INTEGRATION STEPS SKIPPED : 0
NUMBER OF INTEGRATION STEPS REJECTED : 5
T = -1
RELATIVE AND ABSOLUTE ERRORS IN X, Y AND Z :
RE(X) RE(Y) RE(Z) AE(X) AE(Y) AE(Z)
.22"-6 .52"-7 .19"-6 .75"-7 .55"-7 .77"-7
NUMBER OF INTEGRATION STEPS PERFORMED : 10
NUMBER OF INTEGRATION STEPS SKIPPED : 0
NUMBER OF INTEGRATION STEPS REJECTED : 7
SOURCE TEXT(S):
"CODE" 33033;