zunt01()

subroutine zunt01	(	character	rowcol,
		integer	m,
		integer	n,
		complex*16, dimension( ldu, * )	u,
		integer	ldu,
		complex*16, dimension( * )	work,
		integer	lwork,
		double precision, dimension( * )	rwork,
		double precision	resid )

ZUNT01

Purpose:

!>
!> ZUNT01 checks that the matrix U is unitary by computing the ratio
!>
!>    RESID = norm( I - U*U' ) / ( n * EPS ), if ROWCOL = 'R',
!> or
!>    RESID = norm( I - U'*U ) / ( m * EPS ), if ROWCOL = 'C'.
!>
!> Alternatively, if there isn't sufficient workspace to form
!> I - U*U' or I - U'*U, the ratio is computed as
!>
!>    RESID = abs( I - U*U' ) / ( n * EPS ), if ROWCOL = 'R',
!> or
!>    RESID = abs( I - U'*U ) / ( m * EPS ), if ROWCOL = 'C'.
!>
!> where EPS is the machine precision.  ROWCOL is used only if m = n;
!> if m > n, ROWCOL is assumed to be 'C', and if m < n, ROWCOL is
!> assumed to be 'R'.
!> 

Parameters

[in]	ROWCOL	!> ROWCOL is CHARACTER !> Specifies whether the rows or columns of U should be checked !> for orthogonality. Used only if M = N. !> = 'R': Check for orthogonal rows of U !> = 'C': Check for orthogonal columns of U !>
[in]	M	!> M is INTEGER !> The number of rows of the matrix U. !>
[in]	N	!> N is INTEGER !> The number of columns of the matrix U. !>
[in]	U	!> U is COMPLEX*16 array, dimension (LDU,N) !> The unitary matrix U. U is checked for orthogonal columns !> if m > n or if m = n and ROWCOL = 'C'. U is checked for !> orthogonal rows if m < n or if m = n and ROWCOL = 'R'. !>
[in]	LDU	!> LDU is INTEGER !> The leading dimension of the array U. LDU >= max(1,M). !>
[out]	WORK	!> WORK is COMPLEX*16 array, dimension (LWORK) !>
[in]	LWORK	!> LWORK is INTEGER !> The length of the array WORK. For best performance, LWORK !> should be at least N*N if ROWCOL = 'C' or M*M if !> ROWCOL = 'R', but the test will be done even if LWORK is 0. !>
[out]	RWORK	!> RWORK is DOUBLE PRECISION array, dimension (min(M,N)) !> Used only if LWORK is large enough to use the Level 3 BLAS !> code. !>
[out]	RESID	!> RESID is DOUBLE PRECISION !> RESID = norm( I - U * U' ) / ( n * EPS ), if ROWCOL = 'R', or !> RESID = norm( I - U' * U ) / ( m * EPS ), if ROWCOL = 'C'. !>

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Definition at line 124 of file zunt01.f.

*
*  -- LAPACK test routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      CHARACTER          ROWCOL
      INTEGER            LDU, LWORK, M, N
      DOUBLE PRECISION   RESID
*     ..
*     .. Array Arguments ..
      DOUBLE PRECISION   RWORK( * )
      COMPLEX*16         U( LDU, * ), WORK( * )
*     ..
*
*  =====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   ZERO, ONE
      parameter( zero = 0.0d+0, one = 1.0d+0 )
*     ..
*     .. Local Scalars ..
      CHARACTER          TRANSU
      INTEGER            I, J, K, LDWORK, MNMIN
      DOUBLE PRECISION   EPS
      COMPLEX*16         TMP, ZDUM
*     ..
*     .. External Functions ..
      LOGICAL            LSAME
      DOUBLE PRECISION   DLAMCH, ZLANSY
      COMPLEX*16         ZDOTC
      EXTERNAL           lsame, dlamch, zlansy, zdotc
*     ..
*     .. External Subroutines ..
      EXTERNAL           zherk, zlaset
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, dble, dcmplx, dimag, max, min
*     ..
*     .. Statement Functions ..
      DOUBLE PRECISION   CABS1
*     ..
*     .. Statement Function definitions ..
      cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
*     ..
*     .. Executable Statements ..
*
      resid = zero
*
*     Quick return if possible
*
      IF( m.LE.0 .OR. n.LE.0 )
     $   RETURN
*
      eps = dlamch( 'Precision' )
      IF( m.LT.n .OR. ( m.EQ.n .AND. lsame( rowcol, 'R' ) ) ) THEN
         transu = 'N'
         k = n
      ELSE
         transu = 'C'
         k = m
      END IF
      mnmin = min( m, n )
*
      IF( ( mnmin+1 )*mnmin.LE.lwork ) THEN
         ldwork = mnmin
      ELSE
         ldwork = 0
      END IF
      IF( ldwork.GT.0 ) THEN
*
*        Compute I - U*U' or I - U'*U.
*
         CALL zlaset( 'Upper', mnmin, mnmin, dcmplx( zero ),
     $                dcmplx( one ), work, ldwork )
         CALL zherk( 'Upper', transu, mnmin, k, -one, u, ldu, one, work,
     $               ldwork )
*
*        Compute norm( I - U*U' ) / ( K * EPS ) .
*
         resid = zlansy( '1', 'Upper', mnmin, work, ldwork, rwork )
         resid = ( resid / dble( k ) ) / eps
      ELSE IF( transu.EQ.'C' ) THEN
*
*        Find the maximum element in abs( I - U'*U ) / ( m * EPS )
*
         DO 20 j = 1, n
            DO 10 i = 1, j
               IF( i.NE.j ) THEN
                  tmp = zero
               ELSE
                  tmp = one
               END IF
               tmp = tmp - zdotc( m, u( 1, i ), 1, u( 1, j ), 1 )
               resid = max( resid, cabs1( tmp ) )
   10       CONTINUE
   20    CONTINUE
         resid = ( resid / dble( m ) ) / eps
      ELSE
*
*        Find the maximum element in abs( I - U*U' ) / ( n * EPS )
*
         DO 40 j = 1, m
            DO 30 i = 1, j
               IF( i.NE.j ) THEN
                  tmp = zero
               ELSE
                  tmp = one
               END IF
               tmp = tmp - zdotc( n, u( j, 1 ), ldu, u( i, 1 ), ldu )
               resid = max( resid, cabs1( tmp ) )
   30       CONTINUE
   40    CONTINUE
         resid = ( resid / dble( n ) ) / eps
      END IF
      RETURN
*
*     End of ZUNT01
*

Here is the call graph for this function:

Here is the caller graph for this function: