![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b ZUNT03 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE ZUNT03( RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK, 00012 * RWORK, RESULT, INFO ) 00013 * 00014 * .. Scalar Arguments .. 00015 * CHARACTER*( * ) RC 00016 * INTEGER INFO, K, LDU, LDV, LWORK, MU, MV, N 00017 * DOUBLE PRECISION RESULT 00018 * .. 00019 * .. Array Arguments .. 00020 * DOUBLE PRECISION RWORK( * ) 00021 * COMPLEX*16 U( LDU, * ), V( LDV, * ), WORK( * ) 00022 * .. 00023 * 00024 * 00025 *> \par Purpose: 00026 * ============= 00027 *> 00028 *> \verbatim 00029 *> 00030 *> ZUNT03 compares two unitary matrices U and V to see if their 00031 *> corresponding rows or columns span the same spaces. The rows are 00032 *> checked if RC = 'R', and the columns are checked if RC = 'C'. 00033 *> 00034 *> RESULT is the maximum of 00035 *> 00036 *> | V*V' - I | / ( MV ulp ), if RC = 'R', or 00037 *> 00038 *> | V'*V - I | / ( MV ulp ), if RC = 'C', 00039 *> 00040 *> and the maximum over rows (or columns) 1 to K of 00041 *> 00042 *> | U(i) - S*V(i) |/ ( N ulp ) 00043 *> 00044 *> where abs(S) = 1 (chosen to minimize the expression), U(i) is the 00045 *> i-th row (column) of U, and V(i) is the i-th row (column) of V. 00046 *> \endverbatim 00047 * 00048 * Arguments: 00049 * ========== 00050 * 00051 *> \param[in] RC 00052 *> \verbatim 00053 *> RC is CHARACTER*1 00054 *> If RC = 'R' the rows of U and V are to be compared. 00055 *> If RC = 'C' the columns of U and V are to be compared. 00056 *> \endverbatim 00057 *> 00058 *> \param[in] MU 00059 *> \verbatim 00060 *> MU is INTEGER 00061 *> The number of rows of U if RC = 'R', and the number of 00062 *> columns if RC = 'C'. If MU = 0 ZUNT03 does nothing. 00063 *> MU must be at least zero. 00064 *> \endverbatim 00065 *> 00066 *> \param[in] MV 00067 *> \verbatim 00068 *> MV is INTEGER 00069 *> The number of rows of V if RC = 'R', and the number of 00070 *> columns if RC = 'C'. If MV = 0 ZUNT03 does nothing. 00071 *> MV must be at least zero. 00072 *> \endverbatim 00073 *> 00074 *> \param[in] N 00075 *> \verbatim 00076 *> N is INTEGER 00077 *> If RC = 'R', the number of columns in the matrices U and V, 00078 *> and if RC = 'C', the number of rows in U and V. If N = 0 00079 *> ZUNT03 does nothing. N must be at least zero. 00080 *> \endverbatim 00081 *> 00082 *> \param[in] K 00083 *> \verbatim 00084 *> K is INTEGER 00085 *> The number of rows or columns of U and V to compare. 00086 *> 0 <= K <= max(MU,MV). 00087 *> \endverbatim 00088 *> 00089 *> \param[in] U 00090 *> \verbatim 00091 *> U is COMPLEX*16 array, dimension (LDU,N) 00092 *> The first matrix to compare. If RC = 'R', U is MU by N, and 00093 *> if RC = 'C', U is N by MU. 00094 *> \endverbatim 00095 *> 00096 *> \param[in] LDU 00097 *> \verbatim 00098 *> LDU is INTEGER 00099 *> The leading dimension of U. If RC = 'R', LDU >= max(1,MU), 00100 *> and if RC = 'C', LDU >= max(1,N). 00101 *> \endverbatim 00102 *> 00103 *> \param[in] V 00104 *> \verbatim 00105 *> V is COMPLEX*16 array, dimension (LDV,N) 00106 *> The second matrix to compare. If RC = 'R', V is MV by N, and 00107 *> if RC = 'C', V is N by MV. 00108 *> \endverbatim 00109 *> 00110 *> \param[in] LDV 00111 *> \verbatim 00112 *> LDV is INTEGER 00113 *> The leading dimension of V. If RC = 'R', LDV >= max(1,MV), 00114 *> and if RC = 'C', LDV >= max(1,N). 00115 *> \endverbatim 00116 *> 00117 *> \param[out] WORK 00118 *> \verbatim 00119 *> WORK is COMPLEX*16 array, dimension (LWORK) 00120 *> \endverbatim 00121 *> 00122 *> \param[in] LWORK 00123 *> \verbatim 00124 *> LWORK is INTEGER 00125 *> The length of the array WORK. For best performance, LWORK 00126 *> should be at least N*N if RC = 'C' or M*M if RC = 'R', but 00127 *> the tests will be done even if LWORK is 0. 00128 *> \endverbatim 00129 *> 00130 *> \param[out] RWORK 00131 *> \verbatim 00132 *> RWORK is DOUBLE PRECISION array, dimension (max(MV,N)) 00133 *> \endverbatim 00134 *> 00135 *> \param[out] RESULT 00136 *> \verbatim 00137 *> RESULT is DOUBLE PRECISION 00138 *> The value computed by the test described above. RESULT is 00139 *> limited to 1/ulp to avoid overflow. 00140 *> \endverbatim 00141 *> 00142 *> \param[out] INFO 00143 *> \verbatim 00144 *> INFO is INTEGER 00145 *> 0 indicates a successful exit 00146 *> -k indicates the k-th parameter had an illegal value 00147 *> \endverbatim 00148 * 00149 * Authors: 00150 * ======== 00151 * 00152 *> \author Univ. of Tennessee 00153 *> \author Univ. of California Berkeley 00154 *> \author Univ. of Colorado Denver 00155 *> \author NAG Ltd. 00156 * 00157 *> \date November 2011 00158 * 00159 *> \ingroup complex16_eig 00160 * 00161 * ===================================================================== 00162 SUBROUTINE ZUNT03( RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK, 00163 $ RWORK, RESULT, INFO ) 00164 * 00165 * -- LAPACK test routine (version 3.4.0) -- 00166 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00167 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00168 * November 2011 00169 * 00170 * .. Scalar Arguments .. 00171 CHARACTER*( * ) RC 00172 INTEGER INFO, K, LDU, LDV, LWORK, MU, MV, N 00173 DOUBLE PRECISION RESULT 00174 * .. 00175 * .. Array Arguments .. 00176 DOUBLE PRECISION RWORK( * ) 00177 COMPLEX*16 U( LDU, * ), V( LDV, * ), WORK( * ) 00178 * .. 00179 * 00180 * ===================================================================== 00181 * 00182 * 00183 * .. Parameters .. 00184 DOUBLE PRECISION ZERO, ONE 00185 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) 00186 * .. 00187 * .. Local Scalars .. 00188 INTEGER I, IRC, J, LMX 00189 DOUBLE PRECISION RES1, RES2, ULP 00190 COMPLEX*16 S, SU, SV 00191 * .. 00192 * .. External Functions .. 00193 LOGICAL LSAME 00194 INTEGER IZAMAX 00195 DOUBLE PRECISION DLAMCH 00196 EXTERNAL LSAME, IZAMAX, DLAMCH 00197 * .. 00198 * .. Intrinsic Functions .. 00199 INTRINSIC ABS, DBLE, DCMPLX, MAX, MIN 00200 * .. 00201 * .. External Subroutines .. 00202 EXTERNAL XERBLA, ZUNT01 00203 * .. 00204 * .. Executable Statements .. 00205 * 00206 * Check inputs 00207 * 00208 INFO = 0 00209 IF( LSAME( RC, 'R' ) ) THEN 00210 IRC = 0 00211 ELSE IF( LSAME( RC, 'C' ) ) THEN 00212 IRC = 1 00213 ELSE 00214 IRC = -1 00215 END IF 00216 IF( IRC.EQ.-1 ) THEN 00217 INFO = -1 00218 ELSE IF( MU.LT.0 ) THEN 00219 INFO = -2 00220 ELSE IF( MV.LT.0 ) THEN 00221 INFO = -3 00222 ELSE IF( N.LT.0 ) THEN 00223 INFO = -4 00224 ELSE IF( K.LT.0 .OR. K.GT.MAX( MU, MV ) ) THEN 00225 INFO = -5 00226 ELSE IF( ( IRC.EQ.0 .AND. LDU.LT.MAX( 1, MU ) ) .OR. 00227 $ ( IRC.EQ.1 .AND. LDU.LT.MAX( 1, N ) ) ) THEN 00228 INFO = -7 00229 ELSE IF( ( IRC.EQ.0 .AND. LDV.LT.MAX( 1, MV ) ) .OR. 00230 $ ( IRC.EQ.1 .AND. LDV.LT.MAX( 1, N ) ) ) THEN 00231 INFO = -9 00232 END IF 00233 IF( INFO.NE.0 ) THEN 00234 CALL XERBLA( 'ZUNT03', -INFO ) 00235 RETURN 00236 END IF 00237 * 00238 * Initialize result 00239 * 00240 RESULT = ZERO 00241 IF( MU.EQ.0 .OR. MV.EQ.0 .OR. N.EQ.0 ) 00242 $ RETURN 00243 * 00244 * Machine constants 00245 * 00246 ULP = DLAMCH( 'Precision' ) 00247 * 00248 IF( IRC.EQ.0 ) THEN 00249 * 00250 * Compare rows 00251 * 00252 RES1 = ZERO 00253 DO 20 I = 1, K 00254 LMX = IZAMAX( N, U( I, 1 ), LDU ) 00255 IF( V( I, LMX ).EQ.DCMPLX( ZERO ) ) THEN 00256 SV = ONE 00257 ELSE 00258 SV = ABS( V( I, LMX ) ) / V( I, LMX ) 00259 END IF 00260 IF( U( I, LMX ).EQ.DCMPLX( ZERO ) ) THEN 00261 SU = ONE 00262 ELSE 00263 SU = ABS( U( I, LMX ) ) / U( I, LMX ) 00264 END IF 00265 S = SV / SU 00266 DO 10 J = 1, N 00267 RES1 = MAX( RES1, ABS( U( I, J )-S*V( I, J ) ) ) 00268 10 CONTINUE 00269 20 CONTINUE 00270 RES1 = RES1 / ( DBLE( N )*ULP ) 00271 * 00272 * Compute orthogonality of rows of V. 00273 * 00274 CALL ZUNT01( 'Rows', MV, N, V, LDV, WORK, LWORK, RWORK, RES2 ) 00275 * 00276 ELSE 00277 * 00278 * Compare columns 00279 * 00280 RES1 = ZERO 00281 DO 40 I = 1, K 00282 LMX = IZAMAX( N, U( 1, I ), 1 ) 00283 IF( V( LMX, I ).EQ.DCMPLX( ZERO ) ) THEN 00284 SV = ONE 00285 ELSE 00286 SV = ABS( V( LMX, I ) ) / V( LMX, I ) 00287 END IF 00288 IF( U( LMX, I ).EQ.DCMPLX( ZERO ) ) THEN 00289 SU = ONE 00290 ELSE 00291 SU = ABS( U( LMX, I ) ) / U( LMX, I ) 00292 END IF 00293 S = SV / SU 00294 DO 30 J = 1, N 00295 RES1 = MAX( RES1, ABS( U( J, I )-S*V( J, I ) ) ) 00296 30 CONTINUE 00297 40 CONTINUE 00298 RES1 = RES1 / ( DBLE( N )*ULP ) 00299 * 00300 * Compute orthogonality of columns of V. 00301 * 00302 CALL ZUNT01( 'Columns', N, MV, V, LDV, WORK, LWORK, RWORK, 00303 $ RES2 ) 00304 END IF 00305 * 00306 RESULT = MIN( MAX( RES1, RES2 ), ONE / ULP ) 00307 RETURN 00308 * 00309 * End of ZUNT03 00310 * 00311 END