NUMAL Section 3.3.2.2.2

BEGIN SECTION : 3.3.2.2.2 (July, 1974)

AUTHOR   : C.G. VAN DER LAAN.

CONTRIBUTORS : H.FIOLET, C.G. VAN DER LAAN.

INSTITUTE : MATHEMATICAL CENTRE.

RECEIVED: 731016.

BRIEF DESCRIPTION :

    THIS SECTION CONTAINS THE PROCEDURES EIGVALCOM AND EIGCOM.
    EIGVALCOM CALCULATES THE EIGENVALUES OF A COMPLEX MATRIX AND
    EIGCOM CALCULATES THE EIGENVECTORS AS WELL.

KEYWORDS :

    EIGENVALUES,
    EIGENVECTORS,
    COMPLEX MATRICES,
    EQUILIBRATION,
    REDUCTION HESSENBERG FORM,
    HOUSEHOLDER TRANSFORMATION,
    QR-ITERATION.


SUBSECTION: EIGVALCOM.

CALLING SEQUENCE:

    THE HEADING OF THE PROCEDURE READS:
    "INTEGER" "PROCEDURE" EIGVALCOM(AR, AI, N, EM, VALR, VALI);
    "VALUE" N; "INTEGER" N; "ARRAY" AR, AI, EM, VALR, VALI;
    "CODE" 34374;

    THE MEANING OF THE FORMAL PARAMETERS IS:
    AR,AI:     <ARRAY IDENTIFIER>;
               "ARRAY" AR,AI[1:N,1:N];
               ENTRY:
               THE REAL PART AND THE IMAGINARY PART OF THE MATRIX MUST
               BE GIVEN IN THE ARRAYS AR AND AI, RESPECTIVELY;
               THE ELEMENTS OF THE ARRAYS AR AND AI ARE ALTERED;
    N:         <ARITHMETIC EXPRESSION>;
               THE ORDER OF THE GIVEN MATRIX;
    EM:        <ARRAY IDENTIFIER>;
               "ARRAY" EM[0:7];
               ENTRY:
               EM[0]: THE MACHINE PRECISION;
               EM[2]: THE RELATIVE TOLERANCE FOR THE QR-ITERATION;
                      (E.G. THE MACHINE PRECISION);
               EM[4]: THE MAXIMUM ALLOWED NUMBER OF QR-ITERATIONS;
                      (E.G. 10 * N);
               EM[6]: THE MAXIMUM ALLOWED NUMBER OF ITERATIONS FOR
                      EQUILIBRATING THE ORIGINAL MATRIX (E.G. N**2/2);
               EXIT:
               EM[1]: THE EUCLIDEAN NORM OF THE EQUILIBRATED MATRIX;
               EM[3]: THE MAXIMUM ABSOLUTE VALUE OF THE SUBDIAGONAL
                      ELEMENTS NEGLECTED IN THE QR-ITERATION;
               EM[5]: THE NUMBER OF QR-ITERATIONS PERFORMED;
                      EM[5]:=EM[4]+1 IN THE CASE EIGVALCOM^=0;
               EM[7]: THE NUMBER OF ITERATIONS PERFORMED FOR
                      EQUILIBRATING THE ORIGINAL MATRIX;
    VALR,VALI: <ARRAY IDENTIFIER>;
               "ARRAY" VALR,VALI[1:N];
               EXIT:
               THE REAL PART AND THE IMAGINARY PART OF THE CALCULATED
               EIGENVALUES ARE DELIVERED IN THE ARRAYS VALR AND VALI,
               RESPECTIVELY;

    EIGVALCOM:=0, PROVIDED THE QR-ITERATION IS COMPLETED WITHIN EM[4]
    ITERATIONS; OTHERWISE, EIGVALCOM:= THE NUMBER, K, OF EIGENVALUES
    NOT CALCULATED AND ONLY THE LAST N-K ELEMENTS OF THE ARRAYS VALR
    AND VALI ARE APPROXIMATE EIGENVALUES OF THE ORIGINAL MATRIX.

PROCEDURES USED:

    EQILBRCOM = CP34361,
    COMEUCNRM = CP34359,
    HSHCOMHES = CP34366,
    VALQRICOM = CP34372.

REQUIRED CENTRAL MEMORY: FIVE  REAL ARRAYS OF ORDER N AND ONE INTEGER
    ARRAY OF ORDER N ARE DECLARED.

RUNNING TIME: PROPORTIONAL TO N ** 2 * MAX(N,NUMBER OF ITERATIONS).

LANGUAGE: ALGOL 60.

METHOD AND PERFORMANCE: SEE EIGCOM (THIS SECTION).

EXAMPLE OF USE: SEE EIGCOM (THIS SECTION).


SUBSECTION: EIGCOM.

CALLING SEQUENCE:

    THE HEADING OF THE PROCEDURE READS:
    "INTEGER" "PROCEDURE" EIGCOM(AR, AI, N, EM, VALR, VALI, VR, VI);
    "VALUE" N; "INTEGER" N; "ARRAY" AR, AI, EM, VALR, VALI, VR, VI;
    "CODE" 34375;

    THE MEANING OF THE FORMAL PARAMETERS IS:
    AR,AI:     <ARRAY IDENTIFIER>;
               "ARRAY" AR,AI[1:N,1:N];
               ENTRY:
               THE REAL PART AND THE IMAGINARY PART OF THE MATRIX MUST
               BE GIVEN IN THE ARRAYS AR AND AI, RESPECTIVELY;
               THE ELEMENTS OF THE ARRAYS AR AND AI ARE ALTERED;
    N:         <ARITHMETIC EXPRESSION>;
               THE ORDER OF THE GIVEN MATRIX;
    EM:        <ARRAY IDENTIFIER>;
               "ARRAY" EM[0:7];
               ENTRY:
               EM[0]: THE MACHINE PRECISION;
               EM[2]: THE RELATIVE TOLERANCE FOR THE QR-ITERATION;
                      (E.G. THE MACHINE PRECISION);
               EM[4]: THE MAXIMUM ALLOWED NUMBER OF QR-ITERATIONS;
                      (E.G. 10 * N);
               EM[6]: THE MAXIMUM ALLOWED NUMBER OF ITERATIONS FOR
                      EQUILIBRATING THE ORIGINAL MATRIX (E.G. N**2/2);
               EXIT:
               EM[1]: THE EUCLIDEAN NORM OF THE EQUILIBRATED MATRIX;
               EM[3]: THE MAXIMUM ABSOLUTE VALUE OF THE SUBDIAGONAL
                      ELEMENTS NEGLECTED IN THE QR-ITERATION;
               EM[5]: THE NUMBER OF QR-ITERATIONS PERFORMED;
                      EM[5]:=EM[4]+1 IN THE CASE EIGCOM^=0;
               EM[7]: THE NUMBER OF ITERATIONS PERFORMED FOR
                      EQUILIBRATING THE ORIGINAL MATRIX;

    VALR,VALI: <ARRAY IDENTIFIER>;
               "ARRAY" VALR,VALI[1:N];
               EXIT:
               THE REAL PART AND THE IMAGINARY PART OF THE CALCULATED
               EIGENVALUES ARE DELIVERED IN THE ARRAYS VALR AND VALI,
               RESPECTIVELY;
    VR,VI:     <ARRAY IDENTIFIER>;
               "ARRAY" VR,VI[1:N,1:N];
               EXIT:
               THE EIGENVECTORS OF THE MATRIX;
               THE NORMALIZED EIGENVECTOR WITH REAL PART VR[1:N,J] AND
               IMAGINARY PART VI[1:N,J] CORRESPONDS TO THE EIGENVALUE
               VALR[J] + VALI[J] * I, J=1,...,N;

    EIGCOM:=0, PROVIDED THE QR-ITERATION IS COMPLETED WITHIN EM[4]
    ITERATIONS; OTHERWISE, EIGCOM:= THE NUMBER, K, OF EIGENVALUES
    NOT CALCULATED AND ONLY THE LAST N-K ELEMENTS OF THE ARRAYS VALR
    AND VALI ARE APPROXIMATE EIGENVALUES OF THE ORIGINAL MATRIX AND NO
    USEFUL EIGENVECTORS ARE DELIVERED.

PROCEDURES USED:

    EQILBRCOM = CP34361,
    COMEUCNRM = CP34359,
    HSHCOMHES = CP34366,
    QRICOM    = CP34373,
    BAKCOMHES = CP34367,
    BAKLBRCOM = CP34362,
    SCLCOM    = CP34360.

REQUIRED CENTRAL MEMORY: FIVE  REAL ARRAYS OF ORDER N AND ONE INTEGER
    ARRAY OF ORDER N ARE DECLARED.

RUNNING TIME: PROPORTIONAL TO N**2 * MAX(N, NUMBER OF ITERATIONS).

LANGUAGE : ALGOL 60.

THE FOLLOWING HOLDS FOR BOTH PROCEDURES:

METHOD AND PERFORMANCE:

    FOR CALCULATING THE EIGENVALUES AND EIGENVECTORS OF A COMPLEX
    MATRIX WE DISTINGUISH THE FOLLOWING STEPS:
    1) THE MATRIX IS EQUILIBRATED (SEE ALSO SECTION 3.2.1.1.2.).
    2) THE EQUILIBRATED MATRIX IS TRANSFORMED INTO HESSENBERG FORM BY
       MEANS OF HOUSEHOLDER MATRICES (SEE ALSO SECTION 3.2.1.2.2.2.).
    3) THE HESSENBERG MATRIX IS TRANSFORMED INTO AN UPPER TRIANGULAR
       MATRIX BY MEANS OF QR-ITERATION WITH SHIFT OF ORIGIN AND
       DEFLATION (SEE ALSO SECTION 3.3.2.2.1.).
    THE DIAGONAL ELEMENTS OF THE UPPER TRIANGULAR MATRIX ARE THE
    EIGENVALUES OF THE ORIGINAL MATRIX.
    THE EIGENVECTORS OF THE ORIGINAL MATRIX ARE OBTAINED BY CALCULATING
    THE EIGENVECTORS OF THE UPPER TRIANGULAR MATRIX (3) FOLLOWED BY
    BACKTRANSFORMATIONS CORRESPONDING TO (2) AND (1).

REFERENCES:

    DEKKER, T.J. (1968),
    ALGOL 60 PROCEDURES IN NUMERICAL ALGEBRA, PART 1,
    MATH. CENTRE TRACTS 22, MATHEMATISCH CENTRUM;

    DEKKER, T.J. AND W.HOFFMANN (1968),
    ALGOL 60 PROCEDURES IN NUMERICAL ALGEBRA, PART 2,
    MATH. CENTRE TRACTS 23, MATHEMATISCH CENTRUM;

    FRANCIS, J.G.F. (1961),
    THE QR-TRANSFORMATION, PART 1 AND 2,
    COMP.J. 4, P.265-271 AND P.332-345;

    MUELLER, D.J. (1966),
    HOUSEHOLDER,S METHOD  FOR  COMPLEX MATRICES AND EIGENSYSTEMS   OF
    HERMITIAN MATRICES,
    NUMER.MATH., 8, P.72-92;

    OSBORNE, E.E. (1960),
    ON PRECONDITIONING OF MATRICES,
    JACM., 7, P.338-354;

    PARLETT, B.N. AND C.REINSCH (1969),
    BALANCING A MATRIX FOR CALCULATION OF EIGENVALUES AND
    EIGENVECTORS,
    NUM. MATH., 13, P.293-304;

    WILKINSON, J.H. (1965),
    THE ALGEBRAIC EIGENVALUE PROBLEM,
    CLARENDOM PRESS, OXFORD.

EXAMPLE OF USE:

    EIGCOM CALCULATES THE EIGENVALUES AND THE EIGENVECTORS OF THE
    FOLLOWING MATRIX:
    (SEE WILKINSON AND REINSCH, 1971, CONTRIBUTION II/15)
    1+3*I  2+1*I  3+2*I  1+1*I
    3+4*I  1+2*I  2+1*I  4+3*I
    2+3*I  1+5*I  3+1*I  5+2*I
    1+2*I  3+1*I  1+4*I  5+3*I

    ONLY THE EIGENVECTOR CORRESPONDING TO VALR[1] + VALI[1] * I IS
    PRINTED BY THE FOLLOWING PROGRAM.

    "BEGIN"
    "REAL" "ARRAY" AR,AI,VR,VI[1:4,1:4],EM[0:7],VALR,VALI[1:4];
    "INTEGER" I;
    AR[1,1]:=AR[1,4]:=AR[2,2]:=AR[3,2]:=AR[4,1]:=AR[4,3]:=
    AI[1,2]:=AI[1,4]:=AI[2,3]:=AI[3,3]:=AI[4,2]:=1;
    AR[1,2]:=AR[2,3]:=AR[3,1]:=AI[1,3]:=AI[2,2]:=AI[3,4]:=AI[4,1]:=2;
    AR[1,3]:=AR[2,1]:=AR[3,3]:=AR[4,2]:=
    AI[1,1]:=AI[2,4]:=AI[3,1]:=AI[4,4]:=3;
    AR[2,4]:=AI[2,1]:=AI[4,3]:=4;
    AR[3,4]:=AR[4,4]:=AI[3,2]:=5;
    EM[0]:=5"-14;EM[2]:="-12;EM[4]:=10;EM[6]:=10;
    OUTPUT(61,"(""("EIGCOM: ")",D")",
              EIGCOM(AR,AI,4,EM,VALR,VALI,VR,VI));
    OUTPUT(61,"("/,"("EIGENVALUES:")",/")");
    "FOR" I:=1,2,3,4 "DO" OUTPUT(61,"("2(+D.4D),"(" * I")",/")",
                            VALR[I],VALI[I]);
    OUTPUT(61,"(""("FIRST EIGENVECTOR:")",/")");
    "FOR" I:=1,2,3,4 "DO" OUTPUT(61,"("2(+D.4D),"(" * I")",/")",
                            VR[I,1],VI[I,1]);
    OUTPUT(61,"(""("EM[1]: ")",+DD.DD/,"("EM[3]: ")",+D.D"+DD/,
      "("EM[5]: ")",+ZD/,"("EM[7]: ")",+ZD")",EM[1],EM[3],EM[5],EM[7]);
    "END"

    OUTPUT:

    EIGCOM: 0
    EIGENVALUES:
    -3.3710-0.7705 * I
    +9.7837+9.3225 * I
    +1.3657-1.4011 * I
    +2.2217+1.8490 * I
    FIRST EIGENVECTOR:
    -0.5061+0.5835 * I
    +1.0000+0.0000 * I
    +0.5183-0.7147 * I
    -0.5535+0.0188 * I
    EM[1]: +15.30
    EM[3]: +6.0"-12
    EM[5]:  +7
    EM[7]:  +4

SOURCE TEXT(S) :
"CODE" 34374;
"CODE" 34375;