SUBROUTINE ZUNGHR( N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO ) * * -- LAPACK routine (version 3.1) -- * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. * November 2006 * * .. Scalar Arguments .. INTEGER IHI, ILO, INFO, LDA, LWORK, N * .. * .. Array Arguments .. COMPLEX*16 A( LDA, * ), TAU( * ), WORK( * ) * .. * * Purpose * ======= * * ZUNGHR generates a complex unitary matrix Q which is defined as the * product of IHI-ILO elementary reflectors of order N, as returned by * ZGEHRD: * * Q = H(ilo) H(ilo+1) . . . H(ihi-1). * * Arguments * ========= * * N (input) INTEGER * The order of the matrix Q. N >= 0. * * ILO (input) INTEGER * IHI (input) INTEGER * ILO and IHI must have the same values as in the previous call * of ZGEHRD. Q is equal to the unit matrix except in the * submatrix Q(ilo+1:ihi,ilo+1:ihi). * 1 <= ILO <= IHI <= N, if N > 0; ILO=1 and IHI=0, if N=0. * * A (input/output) COMPLEX*16 array, dimension (LDA,N) * On entry, the vectors which define the elementary reflectors, * as returned by ZGEHRD. * On exit, the N-by-N unitary matrix Q. * * LDA (input) INTEGER * The leading dimension of the array A. LDA >= max(1,N). * * TAU (input) COMPLEX*16 array, dimension (N-1) * TAU(i) must contain the scalar factor of the elementary * reflector H(i), as returned by ZGEHRD. * * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK)) * On exit, if INFO = 0, WORK(1) returns the optimal LWORK. * * LWORK (input) INTEGER * The dimension of the array WORK. LWORK >= IHI-ILO. * For optimum performance LWORK >= (IHI-ILO)*NB, where NB is * the optimal blocksize. * * If LWORK = -1, then a workspace query is assumed; the routine * only calculates the optimal size of the WORK array, returns * this value as the first entry of the WORK array, and no error * message related to LWORK is issued by XERBLA. * * INFO (output) INTEGER * = 0: successful exit * < 0: if INFO = -i, the i-th argument had an illegal value * * ===================================================================== * * .. Parameters .. COMPLEX*16 ZERO, ONE PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ), $ ONE = ( 1.0D+0, 0.0D+0 ) ) * .. * .. Local Scalars .. LOGICAL LQUERY INTEGER I, IINFO, J, LWKOPT, NB, NH * .. * .. External Subroutines .. EXTERNAL XERBLA, ZUNGQR * .. * .. External Functions .. INTEGER ILAENV EXTERNAL ILAENV * .. * .. Intrinsic Functions .. INTRINSIC MAX, MIN * .. * .. Executable Statements .. * * Test the input arguments * INFO = 0 NH = IHI - ILO LQUERY = ( LWORK.EQ.-1 ) IF( N.LT.0 ) THEN INFO = -1 ELSE IF( ILO.LT.1 .OR. ILO.GT.MAX( 1, N ) ) THEN INFO = -2 ELSE IF( IHI.LT.MIN( ILO, N ) .OR. IHI.GT.N ) THEN INFO = -3 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN INFO = -5 ELSE IF( LWORK.LT.MAX( 1, NH ) .AND. .NOT.LQUERY ) THEN INFO = -8 END IF * IF( INFO.EQ.0 ) THEN NB = ILAENV( 1, 'ZUNGQR', ' ', NH, NH, NH, -1 ) LWKOPT = MAX( 1, NH )*NB WORK( 1 ) = LWKOPT END IF * IF( INFO.NE.0 ) THEN CALL XERBLA( 'ZUNGHR', -INFO ) RETURN ELSE IF( LQUERY ) THEN RETURN END IF * * Quick return if possible * IF( N.EQ.0 ) THEN WORK( 1 ) = 1 RETURN END IF * * Shift the vectors which define the elementary reflectors one * column to the right, and set the first ilo and the last n-ihi * rows and columns to those of the unit matrix * DO 40 J = IHI, ILO + 1, -1 DO 10 I = 1, J - 1 A( I, J ) = ZERO 10 CONTINUE DO 20 I = J + 1, IHI A( I, J ) = A( I, J-1 ) 20 CONTINUE DO 30 I = IHI + 1, N A( I, J ) = ZERO 30 CONTINUE 40 CONTINUE DO 60 J = 1, ILO DO 50 I = 1, N A( I, J ) = ZERO 50 CONTINUE A( J, J ) = ONE 60 CONTINUE DO 80 J = IHI + 1, N DO 70 I = 1, N A( I, J ) = ZERO 70 CONTINUE A( J, J ) = ONE 80 CONTINUE * IF( NH.GT.0 ) THEN * * Generate Q(ilo+1:ihi,ilo+1:ihi) * CALL ZUNGQR( NH, NH, NH, A( ILO+1, ILO+1 ), LDA, TAU( ILO ), $ WORK, LWORK, IINFO ) END IF WORK( 1 ) = LWKOPT RETURN * * End of ZUNGHR * END