NUMAL Section 3.3.2.1

BEGIN SECTION : 3.3.2.1 (July, 1974)

AUTHOR   : C.G. VAN DER LAAN.

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

INSTITUTE : MATHEMATICAL CENTRE.

RECEIVED: 730917.

BRIEF DESCRIPTION:

    THIS SECTION CONTAINS FOUR PROCEDURES FOR CALCULATING THE
    EIGENVALUES OR THE EIGENVALUES AND EIGENVECTORS OF COMPLEX
    HERMITIAN MATRICES.
    EIGVALHRM CALCULATES THE EIGENVALUES OF A HERMITIAN MATRIX.
    EIGHRM CALCULATES THE EIGENVALUES AND EIGENVECTORS OF A HERMITIAN
    MATRIX.
    QRIVALHRM CALCULATES THE EIGENVALUES OF A HERMITIAN MATRIX.
    QRIHRM CALCULATES THE EIGENVALUES AND EIGENVECTORS OF A HERMITIAN
    MATRIX.
    WHEN A SMALL NUMBER OF EIGENVALUES OR EIGENVALUES AND EIGENVECTORS
    IS REQUIRED, THE USE OF EIGVALHRM OR EIGHRM IS RECOMMANDED; WHEN
    MORE THAN, SAY, 25 PERCENT OF THE EIGENSYSTEM IS REQUIRED. THE
    PROCEDURES QRIVALHRM OR QRIHRM ARE TO BE USED.


SUBSECTION: EIGVALHRM.

CALLING SEQUENCE :

    THE HEADING OF THE PROCEDURE READS :
    "PROCEDURE" EIGVALHRM(A, N, NUMVAL, VAL, EM); "VALUE" N, NUMVAL;
    "INTEGER" N, NUMVAL; "ARRAY" A, VAL, EM;
    "CODE" 34368;

    THE MEANING OF THE FORMAL PARAMETERS IS :
    A:         <ARRAY IDENTIFIER>;
               "ARRAY" A[1:N,1:N];
               ENTRY: THE  REAL  PART  OF  THE  UPPER  TRIANGLE OF  THE
                      HERMITIAN   MATRIX  MUST  BE  GIVEN IN THE  UPPER
                      TRIANGULAR PART OF A (THE ELEMENTS A[I,J], I<=J);
                      THE  IMAGINARY  PART OF THE STRICT LOWER TRIANGLE
                      OF  THE   HERMITIAN  MATRIX  MUST BE GIVEN IN THE
                      STRICT LOWER PART OF A (THE ELEMENTS A[I,J],I>J);
               THE ELEMENTS AF A ARE ALTERED;
    N:         <ARITHMETIC EXPRESSION>;
               THE ORDER OF THE GIVEN MATRIX;
    NUMVAL:    <ARITHMETIC EXPRESSION>;
               EIGVALHRM  CALCULATES  THE LARGEST NUMVAL EIGENVALUES OF
               THE HERMITIAN MATRIX;
    VAL:       <ARRAY IDENTIFIER>;
               "ARRAY" VAL[1:NUMVAL];
               EXIT:
               IN ARRAY VAL THE LARGEST NUMVAL EIGENVALUES ARE
               DELIVERED IN MONOTONICALLY NONINCREASING ORDER;
    EM:        <ARRAY IDENTIFIER>;
               "ARRAY" EM[0:3];
               ENTRY:
               EM[0]: THE MACHINE PRECISION;
               EM[2]: THE RELATIVE TOLERANCE FOR THE EIGENVALUES;
                      MORE PRECISELY: THE TOLERANCE FOR EACH EIGENVALUE
                      LAMBDA, IS ABS(LAMBDA)*EM[2]+EM[1]*EM[0];
               EXIT:
               EM[1]: AN ESTIMATE OF A NORM OF THE ORIGINAL MATRIX;
               EM[3]: THE NUMBER OF ITERATIONS PERFORMED.

PROCEDURES USED:

    HSHHRMTRIVAL = CP34364,
    VALSYMTRI    = CP34151.

REQUIRED CENTRAL MEMORY:
    TWO AUXILIARY ARRAYS OF ORDER N AND N - 1 RESPECTIVELY ARE DECLARED

RUNNING TIME: PROPORTIONAL TO N CUBED.

LANGUAGE: ALGOL 60.

METHOD AND PERFORMANCE: SEE QRIHRM (THIS SECTION).

EXAMPLE OF USE:  SEE EIGHRM (THIS SECTION).


SUBSECTION: EIGHRM.

CALLING SEQUENCE :

    THE HEADING OF THE PROCEDURE READS :
    "PROCEDURE" EIGHRM(A, N, NUMVAL, VAL, VECR, VECI, EM);
    "VALUE" N, NUMVAL; "INTEGER" N, NUMVAL;
    "ARRAY" A, VAL, VECR, VECI, EM;
    "CODE" 34369;

    THE MEANING OF THE FORMAL PARAMETERS IS :
    A:         <ARRAY IDENTIFIER>;
               "ARRAY" A[1:N,1:N];
               ENTRY: THE  REAL  PART  OF  THE  UPPER  TRIANGLE OF  THE
                      HERMITIAN   MATRIX  MUST  BE  GIVEN IN THE  UPPER
                      TRIANGULAR PART OF A (THE ELEMENTS A[I,J], I<=J);
                      THE  IMAGINARY  PART OF THE STRICT LOWER TRIANGLE
                      OF  THE   HERMITIAN  MATRIX  MUST BE GIVEN IN THE
                      STRICT LOWER PART OF A (THE ELEMENTS A[I,J],I>J);
               THE ELEMENTS AF A ARE ALTERED;
    N:         <ARITHMETIC EXPRESSION>;
               THE ORDER OF THE GIVEN MATRIX;
    NUMVAL:    <ARITHMETIC EXPRESSION>;
               EIGHRM  CALCULATES THE LARGEST NUMVAL EIGENVALUES OF THE
               HERMITIAN MATRIX;
    VAL:       <ARRAY IDENTIFIER>;
               "ARRAY" VAL[1:NUMVAL];
               EXIT:
               IN ARRAY VAL THE LARGEST NUMVAL EIGENVALUES ARE
               DELIVERED IN MONOTONICALLY NONINCREASING ORDER;
    VECR,VECI: <ARRAY IDENTIFIER>;
               "ARRAY" VECR,VECI[1:N,1:NUMVAL];
               EXIT:
               THE CALCULATED EIGENVECTORS;
               THE  COMPLEX EIGENVECTOR WITH REAL PART VECR[1:N,I]  AND
               IMAGINARY PART VECI[1:N,I] CORRESPONDS TO THE EIGENVALUE
               VAL[I], I=1,...,NUMVAL;

    EM:        <ARRAY IDENTIFIER>;
               "ARRAY" EM[0:9];
               ENTRY:
               EM[0]: THE MACHINE PRECISION;
               EM[2]: THE RELATIVE TOLERANCE FOR THE EIGENVALUES;
                      MORE PRECISELY: THE TOLERANCE FOR EACH EIGENVALUE
                      LAMBDA, IS ABS(LAMBDA)*EM[2]+EM[1]*EM[0];
               EM[4]: THE ORTHOGONALIZATION PARAMETER (E.G. .01);
               EM[6]: THE TOLERANCE FOR THE EIGENVECTORS;
               EM[8]: THE  MAXIMUM NUMBER OF INVERSE ITERATIONS ALLOWED
                      FOR THE CALCULATION OF EACH EIGENVECTOR;
               EXIT:
               EM[1]: AN ESTIMATE OF A NORM OF THE ORIGINAL MATRIX;
               EM[3]: THE NUMBER OF ITERATIONS PERFORMED;
               EM[5]: THE  NUMBER  OF EIGENVECTORS INVOLVED IN THE LAST
                      GRAM-SCHMIDT ORTHOGONALIZATION;
               EM[7]: THE MAXIMUM EUCLIDEAN NORM OF THE RESIDUES OF THE
                      CALCULATED EIGENVECTORS;
               EM[9]: THE   LARGEST  NUMBER   OF  INVERSE    ITERATIONS
                      PERFORMED   FOR   THE   CALCULATION   OF     SOME
                      EIGENVECTOR; IF, HOWEVER,  FOR  SOME   CALCULATED
                      EIGENVECTOR, THE  EUCLIDEAN NORM OF THE  RESIDUES
                      REMAINS     GREATER  THAN   EM[1]*EM[6],     THEN
                      EM[9]:=EM[8]+1.

PROCEDURES USED:

    HSHHRMTRI = CP34363,
    VALSYMTRI = CP34151,
    VECSYMTRI = CP34152,
    BAKHRMTRI = CP34365.

REQUIRED CENTRAL MEMORY:
    THREE AUXILIARY ARRAYS OF ORDER N-1 AND TWO OF ORDER N ARE DECLARED

RUNNING TIME: PROPORTIONAL TO N CUBED.

LANGUAGE: ALGOL 60.

METHOD AND PERFORMANCE: SEE QRIHRM (THIS SECTION).

EXAMPLE OF USE:

    LET  EIGHRM  CALCULATE THE LARGEST EIGENVALUE AND THE CORRESPONDING
    EIGENVECTOR OF THE FOLLOWING MATRIX:
    (SEE GREGORY AND KARNEY, CHAPTER 6, EXAMPLE 6.6)
         3    1    0   +2I
         1    3   -2I   0
         0   +2I   1    1
        -2I   0    1    1
    THE EIGENVECTORS ARE NORMALIZED BY THE PROCEDURE SCLCOM (SEE
    SECTION 1.2.11. ).

    "BEGIN"
    "COMMENT" GREGORY AND KARNEY,CHAPTER 6, EXAMPLE 6.6;
    "REAL" "ARRAY" A[1:4,1:4],VAL[1:1],VECR,VECI[1:4,1:1],EM[0:9];
    "INTEGER" I;
    INIMAT(1,4,1,4,A,0);
    A[1,1]:=A[2,2]:=3;
    A[1,2]:=A[3,3]:=A[3,4]:=A[4,4]:=1;
    A[3,2]:=2;A[4,1]:=-2;
    EM[0]:=5"-14;EM[2]:="-12;
    EM[4]:=.01;EM[6]:="-10;EM[8]:=5;
    EIGHRM(A,4,1,VAL,VECR,VECI,EM);
    SCLCOM(VECR,VECI,4,1,1);
    OUTPUT(61,"(""("LARGEST EIGENVALUE: ")",N/")",VAL[1]);
    OUTPUT(61,"(""("CORRESPONDING EIGENVECTOR:")",/")");
    "FOR" I:=1,2,3,4 "DO"
    OUTPUT(61,"("+D.D,+D.DDD,"("*I")",/")",VECR[I,1],VECI[I,1]);
    "FOR" I:=1,3,5,7,9 "DO"
    OUTPUT(61,"("/,"("EM[")",D,"("]: ")",+D.DDD"+DD")",I,EM[I]);
    "END"

    DELIVERS:

    LARGEST EIGENVALUE: +4.8284271247462"000
    CORRESPONDING EIGENVECTOR:
    +1.0+0.000*I
    +1.0+0.000*I
    -0.0+0.414*I
    +0.0-0.414*I

    EM[1]: +6.000"+00
    EM[3]: +1.800"+01
    EM[5]: +1.000"+00
    EM[7]: +5.303"-14
    EM[9]: +1.000"+00


SUBSECTION: QRIVALHRM.

CALLING SEQUENCE :

    THE HEADING OF THE PROCEDURE READS :
    "INTEGER" "PROCEDURE" QRIVALHRM(A, N, VAL, EM); "VALUE" N;
    "INTEGER" N; "ARRAY" A, VAL, EM;
    "CODE" 34370;

    THE MEANING OF THE FORMAL PARAMETERS IS:
    A:         <ARRAY IDENTIFIER>;
               "ARRAY" A[1:N,1:N];
               ENTRY: THE  REAL  PART  OF  THE  UPPER  TRIANGLE OF  THE
                      HERMITIAN   MATRIX  MUST  BE  GIVEN IN THE  UPPER
                      TRIANGULAR PART OF A (THE ELEMENTS A[I,J], I<=J);
                      THE  IMAGINARY  PART OF THE STRICT LOWER TRIANGLE
                      OF  THE   HERMITIAN  MATRIX  MUST BE GIVEN IN THE
                      STRICT LOWER PART OF A (THE ELEMENTS A[I,J],I>J);
               THE ELEMENTS AF A ARE ALTERED;
    N:         <ARITHMETIC EXPRESSION>;
               THE ORDER OF THE GIVEN MATRIX;
    VAL:       <ARRAY IDENTIFIER>;
               "ARRAY" VAL[1:N];
               EXIT:
               THE CALCULATED EIGENVALUES;
    EM:        <ARRAY IDENTIFIER>;
               "ARRAY" EM[0:5];
               ENTRY:
               EM[0]: THE MACHINE PRECISION;
               EM[2]: THE RELATIVE TOLERANCE FOR THE QR ITERATION;
               EM[4]: THE MAXIMUM ALLOWED NUMBER OF ITERATIONS;
               EXIT:
               EM[1]: AN ESTIMATE OF A NORM OF THE ORIGINAL MATRIX;
               EM[3]: THE  MAXIMUM   ABSOLUTE VALUE OF THE   CODIAGONAL
                      ELEMENTS NEGLECTED;
               EM[5]: THE NUMBER OF ITERATIONS PERFORMED;
                      EM[5]:= EM[4]+1 IN THE CASE QRIVALHRM^=0;

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

PROCEDURES USED:

    HSHHRMTRIVAL = CP34364,
    QRIVALSYMTRI = CP34160.

REQUIRED CENTRAL MEMORY:
    TWO  AUXILIARY ARRAYS OF ORDER N  ARE DECLARED.

RUNNING TIME: PROPORTIONAL TO N CUBED.

LANGUAGE: ALGOL 60.

METHOD AND PERFORMANCE: SEE QRIHRM (THIS SECTION).

EXAMPLE OF USE: SEE QRIHRM (THIS SECTION).


SUBSECTION: QRIHRM.

CALLING SEQUENCE :

    THE HEADING OF THE PROCEDURE READS :
    "INTEGER" "PROCEDURE" QRIHRM(A, N, VAL, VR, VI, EM); "VALUE" N;
    "INTEGER" N; "ARRAY" A, VAL, VR, VI, EM;
    "CODE" 34371;

    THE MEANING OF THE FORMAL PARAMETERS IS:
    A:         <ARRAY IDENTIFIER>;
               "ARRAY" A[1:N,1:N];
               ENTRY: THE  REAL  PART  OF  THE  UPPER  TRIANGLE OF  THE
                      HERMITIAN   MATRIX  MUST  BE  GIVEN IN THE  UPPER
                      TRIANGULAR PART OF A (THE ELEMENTS A[I,J], I<=J);
                      THE  IMAGINARY  PART OF THE STRICT LOWER TRIANGLE
                      OF  THE   HERMITIAN  MATRIX  MUST BE GIVEN IN THE
                      STRICT LOWER PART OF A (THE ELEMENTS A[I,J],I>J);
               THE ELEMENTS AF A ARE ALTERED;
    N:         <ARITHMETIC EXPRESSION>;
               THE ORDER OF THE GIVEN MATRIX;
    VAL:       <ARRAY IDENTIFIER>;
               "ARRAY" VAL[1:N];
               EXIT:
               THE CALCULATED EIGENVALUES;
    VR,VI:     <ARRAY IDENTIFIER>;
               "ARRAY" VR,VI[1:N,1:N];
               EXIT:
               THE CALCULATED EIGENVECTORS;
               THE  COMPLEX  EIGENVECTOR WITH REAL PART VR[1:N,I]   AND
               IMAGINARY  PART VI[1:N,I] CORRESPONDS TO THE  EIGENVALUE
               VAL[I], I=1,...,N;
    EM:        <ARRAY IDENTIFIER>;
               "ARRAY" EM[0:5];
               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 ITERATIONS;
                      (E.G. 10 * N);
               EXIT:
               EM[1]: AN ESTIMATE OF A NORM OF THE ORIGINAL MATRIX;
               EM[3]: THE  MAXIMUM   ABSOLUTE VALUE OF THE   CODIAGONAL
                      ELEMENTS NEGLECTED;
               EM[5]: THE NUMBER OF ITERATIONS PERFORMED;
                      EM[5]:=EM[4]+1 IN THE CASE QRIHRM^=0;

    QRIHRM:=0, PROVIDED  THE   PROCESS   IS   COMPLETED   WITHIN  EM[4]
    ITERATIONS; OTHERWISE, QRIHRM:= THE  NUMBER OF EIGENVALUES, K, NOT
    CALCULATED AND ONLY THE LAST N-K ELEMENTS OF  VAL ARE APPROXIMATE
    EIGENVALUES AND THE COLUMNS OF THE ARRAYS VR,VI[1:N,N-K:N] ARE
    APPROXIMATE EIGENVECTORS OF THE ORIGINAL HERMITEAN MATRIX .

PROCEDURES USED:

    HSHHRMTRI = CP34363,
    QRISYMTRI = CP34161,
    BAKHRMTRI = CP34365.

REQUIRED CENTRAL MEMORY:
    TWO AUXILIARY ARRAYS OF ORDER N - 1 AND TWO OF ORDER N ARE DECLARED

RUNNING TIME: PROPORTIONAL TO N CUBED.

LANGUAGE: ALGOL 60.

THE FOLLOWING HOLDS FOR THE FOUR PROCEDURES OF THIS SECTION:

METHOD AND PERFORMANCE:

    FOR THE TRANSFORMATION OF THE GIVEN HERMITIAN MATRIX INTO A REAL
    SYMMETRIC TRIDIAGONAL MATRIX, AND FOR THE CORRESPONDING BACK
    TRANSFORMATION, PROCEDURES OF SECTION  3.2.1.2.2.1. ARE USED.
    FOR THE CALCULATION OF THE EIGENVALUES AND EIGENVECTORS OF THE
    RESULTING SYMMETRIC TRIDIAGONAL MATRIX, PROCEDURES OF SECTION
    3.3.1.1.1. ARE USED.

EXAMPLE OF USE:

    QRIHRM CALCULATES THE EIGENVALUES AND EIGENVECTORS OF THE
    FOLLOWING MATRIX:
    (SEE GREGORY AND KARNEY, CHAPTER 6, EXAMPLE 6.6)
         3    1    0   +2I
         1    3   -2I   0
         0   +2I   1    1
        -2I   0    1    1
    THE EIGENVECTORS ARE NORMALIZED BY THE PROCEDURE SCLCOM (SEE
    SECTION 1.2.11. ).
    ONLY THE EIGENVECTORS CORRESPONDING TO VAL[2] AND VAL[3] ARE
    PRINTED BY THE FOLLOWING PROGRAM:

    "BEGIN"
    "COMMENT" GREGORY AND KARNEY,CHAPTER 6, EXAMPLE 6.6;
    "REAL" "ARRAY" A,VR,VI[1:4,1:4],VAL[1:4],EM[0:5];"INTEGER" I;
    INIMAT(1,4,1,4,A,0);
    A[1,1]:=A[2,2]:=3;
    A[3,2]:=2;A[4,1]:=-2;
    A[1,2]:=A[3,3]:=A[3,4]:=A[4,4]:=1;
    EM[0]:=EM[2]:=5"-14;EM[4]:=20;
    OUTPUT(61,"(""("QRIHRM: ")",D/")",QRIHRM(A,4,VAL,VR,VI,EM));
    SCLCOM(VR,VI,4,2,3);
    OUTPUT(61,"(""("EIGENVALUES: ")"")");
    "FOR" I:=1,2,3,4 "DO" OUTPUT(61,"("/,"("VAL[")",D,"("]: ")",
                                    +D.3DBB")",I,VAL[I]);
    OUTPUT(61,"("/,"("EIGENVECTORS CORRESPONDING TO")",/,
                 "(" VAL[2] , VAL[3] ")",/")");
    "FOR" I:=1,2,3,4 "DO"
    OUTPUT(61,"("+D,+D,"("*I  ,  ")",+D,+D,"("*I")",/")",
              VR[I,2],VI[I,2],VR[I,3],VI[I,3]);
    "FOR" I:=1,3,5 "DO"
     OUTPUT(61,"("/,"("EM[")",D,"("]:  ")",+D.DDD"+DD")",I,EM[I])
    "END"

  OUTPUT:
    QRIHRM: 0
    EIGENVALUES:
    VAL[1]: +4.828
    VAL[2]: +4.000
    VAL[3]: -0.000
    VAL[4]: -0.828
    EIGENVECTORS CORRESPONDING TO
     VAL[2] , VAL[3]
    +1+0*I  ,  +0-1*I
    -1+0*I  ,  +0+1*I
    +0-1*I  ,  +1+0*I
    +0-1*I  ,  +1+0*I

    EM[1]:  +6.000"+00
    EM[3]:  +3.804"-22
    EM[5]:  +6.000"+00

SOURCE TEXT(S) :
"CODE" 34368;
"CODE" 34369;
"CODE" 34370;
"CODE" 34371;