NUMAL Section 5.2.1.1.1.1.H
BEGIN SECTION : 5.2.1.1.1.1.H (February, 1979)
PROCEDURE : ARK.
AUTHOR: P.A. BEENTJES.
INSTITUTE: MATHEMATICAL CENTRE.
RECEIVED: 740510.
BRIEF DESCRIPTION:
ARK SOLVES AN INITIAL VALUE PROBLEM, FOR A SYSTEM OF FIRST ORDER
ORDINARY DIFFERENTIAL EQUATIONS . ARK IS RECOMMENDED FOR THE
INTEGRATION OF SEMI-DISCRETE PARABOLIC AND HYPERBOLIC
INITIAL-BOUNDARY PROBLEMS.
KEYWORDS:
INITIAL VALUE PROBLEM,
SEMI-DISCRETE PARABOLIC AND HYPERBOLIC PROBLEM.
CALLING SEQUENCE:
THE HEADING OF THE PROCEDURE READS:
"PROCEDURE" ARK (T, TE, M0, M, U, DERIVATIVE, DATA, OUT);
"INTEGER" M0, M; "REAL" T, TE; "ARRAY" U, DATA;
"PROCEDURE" DERIVATIVE, OUT;
"CODE" 33061;
ARK : INTEGRATES THE SYSTEM OF ORDINARY DIFFERENTIAL EQUATIONS
DU / DT = H(T, U), U = U0 AT T = T0.
THE MEANING OF THE FORMAL PARAMETERS IS:
T: <VARIABLE>;
THE INDEPENDENT VARIABLE T;
ENTRY: THE INITIAL VALUE T0;
EXIT : THE FINAL VALUE TE;
TE: <ARITHMETIC EXPRESSION>;
THE FINAL VALUE OF T (TE >= T);
M0,M: <ARITHMETIC EXPRESSION>;
INDICES OF THE FIRST AND LAST EQUATION OF THE SYSTEM;
U: <ARRAY IDENTIFIER>;
"ARRAY" U[M0 : M];
ENTRY: THE INITIAL VALUES OF THE SOLUTION OF THE SYSTEM OF
DIFFERENTIAL EQUATIONS AT T = T0;
EXIT : THE VALUES OF THE SOLUTION AT T = TE;
DERIVATIVE: <PROCEDURE IDENTIFIER>;
THE HEADING OF THIS PROCEDURE READS:
"PROCEDURE" DERIVATIVE(T, V); "REAL" T; "ARRAY" V;
THIS PROCEDURE PERFORMS AN EVALUATION OF THE RIGHT HAND
SIDE OF THE SYSTEM WITH DEPENDENT VARIABLES V[M0 : M] AND
INDEPENDENT VARIABLE T; UPON COMPLETION OF DERIVATIVE, THE
RIGHT HAND SIDE SHOULD BE OVERWRITTEN ON V[M0 : M];
DATA: <ARRAY IDENTIFIER>;
"ARRAY" DATA[1 : 10 + DATA[1]];
IN ARRAY DATA ONE SHOULD GIVE:
DATA[1]: THE NUMBER OF EVALUATIONS OF H(T, U) PER
INTEGRATION STEP(DATA[1] >= DATA[2]);
DATA[2]: THE ORDER OF ACCURACY OF THE METHOD (DATA[2]<=3);
DATA[3]: STABILITY BOUND(SEE REFERENCE [3]);
DATA[4]: THE SPECTRAL RADIUS OF THE JACOBIAN MATRIX WITH
RESPECT TO THOSE EIGENVALUES, WHICH ARE LOCATED
IN THE NON-POSITIVE HALF PLANE;
DATA[5]: THE MINIMAL STEPSIZE;
DATA[6]: THE ABSOLUTE TOLERANCE;
DATA[7]: THE RELATIVE TOLERANCE;
IF BOTH DATA[6] AND DATA[7] ARE NEGATIVE, THE
INTEGRATION IS PERFORMED WITH A CONSTANT STEP
DATA[5];
DATA[8]: DATA[8] SHOULD BE 0 IF ARK IS CALLED FOR
A FIRST TIME; FOR CONTINUED INTEGRATION
DATA[8] SHOULD NOT BE CHANGED;
DATA[11], ..., DATA[10 + DATA[1]]: POLYNOMIAL COEFFICIENTS
(SEE REFERENCE [3]);
AFTER EACH STEP THE FOLLOWING BY-PRODUCTS ARE DELIVERED:
DATA[8]: THE NUMBER OF INTEGRATION STEPS PERFORMED;
DATA[9]: AN ESTIMATION OF THE LOCAL ERROR LAST MADE;
DATA[10]: INFORMATIVE MESSAGES:
DATA[10] = 0: NO DIFFICULTIES;
DATA[10] = 1: MINIMAL STEPLENGTH EXCEEDS THE
STEPLENGTH PRESCRIBED BY STABILITY
THEORY, I.E. DATA[5] > DATA[3] / DATA[4];
(TERMINATION OF ARK);
DECREASE MINIMAL STEPLENGTH;
IF NECESSARY, DATA[I],I = 4(1)7, CAN BE UPDATED(AFTER EACH
STEP) BY MEANS OF PROCEDURE OUT;
OUT: <PROCEDURE IDENTIFIER>;
THE HEADING OF THIS PROCEDURE READS:
"PROCEDURE" OUT;
AFTER EACH INTEGRATION STEP PERFORMED INFORMATION CAN BE
OBTAINED OR UPDATED BY THIS PROCEDURE, E.G. THE VALUES OF
T, U[M0 : M] AND DATA[I], I = 4(1)10.
DATA AND RESULTS:
FOR THE INDICES M0 AND M THE FOLLOWING REMARKS CAN BE MADE:
WHEN THE METHOD OF LINES IS APPLIED TO HYPERBOLIC DIFFERENTIAL
EQUATIONS THE NUMBER OF RELEVANT ORDINARY DIFFERENTIAL EQUATIONS
DECREASES DURING THE INTEGRATION PROCESS; IN PROCEDURE ARK
THIS MAY BE REALIZED BY INTEGERS M0 AND M, WHICH ARE
DEFINED AS FUNCTIONS OF THE NUMBER OF RIGHT HAND SIDE EVALUATIONS.
A SELECTION OF POSSIBLE ENTRIES FOR ARRAY DATA (DEPENDENT ON THE
KIND OF INITIAL VALUE PROBLEM) IS GIVEN IN REFERENCE [4],SECTION 8.
PROCEDURES USED:
INIVEC = CP31010,
MULVEC = CP31020,
DUPVEC = CP31030,
VECVEC = CP34010,
ELMVEC = CP34020,
DECSOL = CP34301.
REQUIRED CENTRAL MEMORY : CIRCA 75 + 2 * (M - M0) MEMORY PLACES.
METHOD AND PERFORMANCE:
ARK IS AN IMPLEMENTATION OF LOW ORDER STABILIZED RUNGE KUTTA
METHODS (SEE REFERENCE [1]);
AUTOMATIC STEPSIZE CONTROL IS PROVIDED BUT STEP-REJECTION HAS BEEN
EXCLUDED IN ORDER TO SAVE STORAGE;
BECAUSE OF ITS LIMITED STORAGE REQUIREMENTS AND ADAPTIVE STABILITY
FACILITIES THE METHOD IS WELL SUITED FOR THE SOLUTION OF INITIAL
BOUNDARY VALUE PROBLEMS FOR PARTIAL DIFFERENTIAL EQUATIONS;
NUMERICAL RESULTS, OBTAINED WITH A SLIGHTLY DIFFERENT
IMPLEMENTATION CAN BE FOUND IN REFERENCE [2].
REFERENCES:
[1]. P.J. VAN DER HOUWEN.
STABILIZED RUNGE KUTTA METHOD WITH LIMITED
STORAGE REQUIREMENTS.
MATH. CENTR. REPORT TW 124/71;
[2]. P.A. BEENTJES.
AN ALGOL 60 VERSION OF STABILIZED RUNGE KUTTA
METHODS (DUTCH).
MATH. CENTR. REPORT NR 23/72;
[3]. P.J. VAN DER HOUWEN, J. KOK.
NUMERICAL SOLUTION OF A MINIMAX PROBLEM.
MATH. CENTR. REPORT TW 123/71;
[4]. P.J. VAN DER HOUWEN ET AL.
ONE STEP METHODS FOR LINEAR INITIAL VALUE PROBLEMS, I.I.I.
NUMERICAL EXAMPLES,
MATH. CENTR. REPORT TW 130/71.
EXAMPLE OF USE:
THE VALUES OF
1. Y(1) AND Y(2) OF THE INITIAL VALUE PROBLEM
DY / DX = Y - 2 * X / Y, Y(0) = 1
AND
2. U(.6, 0) OF THE CAUCHY PROBLEM (SEE REFERENCE [2]):
DU / DT = .5 * DU / DX, U(0, X) = EXP(-X * X)
MAY BE OBTAINED BY THE FOLLOWING PROGRAM:
"BEGIN" "INTEGER" M0, M, I; "REAL" T, TE, DAT;
"ARRAY" Y[1 : 1], U[-150 : 150], DATA[1 : 14];
"PROCEDURE" DER1(T, V); "REAL" T; "ARRAY" V;
V[1]:= V[1] - 2 * T / V[1];
"PROCEDURE" DER2(T, V); "REAL" T; "ARRAY" V;
"BEGIN" "INTEGER" J; "REAL" V1, V2, V3;
V2:= V[M0]; M0:= M0 + 1; M:= M - 1; V3:= V[M0];
"FOR" J:= M0 "STEP" 1 "UNTIL" M "DO"
"BEGIN" V1:= V2; V2:= V3; V3:= V[J + 1];
V[J]:= 250 * (V3 - V1) / 3
"END"
"END" DER2;
"PROCEDURE" OUT1;
"IF" T = TE "THEN"
"BEGIN" "IF" T = 1 "THEN" OUTPUT(61, "("/, "(" PROBLEM 1")", //,
"(" X NUMBER OF INTEGRATION STEPS Y(COMPUTED) Y(EXACT)")",
//")");
OUTPUT(61, "("ZD, 13ZD,12B,2(-3ZD.7D),"("...")", /")",
T, DATA[8], Y[1], SQRT(2 * T + 1));
TE:= 2
"END" OUT1;
"PROCEDURE" OUT2;
"IF" T = .6 "THEN"
OUTPUT(61, "("//, "(" PROBLEM 2")", //,
"(" NUMBER OF DERIVATIVE CALLS")",
"(" U(.6, 0)COMPUTED U(.6, 0)EXACT")", //, 13ZD,
2(-10Z.7D), "("...")"")", DATA[1] * DATA[8], U[0], EXP(-.09));
I:= 1;
"FOR" DAT:= 3, 3, 1, 1, "-3, "-6, "-6, 0, 0, 0, 1, .5, 1 / 6 "DO"
"BEGIN" DATA[I]:= DAT; I:= I + 1 "END";
T:= 0; Y[1]:= 1; TE:= 1;
ARK(T, TE, 1, 1, Y, DER1, DATA, OUT1);
I:= 1;
"FOR" DAT:= 4, 3, SQRT(8), 500 / 3, DATA[3] / DATA[4], -1, -1,
0, 0, 0, 1, .5, 1 / 6, 1 / 24 "DO"
"BEGIN" DATA[I]:= DAT; I:= I + 1 "END";
M0:= -150; M:= 150; T:= 0; U[0]:= 1;
"FOR" I:= 1 "STEP" 1 "UNTIL" M "DO"
U[I]:= U[-I]:= EXP(-(.003 * I) ** 2);
ARK(T, .6, M0, M, U, DER2, DATA, OUT2)
"END"
THIS PROGRAM DELIVERS:
PROBLEM 1
X NUMBER OF INTEGRATION STEPS Y(COMPUTED) Y(EXACT)
1 38 1.7320535 1.7320508...
2 56 2.2360928 2.2360680...
PROBLEM 2
NUMBER OF DERIVATIVE CALLS U(.6, 0)COMPUTED U(.6, 0)EXACT
144 .9139326 .9139312...
SOURCE TEXT(S):
"CODE" 33061;