116
117
118
119
120
121
122 CHARACTER TRANS
123 INTEGER LDA, LDX, LWORK, M, N, NRHS
124
125
126 COMPLEX A( LDA, * ), WORK( LWORK ), X( LDX, * )
127
128
129
130
131
132 REAL ZERO, ONE
133 parameter( zero = 0.0e0, one = 1.0e0 )
134
135
136 LOGICAL TPSD
137 INTEGER I, INFO, J, LDWORK
138 REAL ANRM, ERR, XNRM
139
140
141 REAL RWORK( 1 )
142
143
144 LOGICAL LSAME
145 REAL CLANGE, SLAMCH
147
148
150
151
152 INTRINSIC abs, conjg, max, min, real
153
154
155
157 IF(
lsame( trans,
'N' ) )
THEN
158 ldwork = m + nrhs
159 tpsd = .false.
160 IF( lwork.LT.( m+nrhs )*( n+2 ) ) THEN
161 CALL xerbla(
'CQRT14', 10 )
162 RETURN
163 ELSE IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
164 RETURN
165 END IF
166 ELSE IF(
lsame( trans,
'C' ) )
THEN
167 ldwork = m
168 tpsd = .true.
169 IF( lwork.LT.( n+nrhs )*( m+2 ) ) THEN
170 CALL xerbla(
'CQRT14', 10 )
171 RETURN
172 ELSE IF( m.LE.0 .OR. nrhs.LE.0 ) THEN
173 RETURN
174 END IF
175 ELSE
176 CALL xerbla(
'CQRT14', 1 )
177 RETURN
178 END IF
179
180
181
182 CALL clacpy(
'All', m, n, a, lda, work, ldwork )
183 anrm =
clange(
'M', m, n, work, ldwork, rwork )
184 IF( anrm.NE.zero )
185 $
CALL clascl(
'G', 0, 0, anrm, one, m, n, work, ldwork, info )
186
187
188
189 IF( tpsd ) THEN
190
191
192
193 CALL clacpy(
'All', m, nrhs, x, ldx, work( n*ldwork+1 ),
194 $ ldwork )
195 xnrm =
clange(
'M', m, nrhs, work( n*ldwork+1 ), ldwork,
196 $ rwork )
197 IF( xnrm.NE.zero )
198 $
CALL clascl(
'G', 0, 0, xnrm, one, m, nrhs,
199 $ work( n*ldwork+1 ), ldwork, info )
200
201
202
203 CALL cgeqr2( m, n+nrhs, work, ldwork,
204 $ work( ldwork*( n+nrhs )+1 ),
205 $ work( ldwork*( n+nrhs )+min( m, n+nrhs )+1 ),
206 $ info )
207
208
209
210
211 err = zero
212 DO 20 j = n + 1, n + nrhs
213 DO 10 i = n + 1, min( m, j )
214 err = max( err, abs( work( i+( j-1 )*m ) ) )
215 10 CONTINUE
216 20 CONTINUE
217
218 ELSE
219
220
221
222 DO 40 i = 1, n
223 DO 30 j = 1, nrhs
224 work( m+j+( i-1 )*ldwork ) = conjg( x( i, j ) )
225 30 CONTINUE
226 40 CONTINUE
227
228 xnrm =
clange(
'M', nrhs, n, work( m+1 ), ldwork, rwork )
229 IF( xnrm.NE.zero )
230 $
CALL clascl(
'G', 0, 0, xnrm, one, nrhs, n, work( m+1 ),
231 $ ldwork, info )
232
233
234
235 CALL cgelq2( ldwork, n, work, ldwork, work( ldwork*n+1 ),
236 $ work( ldwork*( n+1 )+1 ), info )
237
238
239
240
241 err = zero
242 DO 60 j = m + 1, n
243 DO 50 i = j, ldwork
244 err = max( err, abs( work( i+( j-1 )*ldwork ) ) )
245 50 CONTINUE
246 60 CONTINUE
247
248 END IF
249
250 cqrt14 = err / ( real( max( m, n, nrhs ) )*
slamch(
'Epsilon' ) )
251
252 RETURN
253
254
255
subroutine xerbla(srname, info)
real function cqrt14(trans, m, n, nrhs, a, lda, x, ldx, work, lwork)
CQRT14
subroutine cgelq2(m, n, a, lda, tau, work, info)
CGELQ2 computes the LQ factorization of a general rectangular matrix using an unblocked algorithm.
subroutine cgeqr2(m, n, a, lda, tau, work, info)
CGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm.
subroutine clacpy(uplo, m, n, a, lda, b, ldb)
CLACPY copies all or part of one two-dimensional array to another.
real function slamch(cmach)
SLAMCH
real function clange(norm, m, n, a, lda, work)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine clascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
logical function lsame(ca, cb)
LSAME