NUMAL Section 5.2.1.1.1.1.E

BEGIN SECTION : 5.2.1.1.1.1.E (February, 1979)

PROCEDURE : RK5NA.

AUTHOR: J.A.ZONNEVELD.

CONTRIBUTORS: M.BAKKER AND I.BRINK.

INSTITUTE: MATHEMATICAL CENTRE.

RECEIVED: 730715.

BRIEF DESCRIPTION:

    RK5NA SOLVES AN INITIAL VALUE PROBLEM FOR A SYSTEM OF FIRST ORDER
    ORDINARY DIFFERENTIAL EQUATIONS  DY / DX = F(X,Y), OF WHICH THE
    DERIVATIVE COMPONENTS ARE SUPPOSED TO BECOME LARGE, E.G. IN THE
    NEIGHBOURHOOD OF SINGULARITIES. RK5NA INTRODUCES THE ARC LENGTH
    AS INTEGRATION VARIABLE. THE SYSTEM IS SUPPOSED TO BE NON-STIFF.

KEYWORDS:

    INITIAL VALUE PROBLEM.
    SYSTEM OF FIRST ORDER ORDINARY DIFFERENTIAL EQUATIONS.
    LARGE DERIVATIVE COMPONENTS.

CALLING SEQUENCE:

    THE HEADING OF THE PROCEDURE READS:
    "PROCEDURE" RK5NA(X, XA, B, FXJ, J, E, D, FI, N, L, POS);
    "VALUE" FI, N, L, POS; "INTEGER" J, N, L; "BOOLEAN" FI, POS;
    "REAL" B, FXJ; "ARRAY" X, XA, E, D;
    "CODE" 33018;

    THE MEANING OF THE FORMAL PARAMETERS IS:

    X:  <ARRAY IDENTIFIER>;
        "ARRAY" X[0 : N];
        THE DEPENDENT VARIABLES;
        X[0], ..., X[N] CAN BE USED AS JENSEN PARAMETERS;
    XA:  <ARRAY IDENTIFIER>;
        "ARRAY" XA[0 : N];
        ENTRY: THE INITIAL VALUES OF X[J], J = 0, ..., N;
    B:  <ARITHMETIC EXPRESSION>;
        B DEPENDS ON X[0], ..., X[N];
        IF,WITHIN SOME TOLERANCE,B = 0 THEN THE INTEGRATION IS TER-
        MINATED; SEE ALSO THE EXPLANATION OF THE PARAMETER E;
    FXJ:  <ARITHMETIC EXPRESSION>;
        THE RIGHT HAND SIDE OF THE DIFFERENTIAL EQUATION;
        FXJ DEPENDS ON X[0], ..., X[N], J, GIVING THE VALUE OF
        DX[J] / DX[0];
    J:  <VARIABLE>;
        J IS USED AS A JENSEN PARAMETER TO DENOTE,IN THE ACTUAL
        PARAMETER CORRESPONDING TO FXJ,THE NUMBER OF THE FUNCTION
        REQUIRED;
    E:  <ARRAY IDENTIFIER>;
        "ARRAY" E[0 : 2 * N + 3];
        ENTRY:
        E[2 * J] AND E[2 * J + 1] ARE THE RELATIVE AND THE ABSOLUTE
        TOLERANCE,RESPECTIVELY,ASSOCIATED WITH X[J],J = 0,..., N, WHILE
        E[2 * N + 2] AND E[2 * N + 3] ARE THE ONES ASSOCIATED WITH B;
    D:  <ARRAY IDENTIFIER>;
        "ARRAY" D[1 : N + 3];
        AFTER COMPLETION OF EACH STEP WE HAVE:

            ABS(D[1])                   THE ARC LENGTH,
                D[2]                    THE STEP LENGTH,
                D[I + 3]                THE LATEST VALUE OF X[I],
                                        I = 0, ..., N;

    FI:  <BOOLEAN EXPRESSION>;
        IF FI = "TRUE" THEN THE INTEGRATION IS STARTED WITH INITIAL
        CONDITIONS X[I] = XA[I], I = 0, ..., N;
        IF FI = "FALSE" THEN THE INTEGRATION IS CONTINUED WITH
        X[I] = D[I + 3];
    N:  <ARITHMETIC EXPRESSION>;
        THE NUMBER OF EQUATIONS;

    L:  <ARITHMETIC EXPRESSION>;
        1 <= L <= N; SEE THE EXPLANATION OF POS;
    POS:  <BOOLEAN EXPRESSION>;
        IF FI = "TRUE" THEN THE INTEGRATION STARTS IN SUCH A WAY THAT
        X[L] INCREASES IF POS = "TRUE" AND DECREASES IF POS = "FALSE".
        IF FI = "FALSE" THEN POS IS IGNORED.

PROCEDURES USED : ZEROIN = CP34150.

REQUIRED CENTRAL MEMORY : CIRCA 9 * N  MEMORY PLACES.

METHOD AND PERFORMANCE :
    RK5NA IS BASED ON A 5-TH ORDER RUNGE-KUTTA  METHOD AND INTEGRATES
    THE SYSTEM  DX[J] / DX[0] = F(J,X[0],..,X[N]) / F(0,X[0],..,X[N]).
    THE ARC LENGTH IS INTRODUCED AS INTEGRATION VARIABLE.
    THE INTEGRATION PROCESS IS TERMINATED IF SOME CONDITION ON
    X[0], .. ,X[N], TO BE SUPPLIED BY THE USER, IS SATISFIED.
    RK5NA USES STEPLENGTH AND ERROR CONTROL. DETAILS ABOUT THE
    PROCEDURE AND THE UNDERLYING THEORY ARE GIVEN IN [1] ( RK5NA IS
    A SLIGHTLY  ADAPTED VERSION OF RK5N).

REFERENCES:

    [1]. J.A.ZONNEVELD,
         AUTOMATIC NUMERICAL INTEGRATION,
         MATH.CENTRE TRACT 8 (1970) .

EXAMPLE OF USE:

    THE VAN DER POL EQUATION IN THE PHASE PLANE

        DX[1] / DX[0] = (10*(1-X[0]**2)*X[1]-X[0])/X[1]

    CAN BE INTEGRATED BY THE PROCEDURE RK5NA; THE STARTING VALUES ARE
    X[0] = 2, X[1] = 0. THE INTEGRATION PROCEEDS UNTIL THE NEXT ZERO OF
    X[1],THEN IT CONTINUES UNTIL THE NEXT ZERO AND SO ON UNTIL THE
    FOURTH ZERO IS ENCOUNTERED.IN THE EXAMPLE THE OUTPUT IS GIVEN FOR
    THE TOLERANCES E[0]=E[1]=E[2]=E[3]= "-6, E[4]=E[5]= "-10.
    THE PROGRAM READS AS FOLLOWS:

    "BEGIN" "COMMENT" VAN DER POL IN THE PHASE PLANE;
        "INTEGER" J,K; "BOOLEAN" FIRST;
        "ARRAY" E[0:5],XA,X[0:1],D[1:4];
        "PROCEDURE" PRINT(X); "ARRAY" X;
        "BEGIN" OUTPUT(61,"("/B"("X[0]=")"+D.10D,
            10B"("X[1]=")"+D.10D,10B"("S=")"3D.10D")",X[0],
            X[1],ABS(D[1]))
        "END";
        "FOR" K:=0,1,2,3 "DO" E[K]:="-6 ; E[4]:=E[5]:="-10 ;  D[1]:=0;
        XA[0]:=2; XA[1]:=0; J:=0; PRINT(XA); AA:
        FIRST:=J=0;
        RK5NA(X,XA,X[1],"IF" K=0 "THEN" X[1] "ELSE"
        10*(1-X[0]**2)*X[1]-X[0],K,E,D,FIRST,1,1,"FALSE");;J:=J+1;
        PRINT(X); "IF" J<4 "THEN" "GO TO" AA
    "END"

    RESULTS:

    X[0]=+2.0000000000    X[1]=+0.0000000000    S=000.0000000000

    X[0]=-2.0142853657    X[1]=-0.0000000012    S=029.3873834087

    X[0]=+2.0142853659    X[1]=+0.0000000001    S=058.7884331939

    X[0]=-2.0142853659    X[1]=-0.0000000000    S=088.1894829781

    X[0]=+2.0142853659    X[1]=+0.0000000000    S=117.5905327623

SOURCE TEXT(S):
"CODE" 33018;