MB03YA

Annihilation of one or two entries on the subdiagonal of a Hessenberg matrix corresponding to zero elements on the diagonal of a triangular matrix

[Specification] [Arguments] [Method] [References] [Comments] [Example]

Purpose

  To annihilate one or two entries on the subdiagonal of the
  Hessenberg matrix A for dealing with zero elements on the diagonal
  of the triangular matrix B.

  MB03YA is an auxiliary routine called by SLICOT Library routines
  MB03XP and MB03YD.

Specification
      SUBROUTINE MB03YA( WANTT, WANTQ, WANTZ, N, ILO, IHI, ILOQ, IHIQ,
     $                   POS, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO )
C     .. Scalar Arguments ..
      LOGICAL            WANTQ, WANTT, WANTZ
      INTEGER            IHI, IHIQ, ILO, ILOQ, INFO, LDA, LDB, LDQ, LDZ,
     $                   N, POS
C     .. Array Arguments ..
      DOUBLE PRECISION   A(LDA,*), B(LDB,*), Q(LDQ,*), Z(LDZ,*)

Arguments

Mode Parameters

  WANTT   LOGICAL
          Indicates whether the user wishes to compute the full
          Schur form or the eigenvalues only, as follows:
          = .TRUE. :  Compute the full Schur form;
          = .FALSE.:  compute the eigenvalues only.

  WANTQ   LOGICAL
          Indicates whether or not the user wishes to accumulate
          the matrix Q as follows:
          = .TRUE. :  The matrix Q is updated;
          = .FALSE.:  the matrix Q is not required.

  WANTZ   LOGICAL
          Indicates whether or not the user wishes to accumulate
          the matrix Z as follows:
          = .TRUE. :  The matrix Z is updated;
          = .FALSE.:  the matrix Z is not required.

Input/Output Parameters
  N       (input) INTEGER
          The order of the matrices A and B. N >= 0.

  ILO     (input) INTEGER
  IHI     (input) INTEGER
          It is assumed that the matrices A and B are already
          (quasi) upper triangular in rows and columns 1:ILO-1 and
          IHI+1:N. The routine works primarily with the submatrices
          in rows and columns ILO to IHI, but applies the
          transformations to all the rows and columns of the
          matrices A and B, if WANTT = .TRUE..
          1 <= ILO <= max(1,N); min(ILO,N) <= IHI <= N.

  ILOQ    (input) INTEGER
  IHIQ    (input) INTEGER
          Specify the rows of Q and Z to which transformations
          must be applied if WANTQ = .TRUE. and WANTZ = .TRUE.,
          respectively.
          1 <= ILOQ <= ILO; IHI <= IHIQ <= N.

  POS     (input) INTEGER
          The position of the zero element on the diagonal of B.
          ILO <= POS <= IHI.

  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the leading N-by-N part of this array must
          contain the upper Hessenberg matrix A.
          On exit, the leading N-by-N part of this array contains
          the updated matrix A where A(POS,POS-1) = 0, if POS > ILO,
          and A(POS+1,POS) = 0, if POS < IHI.

  LDA     INTEGER
          The leading dimension of the array A.  LDA >= MAX(1,N).

  B       (input/output) DOUBLE PRECISION array, dimension (LDB,N)
          On entry, the leading N-by-N part of this array must
          contain an upper triangular matrix B with B(POS,POS) = 0.
          On exit, the leading N-by-N part of this array contains
          the updated upper triangular matrix B.

  LDB     INTEGER
          The leading dimension of the array B.  LDB >= MAX(1,N).

  Q       (input/output) DOUBLE PRECISION array, dimension (LDQ,N)
          On entry, if WANTQ = .TRUE., then the leading N-by-N part
          of this array must contain the current matrix Q of
          transformations accumulated by MB03XP.
          On exit, if WANTQ = .TRUE., then the leading N-by-N part
          of this array contains the matrix Q updated in the
          submatrix Q(ILOQ:IHIQ,ILO:IHI).
          If WANTQ = .FALSE., Q is not referenced.

  LDQ     INTEGER
          The leading dimension of the array Q.  LDQ >= 1.
          If WANTQ = .TRUE., LDQ >= MAX(1,N).

  Z       (input/output) DOUBLE PRECISION array, dimension (LDZ,N)
          On entry, if WANTZ = .TRUE., then the leading N-by-N part
          of this array must contain the current matrix Z of
          transformations accumulated by MB03XP.
          On exit, if WANTZ = .TRUE., then the leading N-by-N part
          of this array contains the matrix Z updated in the
          submatrix Z(ILOQ:IHIQ,ILO:IHI).
          If WANTZ = .FALSE., Z is not referenced.

  LDZ     INTEGER
          The leading dimension of the array Z.  LDZ >= 1.
          If WANTZ = .TRUE., LDZ >= MAX(1,N).

Error Indicator
  INFO    INTEGER
          = 0:  successful exit;
          < 0:  if INFO = -i, the i-th argument had an illegal
                value.

Method
  The method is illustrated by Wilkinson diagrams for N = 5,
  POS = 3:

        [ x x x x x ]       [ x x x x x ]
        [ x x x x x ]       [ o x x x x ]
    A = [ o x x x x ],  B = [ o o o x x ].
        [ o o x x x ]       [ o o o x x ]
        [ o o o x x ]       [ o o o o x ]

  First, a QR factorization is applied to A(1:3,1:3) and the
  resulting nonzero in the updated matrix B is immediately
  annihilated by a Givens rotation acting on columns 1 and 2:

        [ x x x x x ]       [ x x x x x ]
        [ x x x x x ]       [ o x x x x ]
    A = [ o o x x x ],  B = [ o o o x x ].
        [ o o x x x ]       [ o o o x x ]
        [ o o o x x ]       [ o o o o x ]

  Secondly, an RQ factorization is applied to A(4:5,4:5) and the
  resulting nonzero in the updated matrix B is immediately
  annihilated by a Givens rotation acting on rows 4 and 5:

        [ x x x x x ]       [ x x x x x ]
        [ x x x x x ]       [ o x x x x ]
    A = [ o o x x x ],  B = [ o o o x x ].
        [ o o o x x ]       [ o o o x x ]
        [ o o o x x ]       [ o o o o x ]

References
  [1] Bojanczyk, A.W., Golub, G.H., and Van Dooren, P.
      The periodic Schur decomposition: Algorithms and applications.
      Proc. of the SPIE Conference (F.T. Luk, Ed.), 1770, pp. 31-42,
      1992.

Numerical Aspects
  The algorithm requires O(N**2) floating point operations and is
  backward stable.

Further Comments
  None
Example

Program Text

  None
Program Data
  None
Program Results
  None

Return to Supporting Routines index