 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE ZUNT03( RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK, 00002 $ RWORK, RESULT, INFO ) 00003 * 00004 * -- LAPACK test routine (version 3.1) -- 00005 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 CHARACTER*( * ) RC 00010 INTEGER INFO, K, LDU, LDV, LWORK, MU, MV, N 00011 DOUBLE PRECISION RESULT 00012 * .. 00013 * .. Array Arguments .. 00014 DOUBLE PRECISION RWORK( * ) 00015 COMPLEX*16 U( LDU, * ), V( LDV, * ), WORK( * ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * ZUNT03 compares two unitary matrices U and V to see if their 00022 * corresponding rows or columns span the same spaces. The rows are 00023 * checked if RC = 'R', and the columns are checked if RC = 'C'. 00024 * 00025 * RESULT is the maximum of 00026 * 00027 * | V*V' - I | / ( MV ulp ), if RC = 'R', or 00028 * 00029 * | V'*V - I | / ( MV ulp ), if RC = 'C', 00030 * 00031 * and the maximum over rows (or columns) 1 to K of 00032 * 00033 * | U(i) - S*V(i) |/ ( N ulp ) 00034 * 00035 * where abs(S) = 1 (chosen to minimize the expression), U(i) is the 00036 * i-th row (column) of U, and V(i) is the i-th row (column) of V. 00037 * 00038 * Arguments 00039 * ========== 00040 * 00041 * RC (input) CHARACTER*1 00042 * If RC = 'R' the rows of U and V are to be compared. 00043 * If RC = 'C' the columns of U and V are to be compared. 00044 * 00045 * MU (input) INTEGER 00046 * The number of rows of U if RC = 'R', and the number of 00047 * columns if RC = 'C'. If MU = 0 ZUNT03 does nothing. 00048 * MU must be at least zero. 00049 * 00050 * MV (input) INTEGER 00051 * The number of rows of V if RC = 'R', and the number of 00052 * columns if RC = 'C'. If MV = 0 ZUNT03 does nothing. 00053 * MV must be at least zero. 00054 * 00055 * N (input) INTEGER 00056 * If RC = 'R', the number of columns in the matrices U and V, 00057 * and if RC = 'C', the number of rows in U and V. If N = 0 00058 * ZUNT03 does nothing. N must be at least zero. 00059 * 00060 * K (input) INTEGER 00061 * The number of rows or columns of U and V to compare. 00062 * 0 <= K <= max(MU,MV). 00063 * 00064 * U (input) COMPLEX*16 array, dimension (LDU,N) 00065 * The first matrix to compare. If RC = 'R', U is MU by N, and 00066 * if RC = 'C', U is N by MU. 00067 * 00068 * LDU (input) INTEGER 00069 * The leading dimension of U. If RC = 'R', LDU >= max(1,MU), 00070 * and if RC = 'C', LDU >= max(1,N). 00071 * 00072 * V (input) COMPLEX*16 array, dimension (LDV,N) 00073 * The second matrix to compare. If RC = 'R', V is MV by N, and 00074 * if RC = 'C', V is N by MV. 00075 * 00076 * LDV (input) INTEGER 00077 * The leading dimension of V. If RC = 'R', LDV >= max(1,MV), 00078 * and if RC = 'C', LDV >= max(1,N). 00079 * 00080 * WORK (workspace) COMPLEX*16 array, dimension (LWORK) 00081 * 00082 * LWORK (input) INTEGER 00083 * The length of the array WORK. For best performance, LWORK 00084 * should be at least N*N if RC = 'C' or M*M if RC = 'R', but 00085 * the tests will be done even if LWORK is 0. 00086 * 00087 * RWORK (workspace) DOUBLE PRECISION array, dimension (max(MV,N)) 00088 * 00089 * RESULT (output) DOUBLE PRECISION 00090 * The value computed by the test described above. RESULT is 00091 * limited to 1/ulp to avoid overflow. 00092 * 00093 * INFO (output) INTEGER 00094 * 0 indicates a successful exit 00095 * -k indicates the k-th parameter had an illegal value 00096 * 00097 * ===================================================================== 00098 * 00099 * 00100 * .. Parameters .. 00101 DOUBLE PRECISION ZERO, ONE 00102 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) 00103 * .. 00104 * .. Local Scalars .. 00105 INTEGER I, IRC, J, LMX 00106 DOUBLE PRECISION RES1, RES2, ULP 00107 COMPLEX*16 S, SU, SV 00108 * .. 00109 * .. External Functions .. 00110 LOGICAL LSAME 00111 INTEGER IZAMAX 00112 DOUBLE PRECISION DLAMCH 00113 EXTERNAL LSAME, IZAMAX, DLAMCH 00114 * .. 00115 * .. Intrinsic Functions .. 00116 INTRINSIC ABS, DBLE, DCMPLX, MAX, MIN 00117 * .. 00118 * .. External Subroutines .. 00119 EXTERNAL XERBLA, ZUNT01 00120 * .. 00121 * .. Executable Statements .. 00122 * 00123 * Check inputs 00124 * 00125 INFO = 0 00126 IF( LSAME( RC, 'R' ) ) THEN 00127 IRC = 0 00128 ELSE IF( LSAME( RC, 'C' ) ) THEN 00129 IRC = 1 00130 ELSE 00131 IRC = -1 00132 END IF 00133 IF( IRC.EQ.-1 ) THEN 00134 INFO = -1 00135 ELSE IF( MU.LT.0 ) THEN 00136 INFO = -2 00137 ELSE IF( MV.LT.0 ) THEN 00138 INFO = -3 00139 ELSE IF( N.LT.0 ) THEN 00140 INFO = -4 00141 ELSE IF( K.LT.0 .OR. K.GT.MAX( MU, MV ) ) THEN 00142 INFO = -5 00143 ELSE IF( ( IRC.EQ.0 .AND. LDU.LT.MAX( 1, MU ) ) .OR. 00144 $ ( IRC.EQ.1 .AND. LDU.LT.MAX( 1, N ) ) ) THEN 00145 INFO = -7 00146 ELSE IF( ( IRC.EQ.0 .AND. LDV.LT.MAX( 1, MV ) ) .OR. 00147 $ ( IRC.EQ.1 .AND. LDV.LT.MAX( 1, N ) ) ) THEN 00148 INFO = -9 00149 END IF 00150 IF( INFO.NE.0 ) THEN 00151 CALL XERBLA( 'ZUNT03', -INFO ) 00152 RETURN 00153 END IF 00154 * 00155 * Initialize result 00156 * 00157 RESULT = ZERO 00158 IF( MU.EQ.0 .OR. MV.EQ.0 .OR. N.EQ.0 ) 00159 $ RETURN 00160 * 00161 * Machine constants 00162 * 00163 ULP = DLAMCH( 'Precision' ) 00164 * 00165 IF( IRC.EQ.0 ) THEN 00166 * 00167 * Compare rows 00168 * 00169 RES1 = ZERO 00170 DO 20 I = 1, K 00171 LMX = IZAMAX( N, U( I, 1 ), LDU ) 00172 IF( V( I, LMX ).EQ.DCMPLX( ZERO ) ) THEN 00173 SV = ONE 00174 ELSE 00175 SV = ABS( V( I, LMX ) ) / V( I, LMX ) 00176 END IF 00177 IF( U( I, LMX ).EQ.DCMPLX( ZERO ) ) THEN 00178 SU = ONE 00179 ELSE 00180 SU = ABS( U( I, LMX ) ) / U( I, LMX ) 00181 END IF 00182 S = SV / SU 00183 DO 10 J = 1, N 00184 RES1 = MAX( RES1, ABS( U( I, J )-S*V( I, J ) ) ) 00185 10 CONTINUE 00186 20 CONTINUE 00187 RES1 = RES1 / ( DBLE( N )*ULP ) 00188 * 00189 * Compute orthogonality of rows of V. 00190 * 00191 CALL ZUNT01( 'Rows', MV, N, V, LDV, WORK, LWORK, RWORK, RES2 ) 00192 * 00193 ELSE 00194 * 00195 * Compare columns 00196 * 00197 RES1 = ZERO 00198 DO 40 I = 1, K 00199 LMX = IZAMAX( N, U( 1, I ), 1 ) 00200 IF( V( LMX, I ).EQ.DCMPLX( ZERO ) ) THEN 00201 SV = ONE 00202 ELSE 00203 SV = ABS( V( LMX, I ) ) / V( LMX, I ) 00204 END IF 00205 IF( U( LMX, I ).EQ.DCMPLX( ZERO ) ) THEN 00206 SU = ONE 00207 ELSE 00208 SU = ABS( U( LMX, I ) ) / U( LMX, I ) 00209 END IF 00210 S = SV / SU 00211 DO 30 J = 1, N 00212 RES1 = MAX( RES1, ABS( U( J, I )-S*V( J, I ) ) ) 00213 30 CONTINUE 00214 40 CONTINUE 00215 RES1 = RES1 / ( DBLE( N )*ULP ) 00216 * 00217 * Compute orthogonality of columns of V. 00218 * 00219 CALL ZUNT01( 'Columns', N, MV, V, LDV, WORK, LWORK, RWORK, 00220 $ RES2 ) 00221 END IF 00222 * 00223 RESULT = MIN( MAX( RES1, RES2 ), ONE / ULP ) 00224 RETURN 00225 * 00226 * End of ZUNT03 00227 * 00228 END
1.7.3 
