![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b CQRT14 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * REAL FUNCTION CQRT14( TRANS, M, N, NRHS, A, LDA, X, 00012 * LDX, WORK, LWORK ) 00013 * 00014 * .. Scalar Arguments .. 00015 * CHARACTER TRANS 00016 * INTEGER LDA, LDX, LWORK, M, N, NRHS 00017 * .. 00018 * .. Array Arguments .. 00019 * COMPLEX A( LDA, * ), WORK( LWORK ), X( LDX, * ) 00020 * .. 00021 * 00022 * 00023 *> \par Purpose: 00024 * ============= 00025 *> 00026 *> \verbatim 00027 *> 00028 *> CQRT14 checks whether X is in the row space of A or A'. It does so 00029 *> by scaling both X and A such that their norms are in the range 00030 *> [sqrt(eps), 1/sqrt(eps)], then computing a QR factorization of [A,X] 00031 *> (if TRANS = 'C') or an LQ factorization of [A',X]' (if TRANS = 'N'), 00032 *> and returning the norm of the trailing triangle, scaled by 00033 *> MAX(M,N,NRHS)*eps. 00034 *> \endverbatim 00035 * 00036 * Arguments: 00037 * ========== 00038 * 00039 *> \param[in] TRANS 00040 *> \verbatim 00041 *> TRANS is CHARACTER*1 00042 *> = 'N': No transpose, check for X in the row space of A 00043 *> = 'C': Conjugate transpose, check for X in row space of A'. 00044 *> \endverbatim 00045 *> 00046 *> \param[in] M 00047 *> \verbatim 00048 *> M is INTEGER 00049 *> The number of rows of the matrix A. 00050 *> \endverbatim 00051 *> 00052 *> \param[in] N 00053 *> \verbatim 00054 *> N is INTEGER 00055 *> The number of columns of the matrix A. 00056 *> \endverbatim 00057 *> 00058 *> \param[in] NRHS 00059 *> \verbatim 00060 *> NRHS is INTEGER 00061 *> The number of right hand sides, i.e., the number of columns 00062 *> of X. 00063 *> \endverbatim 00064 *> 00065 *> \param[in] A 00066 *> \verbatim 00067 *> A is COMPLEX array, dimension (LDA,N) 00068 *> The M-by-N matrix A. 00069 *> \endverbatim 00070 *> 00071 *> \param[in] LDA 00072 *> \verbatim 00073 *> LDA is INTEGER 00074 *> The leading dimension of the array A. 00075 *> \endverbatim 00076 *> 00077 *> \param[in] X 00078 *> \verbatim 00079 *> X is COMPLEX array, dimension (LDX,NRHS) 00080 *> If TRANS = 'N', the N-by-NRHS matrix X. 00081 *> IF TRANS = 'C', the M-by-NRHS matrix X. 00082 *> \endverbatim 00083 *> 00084 *> \param[in] LDX 00085 *> \verbatim 00086 *> LDX is INTEGER 00087 *> The leading dimension of the array X. 00088 *> \endverbatim 00089 *> 00090 *> \param[out] WORK 00091 *> \verbatim 00092 *> WORK is COMPLEX array dimension (LWORK) 00093 *> \endverbatim 00094 *> 00095 *> \param[in] LWORK 00096 *> \verbatim 00097 *> LWORK is INTEGER 00098 *> length of workspace array required 00099 *> If TRANS = 'N', LWORK >= (M+NRHS)*(N+2); 00100 *> if TRANS = 'C', LWORK >= (N+NRHS)*(M+2). 00101 *> \endverbatim 00102 * 00103 * Authors: 00104 * ======== 00105 * 00106 *> \author Univ. of Tennessee 00107 *> \author Univ. of California Berkeley 00108 *> \author Univ. of Colorado Denver 00109 *> \author NAG Ltd. 00110 * 00111 *> \date November 2011 00112 * 00113 *> \ingroup complex_lin 00114 * 00115 * ===================================================================== 00116 REAL FUNCTION CQRT14( TRANS, M, N, NRHS, A, LDA, X, 00117 $ LDX, WORK, LWORK ) 00118 * 00119 * -- LAPACK test routine (version 3.4.0) -- 00120 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00121 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00122 * November 2011 00123 * 00124 * .. Scalar Arguments .. 00125 CHARACTER TRANS 00126 INTEGER LDA, LDX, LWORK, M, N, NRHS 00127 * .. 00128 * .. Array Arguments .. 00129 COMPLEX A( LDA, * ), WORK( LWORK ), X( LDX, * ) 00130 * .. 00131 * 00132 * ===================================================================== 00133 * 00134 * .. Parameters .. 00135 REAL ZERO, ONE 00136 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0 ) 00137 * .. 00138 * .. Local Scalars .. 00139 LOGICAL TPSD 00140 INTEGER I, INFO, J, LDWORK 00141 REAL ANRM, ERR, XNRM 00142 * .. 00143 * .. Local Arrays .. 00144 REAL RWORK( 1 ) 00145 * .. 00146 * .. External Functions .. 00147 LOGICAL LSAME 00148 REAL CLANGE, SLAMCH 00149 EXTERNAL LSAME, CLANGE, SLAMCH 00150 * .. 00151 * .. External Subroutines .. 00152 EXTERNAL CGELQ2, CGEQR2, CLACPY, CLASCL, XERBLA 00153 * .. 00154 * .. Intrinsic Functions .. 00155 INTRINSIC ABS, CONJG, MAX, MIN, REAL 00156 * .. 00157 * .. Executable Statements .. 00158 * 00159 CQRT14 = ZERO 00160 IF( LSAME( TRANS, 'N' ) ) THEN 00161 LDWORK = M + NRHS 00162 TPSD = .FALSE. 00163 IF( LWORK.LT.( M+NRHS )*( N+2 ) ) THEN 00164 CALL XERBLA( 'CQRT14', 10 ) 00165 RETURN 00166 ELSE IF( N.LE.0 .OR. NRHS.LE.0 ) THEN 00167 RETURN 00168 END IF 00169 ELSE IF( LSAME( TRANS, 'C' ) ) THEN 00170 LDWORK = M 00171 TPSD = .TRUE. 00172 IF( LWORK.LT.( N+NRHS )*( M+2 ) ) THEN 00173 CALL XERBLA( 'CQRT14', 10 ) 00174 RETURN 00175 ELSE IF( M.LE.0 .OR. NRHS.LE.0 ) THEN 00176 RETURN 00177 END IF 00178 ELSE 00179 CALL XERBLA( 'CQRT14', 1 ) 00180 RETURN 00181 END IF 00182 * 00183 * Copy and scale A 00184 * 00185 CALL CLACPY( 'All', M, N, A, LDA, WORK, LDWORK ) 00186 ANRM = CLANGE( 'M', M, N, WORK, LDWORK, RWORK ) 00187 IF( ANRM.NE.ZERO ) 00188 $ CALL CLASCL( 'G', 0, 0, ANRM, ONE, M, N, WORK, LDWORK, INFO ) 00189 * 00190 * Copy X or X' into the right place and scale it 00191 * 00192 IF( TPSD ) THEN 00193 * 00194 * Copy X into columns n+1:n+nrhs of work 00195 * 00196 CALL CLACPY( 'All', M, NRHS, X, LDX, WORK( N*LDWORK+1 ), 00197 $ LDWORK ) 00198 XNRM = CLANGE( 'M', M, NRHS, WORK( N*LDWORK+1 ), LDWORK, 00199 $ RWORK ) 00200 IF( XNRM.NE.ZERO ) 00201 $ CALL CLASCL( 'G', 0, 0, XNRM, ONE, M, NRHS, 00202 $ WORK( N*LDWORK+1 ), LDWORK, INFO ) 00203 ANRM = CLANGE( 'One-norm', M, N+NRHS, WORK, LDWORK, RWORK ) 00204 * 00205 * Compute QR factorization of X 00206 * 00207 CALL CGEQR2( M, N+NRHS, WORK, LDWORK, 00208 $ WORK( LDWORK*( N+NRHS )+1 ), 00209 $ WORK( LDWORK*( N+NRHS )+MIN( M, N+NRHS )+1 ), 00210 $ INFO ) 00211 * 00212 * Compute largest entry in upper triangle of 00213 * work(n+1:m,n+1:n+nrhs) 00214 * 00215 ERR = ZERO 00216 DO 20 J = N + 1, N + NRHS 00217 DO 10 I = N + 1, MIN( M, J ) 00218 ERR = MAX( ERR, ABS( WORK( I+( J-1 )*M ) ) ) 00219 10 CONTINUE 00220 20 CONTINUE 00221 * 00222 ELSE 00223 * 00224 * Copy X' into rows m+1:m+nrhs of work 00225 * 00226 DO 40 I = 1, N 00227 DO 30 J = 1, NRHS 00228 WORK( M+J+( I-1 )*LDWORK ) = CONJG( X( I, J ) ) 00229 30 CONTINUE 00230 40 CONTINUE 00231 * 00232 XNRM = CLANGE( 'M', NRHS, N, WORK( M+1 ), LDWORK, RWORK ) 00233 IF( XNRM.NE.ZERO ) 00234 $ CALL CLASCL( 'G', 0, 0, XNRM, ONE, NRHS, N, WORK( M+1 ), 00235 $ LDWORK, INFO ) 00236 * 00237 * Compute LQ factorization of work 00238 * 00239 CALL CGELQ2( LDWORK, N, WORK, LDWORK, WORK( LDWORK*N+1 ), 00240 $ WORK( LDWORK*( N+1 )+1 ), INFO ) 00241 * 00242 * Compute largest entry in lower triangle in 00243 * work(m+1:m+nrhs,m+1:n) 00244 * 00245 ERR = ZERO 00246 DO 60 J = M + 1, N 00247 DO 50 I = J, LDWORK 00248 ERR = MAX( ERR, ABS( WORK( I+( J-1 )*LDWORK ) ) ) 00249 50 CONTINUE 00250 60 CONTINUE 00251 * 00252 END IF 00253 * 00254 CQRT14 = ERR / ( REAL( MAX( M, N, NRHS ) )*SLAMCH( 'Epsilon' ) ) 00255 * 00256 RETURN 00257 * 00258 * End of CQRT14 00259 * 00260 END