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

## ◆ zhet21()

 subroutine zhet21 ( integer ITYPE, character UPLO, integer N, integer KBAND, complex*16, dimension( lda, * ) A, integer LDA, double precision, dimension( * ) D, double precision, dimension( * ) E, complex*16, dimension( ldu, * ) U, integer LDU, complex*16, dimension( ldv, * ) V, integer LDV, complex*16, dimension( * ) TAU, complex*16, dimension( * ) WORK, double precision, dimension( * ) RWORK, double precision, dimension( 2 ) RESULT )

ZHET21

Purpose:
``` ZHET21 generally checks a decomposition of the form

A = U S U**H

where **H means conjugate transpose, A is hermitian, U is unitary, and
S is diagonal (if KBAND=0) or (real) symmetric tridiagonal (if
KBAND=1).

If ITYPE=1, then U is represented as a dense matrix; otherwise U is
expressed as a product of Householder transformations, whose vectors
are stored in the array "V" and whose scaling constants are in "TAU".
We shall use the letter "V" to refer to the product of Householder
transformations (which should be equal to U).

Specifically, if ITYPE=1, then:

RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
RESULT(2) = | I - U U**H | / ( n ulp )

If ITYPE=2, then:

RESULT(1) = | A - V S V**H | / ( |A| n ulp )

If ITYPE=3, then:

RESULT(1) = | I - U V**H | / ( n ulp )

For ITYPE > 1, the transformation U is expressed as a product
V = H(1)...H(n-2),  where H(j) = I  -  tau(j) v(j) v(j)**H and each
vector v(j) has its first j elements 0 and the remaining n-j elements
stored in V(j+1:n,j).```
Parameters
 [in] ITYPE ``` ITYPE is INTEGER Specifies the type of tests to be performed. 1: U expressed as a dense unitary matrix: RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and RESULT(2) = | I - U U**H | / ( n ulp ) 2: U expressed as a product V of Housholder transformations: RESULT(1) = | A - V S V**H | / ( |A| n ulp ) 3: U expressed both as a dense unitary matrix and as a product of Housholder transformations: RESULT(1) = | I - U V**H | / ( n ulp )``` [in] UPLO ``` UPLO is CHARACTER If UPLO='U', the upper triangle of A and V will be used and the (strictly) lower triangle will not be referenced. If UPLO='L', the lower triangle of A and V will be used and the (strictly) upper triangle will not be referenced.``` [in] N ``` N is INTEGER The size of the matrix. If it is zero, ZHET21 does nothing. It must be at least zero.``` [in] KBAND ``` KBAND is INTEGER The bandwidth of the matrix. It may only be zero or one. If zero, then S is diagonal, and E is not referenced. If one, then S is symmetric tri-diagonal.``` [in] A ``` A is COMPLEX*16 array, dimension (LDA, N) The original (unfactored) matrix. It is assumed to be hermitian, and only the upper (UPLO='U') or only the lower (UPLO='L') will be referenced.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. It must be at least 1 and at least N.``` [in] D ``` D is DOUBLE PRECISION array, dimension (N) The diagonal of the (symmetric tri-) diagonal matrix.``` [in] E ``` E is DOUBLE PRECISION array, dimension (N-1) The off-diagonal of the (symmetric tri-) diagonal matrix. E(1) is the (1,2) and (2,1) element, E(2) is the (2,3) and (3,2) element, etc. Not referenced if KBAND=0.``` [in] U ``` U is COMPLEX*16 array, dimension (LDU, N) If ITYPE=1 or 3, this contains the unitary matrix in the decomposition, expressed as a dense matrix. If ITYPE=2, then it is not referenced.``` [in] LDU ``` LDU is INTEGER The leading dimension of U. LDU must be at least N and at least 1.``` [in] V ``` V is COMPLEX*16 array, dimension (LDV, N) If ITYPE=2 or 3, the columns of this array contain the Householder vectors used to describe the unitary matrix in the decomposition. If UPLO='L', then the vectors are in the lower triangle, if UPLO='U', then in the upper triangle. *NOTE* If ITYPE=2 or 3, V is modified and restored. The subdiagonal (if UPLO='L') or the superdiagonal (if UPLO='U') is set to one, and later reset to its original value, during the course of the calculation. If ITYPE=1, then it is neither referenced nor modified.``` [in] LDV ``` LDV is INTEGER The leading dimension of V. LDV must be at least N and at least 1.``` [in] TAU ``` TAU is COMPLEX*16 array, dimension (N) If ITYPE >= 2, then TAU(j) is the scalar factor of v(j) v(j)**H in the Householder transformation H(j) of the product U = H(1)...H(n-2) If ITYPE < 2, then TAU is not referenced.``` [out] WORK ` WORK is COMPLEX*16 array, dimension (2*N**2)` [out] RWORK ` RWORK is DOUBLE PRECISION array, dimension (N)` [out] RESULT ``` RESULT is DOUBLE PRECISION array, dimension (2) The values computed by the two tests described above. The values are currently limited to 1/ulp, to avoid overflow. RESULT(1) is always modified. RESULT(2) is modified only if ITYPE=1.```

Definition at line 212 of file zhet21.f.

214*
215* -- LAPACK test routine --
216* -- LAPACK is a software package provided by Univ. of Tennessee, --
217* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
218*
219* .. Scalar Arguments ..
220 CHARACTER UPLO
221 INTEGER ITYPE, KBAND, LDA, LDU, LDV, N
222* ..
223* .. Array Arguments ..
224 DOUBLE PRECISION D( * ), E( * ), RESULT( 2 ), RWORK( * )
225 COMPLEX*16 A( LDA, * ), TAU( * ), U( LDU, * ),
226 \$ V( LDV, * ), WORK( * )
227* ..
228*
229* =====================================================================
230*
231* .. Parameters ..
232 DOUBLE PRECISION ZERO, ONE, TEN
233 parameter( zero = 0.0d+0, one = 1.0d+0, ten = 10.0d+0 )
234 COMPLEX*16 CZERO, CONE
235 parameter( czero = ( 0.0d+0, 0.0d+0 ),
236 \$ cone = ( 1.0d+0, 0.0d+0 ) )
237* ..
238* .. Local Scalars ..
239 LOGICAL LOWER
240 CHARACTER CUPLO
241 INTEGER IINFO, J, JCOL, JR, JROW
242 DOUBLE PRECISION ANORM, ULP, UNFL, WNORM
243 COMPLEX*16 VSAVE
244* ..
245* .. External Functions ..
246 LOGICAL LSAME
247 DOUBLE PRECISION DLAMCH, ZLANGE, ZLANHE
248 EXTERNAL lsame, dlamch, zlange, zlanhe
249* ..
250* .. External Subroutines ..
251 EXTERNAL zgemm, zher, zher2, zlacpy, zlarfy, zlaset,
252 \$ zunm2l, zunm2r
253* ..
254* .. Intrinsic Functions ..
255 INTRINSIC dble, dcmplx, max, min
256* ..
257* .. Executable Statements ..
258*
259 result( 1 ) = zero
260 IF( itype.EQ.1 )
261 \$ result( 2 ) = zero
262 IF( n.LE.0 )
263 \$ RETURN
264*
265 IF( lsame( uplo, 'U' ) ) THEN
266 lower = .false.
267 cuplo = 'U'
268 ELSE
269 lower = .true.
270 cuplo = 'L'
271 END IF
272*
273 unfl = dlamch( 'Safe minimum' )
274 ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
275*
276* Some Error Checks
277*
278 IF( itype.LT.1 .OR. itype.GT.3 ) THEN
279 result( 1 ) = ten / ulp
280 RETURN
281 END IF
282*
283* Do Test 1
284*
285* Norm of A:
286*
287 IF( itype.EQ.3 ) THEN
288 anorm = one
289 ELSE
290 anorm = max( zlanhe( '1', cuplo, n, a, lda, rwork ), unfl )
291 END IF
292*
293* Compute error matrix:
294*
295 IF( itype.EQ.1 ) THEN
296*
297* ITYPE=1: error = A - U S U**H
298*
299 CALL zlaset( 'Full', n, n, czero, czero, work, n )
300 CALL zlacpy( cuplo, n, n, a, lda, work, n )
301*
302 DO 10 j = 1, n
303 CALL zher( cuplo, n, -d( j ), u( 1, j ), 1, work, n )
304 10 CONTINUE
305*
306 IF( n.GT.1 .AND. kband.EQ.1 ) THEN
307 DO 20 j = 2, n - 1
308 CALL zher2( cuplo, n, -dcmplx( e( j ) ), u( 1, j ), 1,
309 \$ u( 1, j-1 ), 1, work, n )
310 20 CONTINUE
311 END IF
312 wnorm = zlanhe( '1', cuplo, n, work, n, rwork )
313*
314 ELSE IF( itype.EQ.2 ) THEN
315*
316* ITYPE=2: error = V S V**H - A
317*
318 CALL zlaset( 'Full', n, n, czero, czero, work, n )
319*
320 IF( lower ) THEN
321 work( n**2 ) = d( n )
322 DO 40 j = n - 1, 1, -1
323 IF( kband.EQ.1 ) THEN
324 work( ( n+1 )*( j-1 )+2 ) = ( cone-tau( j ) )*e( j )
325 DO 30 jr = j + 2, n
326 work( ( j-1 )*n+jr ) = -tau( j )*e( j )*v( jr, j )
327 30 CONTINUE
328 END IF
329*
330 vsave = v( j+1, j )
331 v( j+1, j ) = one
332 CALL zlarfy( 'L', n-j, v( j+1, j ), 1, tau( j ),
333 \$ work( ( n+1 )*j+1 ), n, work( n**2+1 ) )
334 v( j+1, j ) = vsave
335 work( ( n+1 )*( j-1 )+1 ) = d( j )
336 40 CONTINUE
337 ELSE
338 work( 1 ) = d( 1 )
339 DO 60 j = 1, n - 1
340 IF( kband.EQ.1 ) THEN
341 work( ( n+1 )*j ) = ( cone-tau( j ) )*e( j )
342 DO 50 jr = 1, j - 1
343 work( j*n+jr ) = -tau( j )*e( j )*v( jr, j+1 )
344 50 CONTINUE
345 END IF
346*
347 vsave = v( j, j+1 )
348 v( j, j+1 ) = one
349 CALL zlarfy( 'U', j, v( 1, j+1 ), 1, tau( j ), work, n,
350 \$ work( n**2+1 ) )
351 v( j, j+1 ) = vsave
352 work( ( n+1 )*j+1 ) = d( j+1 )
353 60 CONTINUE
354 END IF
355*
356 DO 90 jcol = 1, n
357 IF( lower ) THEN
358 DO 70 jrow = jcol, n
359 work( jrow+n*( jcol-1 ) ) = work( jrow+n*( jcol-1 ) )
360 \$ - a( jrow, jcol )
361 70 CONTINUE
362 ELSE
363 DO 80 jrow = 1, jcol
364 work( jrow+n*( jcol-1 ) ) = work( jrow+n*( jcol-1 ) )
365 \$ - a( jrow, jcol )
366 80 CONTINUE
367 END IF
368 90 CONTINUE
369 wnorm = zlanhe( '1', cuplo, n, work, n, rwork )
370*
371 ELSE IF( itype.EQ.3 ) THEN
372*
373* ITYPE=3: error = U V**H - I
374*
375 IF( n.LT.2 )
376 \$ RETURN
377 CALL zlacpy( ' ', n, n, u, ldu, work, n )
378 IF( lower ) THEN
379 CALL zunm2r( 'R', 'C', n, n-1, n-1, v( 2, 1 ), ldv, tau,
380 \$ work( n+1 ), n, work( n**2+1 ), iinfo )
381 ELSE
382 CALL zunm2l( 'R', 'C', n, n-1, n-1, v( 1, 2 ), ldv, tau,
383 \$ work, n, work( n**2+1 ), iinfo )
384 END IF
385 IF( iinfo.NE.0 ) THEN
386 result( 1 ) = ten / ulp
387 RETURN
388 END IF
389*
390 DO 100 j = 1, n
391 work( ( n+1 )*( j-1 )+1 ) = work( ( n+1 )*( j-1 )+1 ) - cone
392 100 CONTINUE
393*
394 wnorm = zlange( '1', n, n, work, n, rwork )
395 END IF
396*
397 IF( anorm.GT.wnorm ) THEN
398 result( 1 ) = ( wnorm / anorm ) / ( n*ulp )
399 ELSE
400 IF( anorm.LT.one ) THEN
401 result( 1 ) = ( min( wnorm, n*anorm ) / anorm ) / ( n*ulp )
402 ELSE
403 result( 1 ) = min( wnorm / anorm, dble( n ) ) / ( n*ulp )
404 END IF
405 END IF
406*
407* Do Test 2
408*
409* Compute U U**H - I
410*
411 IF( itype.EQ.1 ) THEN
412 CALL zgemm( 'N', 'C', n, n, n, cone, u, ldu, u, ldu, czero,
413 \$ work, n )
414*
415 DO 110 j = 1, n
416 work( ( n+1 )*( j-1 )+1 ) = work( ( n+1 )*( j-1 )+1 ) - cone
417 110 CONTINUE
418*
419 result( 2 ) = min( zlange( '1', n, n, work, n, rwork ),
420 \$ dble( n ) ) / ( n*ulp )
421 END IF
422*
423 RETURN
424*
425* End of ZHET21
426*
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine zher(UPLO, N, ALPHA, X, INCX, A, LDA)
ZHER
Definition: zher.f:135
subroutine zher2(UPLO, N, ALPHA, X, INCX, Y, INCY, A, LDA)
ZHER2
Definition: zher2.f:150
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:187
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: zlange.f:115
double precision function zlanhe(NORM, UPLO, N, A, LDA, WORK)
ZLANHE returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: zlanhe.f:124
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103
subroutine zlarfy(UPLO, N, V, INCV, TAU, C, LDC, WORK)
ZLARFY
Definition: zlarfy.f:108
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: zlaset.f:106
subroutine zunm2l(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
ZUNM2L multiplies a general matrix by the unitary matrix from a QL factorization determined by cgeqlf...
Definition: zunm2l.f:159
subroutine zunm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
ZUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
Definition: zunm2r.f:159
Here is the call graph for this function:
Here is the caller graph for this function: