LAPACK 3.11.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ dgehrd()

 subroutine dgehrd ( integer N, integer ILO, integer IHI, double precision, dimension( lda, * ) A, integer LDA, double precision, dimension( * ) TAU, double precision, dimension( * ) WORK, integer LWORK, integer INFO )

DGEHRD

Purpose:
``` DGEHRD reduces a real general matrix A to upper Hessenberg form H by
an orthogonal similarity transformation:  Q**T * A * Q = H .```
Parameters
 [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in] ILO ` ILO is INTEGER` [in] IHI ``` IHI is INTEGER It is assumed that A is already upper triangular in rows and columns 1:ILO-1 and IHI+1:N. ILO and IHI are normally set by a previous call to DGEBAL; otherwise they should be set to 1 and N respectively. See Further Details. 1 <= ILO <= IHI <= N, if N > 0; ILO=1 and IHI=0, if N=0.``` [in,out] A ``` A is DOUBLE PRECISION array, dimension (LDA,N) On entry, the N-by-N general matrix to be reduced. On exit, the upper triangle and the first subdiagonal of A are overwritten with the upper Hessenberg matrix H, and the elements below the first subdiagonal, with the array TAU, represent the orthogonal matrix Q as a product of elementary reflectors. See Further Details.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,N).``` [out] TAU ``` TAU is DOUBLE PRECISION array, dimension (N-1) The scalar factors of the elementary reflectors (see Further Details). Elements 1:ILO-1 and IHI:N-1 of TAU are set to zero.``` [out] WORK ``` WORK is DOUBLE PRECISION array, dimension (LWORK) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The length of the array WORK. LWORK >= max(1,N). For good performance, LWORK should generally be larger. 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.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value.```
Further Details:
```  The matrix Q is represented as a product of (ihi-ilo) elementary
reflectors

Q = H(ilo) H(ilo+1) . . . H(ihi-1).

Each H(i) has the form

H(i) = I - tau * v * v**T

where tau is a real scalar, and v is a real vector with
v(1:i) = 0, v(i+1) = 1 and v(ihi+1:n) = 0; v(i+2:ihi) is stored on
exit in A(i+2:ihi,i), and tau in TAU(i).

The contents of A are illustrated by the following example, with
n = 7, ilo = 2 and ihi = 6:

on entry,                        on exit,

( a   a   a   a   a   a   a )    (  a   a   h   h   h   h   a )
(     a   a   a   a   a   a )    (      a   h   h   h   h   a )
(     a   a   a   a   a   a )    (      h   h   h   h   h   h )
(     a   a   a   a   a   a )    (      v2  h   h   h   h   h )
(     a   a   a   a   a   a )    (      v2  v3  h   h   h   h )
(     a   a   a   a   a   a )    (      v2  v3  v4  h   h   h )
(                         a )    (                          a )

where a denotes an element of the original matrix A, h denotes a
modified element of the upper Hessenberg matrix H, and vi denotes an
element of the vector defining H(i).

This file is a slight modification of LAPACK-3.0's DGEHRD
subroutine incorporating improvements proposed by Quintana-Orti and
Van de Geijn (2006). (See DLAHR2.)```

Definition at line 166 of file dgehrd.f.

167*
168* -- LAPACK computational routine --
169* -- LAPACK is a software package provided by Univ. of Tennessee, --
170* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
171*
172* .. Scalar Arguments ..
173 INTEGER IHI, ILO, INFO, LDA, LWORK, N
174* ..
175* .. Array Arguments ..
176 DOUBLE PRECISION A( LDA, * ), TAU( * ), WORK( * )
177* ..
178*
179* =====================================================================
180*
181* .. Parameters ..
182 INTEGER NBMAX, LDT, TSIZE
183 parameter( nbmax = 64, ldt = nbmax+1,
184 \$ tsize = ldt*nbmax )
185 DOUBLE PRECISION ZERO, ONE
186 parameter( zero = 0.0d+0,
187 \$ one = 1.0d+0 )
188* ..
189* .. Local Scalars ..
190 LOGICAL LQUERY
191 INTEGER I, IB, IINFO, IWT, J, LDWORK, LWKOPT, NB,
192 \$ NBMIN, NH, NX
193 DOUBLE PRECISION EI
194* ..
195* .. External Subroutines ..
196 EXTERNAL daxpy, dgehd2, dgemm, dlahr2, dlarfb, dtrmm,
197 \$ xerbla
198* ..
199* .. Intrinsic Functions ..
200 INTRINSIC max, min
201* ..
202* .. External Functions ..
203 INTEGER ILAENV
204 EXTERNAL ilaenv
205* ..
206* .. Executable Statements ..
207*
208* Test the input parameters
209*
210 info = 0
211 lquery = ( lwork.EQ.-1 )
212 IF( n.LT.0 ) THEN
213 info = -1
214 ELSE IF( ilo.LT.1 .OR. ilo.GT.max( 1, n ) ) THEN
215 info = -2
216 ELSE IF( ihi.LT.min( ilo, n ) .OR. ihi.GT.n ) THEN
217 info = -3
218 ELSE IF( lda.LT.max( 1, n ) ) THEN
219 info = -5
220 ELSE IF( lwork.LT.max( 1, n ) .AND. .NOT.lquery ) THEN
221 info = -8
222 END IF
223*
224 IF( info.EQ.0 ) THEN
225*
226* Compute the workspace requirements
227*
228 nb = min( nbmax, ilaenv( 1, 'DGEHRD', ' ', n, ilo, ihi, -1 ) )
229 lwkopt = n*nb + tsize
230 work( 1 ) = lwkopt
231 END IF
232*
233 IF( info.NE.0 ) THEN
234 CALL xerbla( 'DGEHRD', -info )
235 RETURN
236 ELSE IF( lquery ) THEN
237 RETURN
238 END IF
239*
240* Set elements 1:ILO-1 and IHI:N-1 of TAU to zero
241*
242 DO 10 i = 1, ilo - 1
243 tau( i ) = zero
244 10 CONTINUE
245 DO 20 i = max( 1, ihi ), n - 1
246 tau( i ) = zero
247 20 CONTINUE
248*
249* Quick return if possible
250*
251 nh = ihi - ilo + 1
252 IF( nh.LE.1 ) THEN
253 work( 1 ) = 1
254 RETURN
255 END IF
256*
257* Determine the block size
258*
259 nb = min( nbmax, ilaenv( 1, 'DGEHRD', ' ', n, ilo, ihi, -1 ) )
260 nbmin = 2
261 IF( nb.GT.1 .AND. nb.LT.nh ) THEN
262*
263* Determine when to cross over from blocked to unblocked code
264* (last block is always handled by unblocked code)
265*
266 nx = max( nb, ilaenv( 3, 'DGEHRD', ' ', n, ilo, ihi, -1 ) )
267 IF( nx.LT.nh ) THEN
268*
269* Determine if workspace is large enough for blocked code
270*
271 IF( lwork.LT.n*nb+tsize ) THEN
272*
273* Not enough workspace to use optimal NB: determine the
274* minimum value of NB, and reduce NB or force use of
275* unblocked code
276*
277 nbmin = max( 2, ilaenv( 2, 'DGEHRD', ' ', n, ilo, ihi,
278 \$ -1 ) )
279 IF( lwork.GE.(n*nbmin + tsize) ) THEN
280 nb = (lwork-tsize) / n
281 ELSE
282 nb = 1
283 END IF
284 END IF
285 END IF
286 END IF
287 ldwork = n
288*
289 IF( nb.LT.nbmin .OR. nb.GE.nh ) THEN
290*
291* Use unblocked code below
292*
293 i = ilo
294*
295 ELSE
296*
297* Use blocked code
298*
299 iwt = 1 + n*nb
300 DO 40 i = ilo, ihi - 1 - nx, nb
301 ib = min( nb, ihi-i )
302*
303* Reduce columns i:i+ib-1 to Hessenberg form, returning the
304* matrices V and T of the block reflector H = I - V*T*V**T
305* which performs the reduction, and also the matrix Y = A*V*T
306*
307 CALL dlahr2( ihi, i, ib, a( 1, i ), lda, tau( i ),
308 \$ work( iwt ), ldt, work, ldwork )
309*
310* Apply the block reflector H to A(1:ihi,i+ib:ihi) from the
311* right, computing A := A - Y * V**T. V(i+ib,ib-1) must be set
312* to 1
313*
314 ei = a( i+ib, i+ib-1 )
315 a( i+ib, i+ib-1 ) = one
316 CALL dgemm( 'No transpose', 'Transpose',
317 \$ ihi, ihi-i-ib+1,
318 \$ ib, -one, work, ldwork, a( i+ib, i ), lda, one,
319 \$ a( 1, i+ib ), lda )
320 a( i+ib, i+ib-1 ) = ei
321*
322* Apply the block reflector H to A(1:i,i+1:i+ib-1) from the
323* right
324*
325 CALL dtrmm( 'Right', 'Lower', 'Transpose',
326 \$ 'Unit', i, ib-1,
327 \$ one, a( i+1, i ), lda, work, ldwork )
328 DO 30 j = 0, ib-2
329 CALL daxpy( i, -one, work( ldwork*j+1 ), 1,
330 \$ a( 1, i+j+1 ), 1 )
331 30 CONTINUE
332*
333* Apply the block reflector H to A(i+1:ihi,i+ib:n) from the
334* left
335*
336 CALL dlarfb( 'Left', 'Transpose', 'Forward',
337 \$ 'Columnwise',
338 \$ ihi-i, n-i-ib+1, ib, a( i+1, i ), lda,
339 \$ work( iwt ), ldt, a( i+1, i+ib ), lda,
340 \$ work, ldwork )
341 40 CONTINUE
342 END IF
343*
344* Use unblocked code to reduce the rest of the matrix
345*
346 CALL dgehd2( n, i, ihi, a, lda, tau, work, iinfo )
347 work( 1 ) = lwkopt
348*
349 RETURN
350*
351* End of DGEHRD
352*
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine daxpy(N, DA, DX, INCX, DY, INCY)
DAXPY
Definition: daxpy.f:89
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:187
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM
Definition: dtrmm.f:177
subroutine dgehd2(N, ILO, IHI, A, LDA, TAU, WORK, INFO)
DGEHD2 reduces a general square matrix to upper Hessenberg form using an unblocked algorithm.
Definition: dgehd2.f:149
subroutine dlarfb(SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV, T, LDT, C, LDC, WORK, LDWORK)
DLARFB applies a block reflector or its transpose to a general rectangular matrix.
Definition: dlarfb.f:197
subroutine dlahr2(N, K, NB, A, LDA, TAU, T, LDT, Y, LDY)
DLAHR2 reduces the specified number of first columns of a general rectangular matrix A so that elemen...
Definition: dlahr2.f:181
Here is the call graph for this function:
Here is the caller graph for this function: