LAPACK 3.3.0
|
00001 SUBROUTINE DGEHD2( N, ILO, IHI, A, LDA, TAU, WORK, INFO ) 00002 * 00003 * -- LAPACK routine (version 3.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 INTEGER IHI, ILO, INFO, LDA, N 00010 * .. 00011 * .. Array Arguments .. 00012 DOUBLE PRECISION A( LDA, * ), TAU( * ), WORK( * ) 00013 * .. 00014 * 00015 * Purpose 00016 * ======= 00017 * 00018 * DGEHD2 reduces a real general matrix A to upper Hessenberg form H by 00019 * an orthogonal similarity transformation: Q' * A * Q = H . 00020 * 00021 * Arguments 00022 * ========= 00023 * 00024 * N (input) INTEGER 00025 * The order of the matrix A. N >= 0. 00026 * 00027 * ILO (input) INTEGER 00028 * IHI (input) INTEGER 00029 * It is assumed that A is already upper triangular in rows 00030 * and columns 1:ILO-1 and IHI+1:N. ILO and IHI are normally 00031 * set by a previous call to DGEBAL; otherwise they should be 00032 * set to 1 and N respectively. See Further Details. 00033 * 1 <= ILO <= IHI <= max(1,N). 00034 * 00035 * A (input/output) DOUBLE PRECISION array, dimension (LDA,N) 00036 * On entry, the n by n general matrix to be reduced. 00037 * On exit, the upper triangle and the first subdiagonal of A 00038 * are overwritten with the upper Hessenberg matrix H, and the 00039 * elements below the first subdiagonal, with the array TAU, 00040 * represent the orthogonal matrix Q as a product of elementary 00041 * reflectors. See Further Details. 00042 * 00043 * LDA (input) INTEGER 00044 * The leading dimension of the array A. LDA >= max(1,N). 00045 * 00046 * TAU (output) DOUBLE PRECISION array, dimension (N-1) 00047 * The scalar factors of the elementary reflectors (see Further 00048 * Details). 00049 * 00050 * WORK (workspace) DOUBLE PRECISION array, dimension (N) 00051 * 00052 * INFO (output) INTEGER 00053 * = 0: successful exit. 00054 * < 0: if INFO = -i, the i-th argument had an illegal value. 00055 * 00056 * Further Details 00057 * =============== 00058 * 00059 * The matrix Q is represented as a product of (ihi-ilo) elementary 00060 * reflectors 00061 * 00062 * Q = H(ilo) H(ilo+1) . . . H(ihi-1). 00063 * 00064 * Each H(i) has the form 00065 * 00066 * H(i) = I - tau * v * v' 00067 * 00068 * where tau is a real scalar, and v is a real vector with 00069 * v(1:i) = 0, v(i+1) = 1 and v(ihi+1:n) = 0; v(i+2:ihi) is stored on 00070 * exit in A(i+2:ihi,i), and tau in TAU(i). 00071 * 00072 * The contents of A are illustrated by the following example, with 00073 * n = 7, ilo = 2 and ihi = 6: 00074 * 00075 * on entry, on exit, 00076 * 00077 * ( a a a a a a a ) ( a a h h h h a ) 00078 * ( a a a a a a ) ( a h h h h a ) 00079 * ( a a a a a a ) ( h h h h h h ) 00080 * ( a a a a a a ) ( v2 h h h h h ) 00081 * ( a a a a a a ) ( v2 v3 h h h h ) 00082 * ( a a a a a a ) ( v2 v3 v4 h h h ) 00083 * ( a ) ( a ) 00084 * 00085 * where a denotes an element of the original matrix A, h denotes a 00086 * modified element of the upper Hessenberg matrix H, and vi denotes an 00087 * element of the vector defining H(i). 00088 * 00089 * ===================================================================== 00090 * 00091 * .. Parameters .. 00092 DOUBLE PRECISION ONE 00093 PARAMETER ( ONE = 1.0D+0 ) 00094 * .. 00095 * .. Local Scalars .. 00096 INTEGER I 00097 DOUBLE PRECISION AII 00098 * .. 00099 * .. External Subroutines .. 00100 EXTERNAL DLARF, DLARFG, XERBLA 00101 * .. 00102 * .. Intrinsic Functions .. 00103 INTRINSIC MAX, MIN 00104 * .. 00105 * .. Executable Statements .. 00106 * 00107 * Test the input parameters 00108 * 00109 INFO = 0 00110 IF( N.LT.0 ) THEN 00111 INFO = -1 00112 ELSE IF( ILO.LT.1 .OR. ILO.GT.MAX( 1, N ) ) THEN 00113 INFO = -2 00114 ELSE IF( IHI.LT.MIN( ILO, N ) .OR. IHI.GT.N ) THEN 00115 INFO = -3 00116 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 00117 INFO = -5 00118 END IF 00119 IF( INFO.NE.0 ) THEN 00120 CALL XERBLA( 'DGEHD2', -INFO ) 00121 RETURN 00122 END IF 00123 * 00124 DO 10 I = ILO, IHI - 1 00125 * 00126 * Compute elementary reflector H(i) to annihilate A(i+2:ihi,i) 00127 * 00128 CALL DLARFG( IHI-I, A( I+1, I ), A( MIN( I+2, N ), I ), 1, 00129 $ TAU( I ) ) 00130 AII = A( I+1, I ) 00131 A( I+1, I ) = ONE 00132 * 00133 * Apply H(i) to A(1:ihi,i+1:ihi) from the right 00134 * 00135 CALL DLARF( 'Right', IHI, IHI-I, A( I+1, I ), 1, TAU( I ), 00136 $ A( 1, I+1 ), LDA, WORK ) 00137 * 00138 * Apply H(i) to A(i+1:ihi,i+1:n) from the left 00139 * 00140 CALL DLARF( 'Left', IHI-I, N-I, A( I+1, I ), 1, TAU( I ), 00141 $ A( I+1, I+1 ), LDA, WORK ) 00142 * 00143 A( I+1, I ) = AII 00144 10 CONTINUE 00145 * 00146 RETURN 00147 * 00148 * End of DGEHD2 00149 * 00150 END