NUMAL Section 5.1.1.2.2

BEGIN SECTION : 5.1.1.2.2 (October, 1974)

AUTHOR: J.C.P.BUS.

INSTITUTE: MATHEMATICAL CENTRE.

RECEIVED: 740218.

BRIEF DESCRIPTION:

    THIS  SECTION  CONTAINS  TWO  PROCEDURES  FOR  SOLVING  SYSTEMS  OF
    NON-LINEAR EQUATIONS, OF WHICH THE JACOBIAN IS KNOWN TO BE A BAND
    MATRIX.
    QUANEWBND  ASKS FOR AN APPROXIMATION OF THE JACOBIAN AT THE INITIAL
    GUESS.
    QUANEWBND1  CALCULATES AN APPROXIMATION  OF THE JACOBIAN AT THE
    INITIAL GUESS, USING FORWARD DIFFERENCES.

KEYWORDS:

    NON-LINEAR EQUATIONS,
    SYSTEMS OF EQUATIONS,
    NO EXPLICIT JACOBIAN.


SUBSECTION: QUANEWBND.

CALLING SEQUENCE:

    THE HEADING OF THIS PROCEDURE READS :

    "PROCEDURE" QUANEWBND(N, LW, RW, X, F, JAC, FUNCT, IN, OUT);
    "VALUE" N, LW, RW; "INTEGER" N, LW, RW;
    "ARRAY" X, F, JAC, IN, OUT; "BOOLEAN" "PROCEDURE" FUNCT;
    "CODE" 34430;

    THE MEANING OF THE FORMAL PARAMETERS IS:
    N:      <ARITHMETIC EXPRESSION>;
            THE  NUMBER  OF   INDEPENDENT  VARIABLES;   THE  NUMBER  OF
            EQUATIONS SHOULD ALSO BE EQUAL TO N;
    LW:     <ARITHMETIC EXPRESSION>;
            THE NUMBER OF CODIAGONALS  TO THE LEFT OF THE MAIN DIAGONAL
            OF THE JACOBIAN;
    RW:     <ARITHMETIC EXPRESSION>;
            THE NUMBER OF CODIAGONALS TO THE RIGHT OF THE MAIN DIAGONAL
            OF THE JACOBIAN;
    X:      <ARRAY IDENTIFIER>;
            "ARRAY" X[1:N];
            ENTRY:  AN INITIAL ESTIMATE  OF THE SOLUTION  OF THE SYSTEM
                    THAT HAS TO BE SOLVED;
            EXIT:   THE CALCULATED SOLUTION OF THE SYSTEM;
    F:      <ARRAY IDENTIFIER>;
            "ARRAY" F[1:N];
            ENTRY:  THE  VALUES  OF  THE  FUNCTION  COMPONENTS  AT  THE
                    INITIAL GUESS;
            EXIT:   THE  VALUES  OF  THE  FUNCTION  COMPONENTS  AT  THE
                    CALCULATED SOLUTION;
    JAC:    <ARRAY IDENTIFIER>;
            "ARRAY" JAC[1:(LW + RW) * (N - 1) + N];
            ENTRY:  AN APPROXIMATION  OF  THE JACOBIAN  AT  THE INITIAL
                    ESTIMATE OF THE SOLUTION;
            EXIT:   AN APPROXIMATION OF THE JACOBIAN  AT THE CALCULATED
                    SOLUTION;
            AN APPROXIMATION OF THE (I, J)-TH  ELEMENT  OF THE JACOBIAN
            IS GIVEN IN JAC[(LW + RW) * (I - 1) + J], FOR I = 1, ..., N
            AND J = MAX(1, I - LW), ..., MIN(N, I + RW);

    FUNCT:  <PROCEDURE IDENTIFIER>;
            THE HEADING OF THIS PROCEDURE READS :
            "BOOLEAN" "PROCEDURE" FUNCT(N, L, U, X, F);
            "VALUE" N, L, U; "INTEGER" N, L, U; "ARRAY" X, F;
            THE MEANING OF THE FORMAL PARAMETERS IS :
            N:  <ARITHMETIC EXPRESSION>;
                THE NUMBER OF INDEPENDENT VARIABLES OF THE FUNCTION  F;
            L, U: <ARITHMETIC EXPRESSION>;
                LOWER  AND  UPPER  BOUND   OF  THE  FUNCTION  COMPONENT
                SUBSCRIPT;
            X:  <ARRAY IDENTIFIER>;
                THE INDEPENDENT VARIABLES ARE GIVEN IN X[1:N];
            F:  <ARRAY IDENTIFIER>;
                AFTER  A CALL OF  FUNCT  THE FUNCTION COMPONENTS  F[I],
                I = L, ..., U, SHOULD BE GIVEN IN F[L:U];
            IF  THE VALUE  OF  THE PROCEDURE IDENTIFIER  EQUALS  FALSE,
            THEN THE EXECUTION OF  QUANEWBND  WILL BE TERMINATED, WHILE
            THE VALUE OF OUT[5] IS SET EQUAL TO 2;
    IN:     <ARRAY IDENTIFIER>;
            "ARRAY" IN[0:4];
            ENTRY :
            IN  THIS  AUXILIARY  ARRAY  ONE  SHOULD  GIVE THE FOLLOWING
            VALUES FOR CONTROLLING THE PROCESS:
            IN[0]:  THE MACHINE PRECISION;
            IN[1]:  THE RELATIVE PRECISION ASKED FOR;
            IN[2]:  THE ABSOLUTE PRECISION  ASKED FOR;  IF  THE  VALUE,
                    DELIVERED IN  OUT[5]  EQUALS ZERO,  THEN  THE  LAST
                    CORRECTION VECTOR  D, SAY,  WHICH  IS A MEASURE FOR
                    THE ERROR IN THE SOLUTION,  SATIFIES THE INEQUALITY
                    NORM(D) <= NORM(X) * IN[1] + IN[2],
                    WHEREBY X DENOTES THE CALCULATED SOLUTION, GIVEN IN
                    ARRAY  X  AND  NORM(.)  DENOTES THE EUCLIDIAN NORM;
                    HOWEVER,WE CAN NOT GUARANTEE THAT THE TRUE ERROR IN
                    THE SOLUTION SATISFIES THIS INEQUALITY,  ESPECIALLY
                    IF  THE  JACOBIAN  IS  (NEARLY)  SINGULAR   AT  THE
                    SOLUTION;
            IN[3]:  THE  MAXIMUM  VALUE  OF  THE NORM  OF  THE RESIDUAL
                    VECTOR ALLOWED; IF  OUT[5] = 0,  THEN THIS RESIDUAL
                    VECTOR F, SAY, SATIFIES: NORM(F) <= IN[3];
            IN[4]:  THE   MAXIMUM   NUMBER    OF   FUNCTION   COMPONENT
                    EVALUATIONS ALLOWED;  L - U + 1  FUNCTION COMPONENT
                    EVALUATIONS    ARE    COUNTED    EACH    CALL    OF
                    FUNCT(N, L, U, X, F); IF OUT[5]=1, THEN THE PROCESS
                    IS TERMINATED,  BECAUSE  THE NUMBER  OF EVALUATIONS
                    EXCEEDED THE VALUE GIVEN IN IN[4];

    OUT:    <ARRAY IDENTIFIER>;
            "ARRAY" OUT[1:5];
            EXIT :
            OUT[1]: THE EUCLIDIAN NORM OF THE LAST STEP ACCEPTED;
            OUT[2]: THE EUCLIDIAN NORM  OF  THE RESIDUAL VECTOR  AT THE
                    CALCULATED SOLUTION;
            OUT[3]: THE  NUMBER   OF   FUNCTION  COMPONENT  EVALUATIONS
                    PERFORMED;
            OUT[4]: THE NUMBER OF ITERATIONS CARRIED OUT;
            OUT[5]: THE INTEGER VALUE DELIVERED IN  OUT[5]  GIVES  SOME
                    INFORMATION  ABOUT THE TERMINATION  OF THE PROCESS;
                    OUT[5] = 0: THE PROCESS  IS TERMINATED  IN A NORMAL
                                WAY; THE LAST STEP AND  THE NORM OF THE
                                RESIDUAL VECTOR SATISFY  THE CONDITIONS
                                (SEE IN[2], IN[3]);
                    IF  OUT[5] ^= 0,  THEN  THE PROCESS  IS  TERMINATED
                    PREMATURALY;
                    OUT[5] = 1: THE   NUMBER   OF   FUNCTION  COMPONENT
                                EVALUATIONS  EXCEEDS THE VALUE GIVEN IN
                                IN[4];
                    OUT[5] = 2: A CALL OF  FUNCT  DELIVERED  THE  VALUE
                                FALSE;
                    OUT[5] = 3: THE  APPROXIMATION   OF   THE  JACOBIAN
                                MATRIX TURNS OUT TO BE SINGULAR.

PROCEDURES USED:

    MULVEC      = CP31020,
    DUPVEC      = CP31030,
    VECVEC      = CP34010,
    ELMVEC      = CP34020,
    DECSOLBND   = CP34322.

EXECUTION FIELD LENGTH:

    THE MAXIMUM NUMBER OF WORDS, NECESSARY  FOR THE ARRAYS  DECLARED IN
    QUANEWBND EQUALS  MAX(N * 3 + (N - 1) * (LW + RW), 4N).

RUNNING TIME: PROPORTIONAL TO  N * LW * ( LW + RW + 1).

LANGUAGE:   ALGOL 60.

METHOD AND PERFORMANCE:

    THE METHOD USED IN QUANEWBND IS THE SAME AS GIVEN IN [1];  THE SAME
    PROBLEMS HAVE BEEN TESTED AND THE RESULTS ARE  THE SAME  OR  BETTER
    THAN THOSE REPORTED IN  [1];  CITING  [1],  WE CAN SAY  THAT  "THIS
    METHOD OFFERS A USEFUL IF MODEST IMPROVEMENT UPON  NEWTON'S METHOD,
    BUT THIS IMPROVEMENT TENDS TO VANISH  AS THE NONLINEARITIES  BECOME
    MORE PRONOUNCED".

EXAMPLE OF USE: SEE QUANEWBND1 (THIS SECTION).

REFERENCES:

    [1] BROYDEN C.G.
        THE CONVERGENCE  OF  AN ALGORITHM  FOR SOLVING SPARSE NONLINEAR
        SYSTEMS.
        MATH. COMP., VOL.25 (1971).


SUBSECTION: QUANEWBND1.

CALLING SEQUENCE:

    THE HEADING OF THIS PROCEDURE READS :

    "PROCEDURE" QUANEWBND1(N, LW, RW, X, F, FUNCT, IN, OUT);
    "VALUE" N, LW, RW; "INTEGER" N, LW, RW;
    "ARRAY" X, F, IN, OUT; "BOOLEAN" "PROCEDURE" FUNCT;
    "CODE" 34431;

    THE MEANING OF THE FORMAL PARAMETERS IS:
    N:      <ARITHMETIC EXPRESSION>;
            THE  NUMBER  OF   INDEPENDENT  VARIABLES;   THE  NUMBER  OF
            EQUATIONS SHOULD ALSO BE EQUAL TO N;
    LW:     <ARITHMETIC EXPRESSION>;
            THE NUMBER OF CODIAGONALS  TO THE LEFT OF THE MAIN DIAGONAL
            OF THE JACOBIAN;
    RW:     <ARITHMETIC EXPRESSION>;
            THE NUMBER OF CODIAGONALS TO THE RIGHT OF THE MAIN DIAGONAL
            OF THE JACOBIAN;
    X:      <ARRAY IDENTIFIER>;
            "ARRAY" X[1:N];
            ENTRY:  AN INITIAL ESTIMATE  OF THE SOLUTION  OF THE SYSTEM
                    THAT HAS TO BE SOLVED;
            EXIT:   THE CALCULATED SOLUTION OF THE SYSTEM;
    F:      <ARRAY IDENTIFIER>;
            "ARRAY" F[1:N];
            EXIT:   THE  VALUES  OF  THE  FUNCTION  COMPONENTS  AT  THE
                    CALCULATED SOLUTION;

    FUNCT:  <PROCEDURE IDENTIFIER>;
            THE HEADING OF THIS PROCEDURE READS :
            "BOOLEAN" "PROCEDURE" FUNCT(N, L, U, X, F);
            "VALUE" N, L, U; "INTEGER" N, L, U; "ARRAY" X, F;
            THE MEANING OF THE FORMAL PARAMETERS IS :
            N:  <ARITHMETIC EXPRESSION>;
                THE NUMBER OF INDEPENDENT VARIABLES OF THE FUNCTION  F;
            L, U: <ARITHMETIC EXPRESSION>;
                LOWER  AND  UPPER  BOUND   OF  THE  FUNCTION  COMPONENT
                SUBSCRIPT;
            X:  <ARRAY IDENTIFIER>;
                THE INDEPENDENT VARIABLES ARE GIVEN IN X[1:N];
            F:  <ARRAY IDENTIFIER>;
                AFTER  A CALL OF  FUNCT  THE FUNCTION COMPONENTS  F[I],
                I = L, ..., U, SHOULD BE GIVEN IN F[L:U];
            IF  THE VALUE  OF  THE PROCEDURE IDENTIFIER  EQUALS  FALSE,
            THEN THE EXECUTION OF  QUANEWBND  WILL BE TERMINATED, WHILE
            THE VALUE OF OUT[5] IS SET EQUAL TO 2;
    IN:     <ARRAY IDENTIFIER>;
            "ARRAY" IN[0:4];
            ENTRY :
            IN  THIS  AUXILIARY  ARRAY  ONE  SHOULD  GIVE THE FOLLOWING
            VALUES FOR CONTROLLING THE PROCESS:
            IN[0]:  THE MACHINE PRECISION;
            IN[1]:  THE RELATIVE PRECISION ASKED FOR;
            IN[2]:  THE ABSOLUTE PRECISION  ASKED FOR;  IF  THE  VALUE,
                    DELIVERED IN  OUT[5]  EQUALS ZERO,  THEN  THE  LAST
                    CORRECTION VECTOR  D, SAY,  WHICH  IS A MEASURE FOR
                    THE ERROR IN THE SOLUTION,  SATIFIES THE INEQUALITY
                    NORM(D) <= NORM(X) * IN[1] + IN[2],
                    WHEREBY X DENOTES THE CALCULATED SOLUTION, GIVEN IN
                    ARRAY  X  AND  NORM(.)  DENOTES THE EUCLIDIAN NORM;
                    HOWEVER,WE CAN NOT GUARANTEE THAT THE TRUE ERROR IN
                    THE SOLUTION SATISFIES THIS INEQUALITY,  ESPECIALLY
                    IF  THE  JACOBIAN  IS  (NEARLY)  SINGULAR   AT  THE
                    SOLUTION;
            IN[3]:  THE  MAXIMUM  VALUE  OF  THE NORM  OF  THE RESIDUAL
                    VECTOR ALLOWED; IF  OUT[5] = 0,  THEN THIS RESIDUAL
                    VECTOR F, SAY, SATIFIES: NORM(F) <= IN[3];
            IN[4]:  THE   MAXIMUM   NUMBER    OF   FUNCTION   COMPONENT
                    EVALUATIONS ALLOWED;  L - U + 1  FUNCTION COMPONENT
                    EVALUATIONS    ARE    COUNTED    EACH    CALL    OF
                    FUNCT(N, L, U, X, F); IF OUT[5]=1, THEN THE PROCESS
                    IS TERMINATED,  BECAUSE  THE NUMBER  OF EVALUATIONS
                    EXCEEDED THE VALUE GIVEN IN IN[4];
            IN[5]:  THE  JACOBIAN  MATRIX  AT  THE  INITIAL  GUESS   IS
                    APPROXIMATED  USING  FORWARD DIFFERENCES,  WITH  AN
                    FIXED INCREMENT TO EACH VARIABLE  THAT  EQUALS  THE
                    VALUE GIVEN IN IN[5];

    OUT:    <ARRAY IDENTIFIER>;
            "ARRAY" OUT[1:5];
            EXIT :
            OUT[1]: THE EUCLIDIAN NORM OF THE LAST STEP ACCEPTED;
            OUT[2]: THE EUCLIDIAN NORM  OF  THE RESIDUAL VECTOR  AT THE
                    CALCULATED SOLUTION;
            OUT[3]: THE  NUMBER   OF   FUNCTION  COMPONENT  EVALUATIONS
                    PERFORMED;
            OUT[4]: THE NUMBER OF ITERATIONS CARRIED OUT;
            OUT[5]: THE INTEGER VALUE DELIVERED IN  OUT[5]  GIVES  SOME
                    INFORMATION  ABOUT THE TERMINATION  OF THE PROCESS;
                    OUT[5] = 0: THE PROCESS  IS TERMINATED  IN A NORMAL
                                WAY; THE LAST STEP AND  THE NORM OF THE
                                RESIDUAL VECTOR SATISFY  THE CONDITIONS
                                (SEE IN[2], IN[3]);
                    IF  OUT[5] ^= 0,  THEN  THE PROCESS  IS  TERMINATED
                    PREMATURALY;
                    OUT[5] = 1: THE   NUMBER   OF   FUNCTION  COMPONENT
                                EVALUATIONS  EXCEEDS THE VALUE GIVEN IN
                                IN[4];
                    OUT[5] = 2: A CALL OF  FUNCT  DELIVERED  THE  VALUE
                                FALSE;
                    OUT[5] = 3: THE  APPROXIMATION   OF   THE  JACOBIAN
                                MATRIX TURNS OUT TO BE SINGULAR.

PROCEDURES USED:

    JACOBNBNDF  = CP34439,
    QUANEWBND   = CP34430.

EXECUTION FIELD LENGTH:

    QUANEWBND1  DECLARES AN AUXILIARY ARRAY OF DIMENSION ONE  AND ORDER
    N + (N - 1) * (LW + RW).

RUNNING TIME: PROPORTIONAL TO  N * LW * ( LW + RW + 1).

LANGUAGE:   ALGOL 60.

METHOD AND PERFORMANCE:

    QUANEWBND1  USES  JACOBNBNDF  (SECTION  4.3.2.1)  TO  CALCULATE  AN
    INITIAL APPROXIMATION OF THE JACOBIAN MATRIX   AT THE INITIAL GUESS
    GIVEN  IN  X[1:N]  AND  SOLVES  THE  NONLINEAR  SYSTEM  BY  CALLING
    QUANEWBND (THIS SECTION).

EXAMPLE OF USE:

    LET THE FUNCTION F BE DEFINED BY (SEE [1]):
    F[1] = (3 - 2 * X[1]) * X[1] + 1 - 2 * X[2],
    F[I] = (3 - 2 * X[I]) * X[I] + 1 - X[I - 1] - 2 * X[I + 1],  I = 2,
    ..., N - 1,
    F[N] = (3 - 2 * X[N]) * X[N] + 1 - X[N - 1];
    LET AN INITIAL ESTIMATE OF THE SOLUTION OF THE SYSTEM  F(X) = 0  BE
    GIVEN BY X[I] = -1, I = 1, ..., N;  THEN THE FOLLOWING PROGRAM  MAY
    SOLVE THIS SYSTEM FOR N = 600 AND PRINTS SOME RESULTS.

    "BEGIN"
        "BOOLEAN" "PROCEDURE" FUN(N, L, U, X, F); "VALUE" N, L, U;
        "INTEGER" N, L, U; "ARRAY" X, F;
        "BEGIN" "INTEGER" I; "REAL" X1, X2, X3;
            X1:= "IF" L = 1 "THEN" 0 "ELSE" X[L - 1]; X2:= X[L];
            X3:= "IF" L = N "THEN" 0 "ELSE" X[L + 1];
            "FOR" I:= L "STEP" 1 "UNTIL" U "DO"
            "BEGIN" F[I]:= (3 - 2 * X2) * X2 + 1 - X1 - X3 * 2;
                X1:= X2; X2:= X3;
                X3:= "IF" I <= N - 2 "THEN" X[I + 2] "ELSE" 0
            "END"; FUN:= "TRUE"
        "END" FUN;

        "INTEGER" I; "ARRAY" X, F[1:600], IN[0:5], OUT[1:5];
        "FOR" I := 1 "STEP" 1 "UNTIL" 600 "DO" X[I]:= -1;
        IN[0]:= "-14; IN[1]:= IN[2]:= IN[3]:= "-6; IN[4]:= 20000;
        IN[5]:= 0.001;
        QUANEWBND1(600, 1, 1, X, F, FUN, IN, OUT);
        OUTPUT(71, "("//,"(" NORM RESIDUALVECTOR: ")"+.15D"+3D,/,
        "(" LENGTH OF LAST STEP: ")"+.15D"+3D,/,
        "(" NUMBER OF FUNCTION COMPONENT EVALUATIONS: ")"5ZD,/,
        "(" NUMBER OF ITERATIONS: ")"4ZD/,"("REPORT: ")"D/")",
        OUT[2], OUT[1], OUT[3], OUT[4], OUT[5])
    "END"

    RESULTS:

    NORM RESIDUALVECTOR: +.221010684482660"-006
    LENGTH OF LAST STEP: +.302712457332660"-006
    NUMBER OF FUNCTION COMPONENT EVALUATIONS:   6598
    NUMBER OF ITERATIONS:     7
    REPORT: 0

REFERENCES:

    [1] BROYDEN C.G.
        THE CONVERGENCE  OF  AN ALGORITHM  FOR SOLVING SPARSE NONLINEAR
        SYSTEMS.
        MATH. COMP., VOL.25 (1971).

SOURCE TEXT(S):
"CODE" 34430;

"CODE" 34431;