149
150
151
152
153
154
155 INTEGER LDA, LDB, LWORK, M, N, NRHS, RANK, RKSEL, SCALE
156 REAL NORMA, NORMB
157
158
159 INTEGER ISEED( 4 )
160 REAL S( * )
161 COMPLEX A( LDA, * ), B( LDB, * ), WORK( LWORK )
162
163
164
165
166
167 REAL ZERO, ONE, TWO, SVMIN
168 parameter( zero = 0.0e+0, one = 1.0e+0, two = 2.0e+0,
169 $ svmin = 0.1e+0 )
170 COMPLEX CZERO, CONE
171 parameter( czero = ( 0.0e+0, 0.0e+0 ),
172 $ cone = ( 1.0e+0, 0.0e+0 ) )
173
174
175 INTEGER INFO, J, MN
176 REAL BIGNUM, EPS, SMLNUM, TEMP
177
178
179 REAL DUMMY( 1 )
180
181
182 REAL CLANGE, SASUM, SCNRM2, SLAMCH, SLARND
184
185
188
189
190 INTRINSIC abs, cmplx, max, min
191
192
193
194 mn = min( m, n )
195 IF( lwork.LT.max( m+mn, mn*nrhs, 2*n+m ) ) THEN
196 CALL xerbla(
'CQRT15', 16 )
197 RETURN
198 END IF
199
200 smlnum =
slamch(
'Safe minimum' )
201 bignum = one / smlnum
203 smlnum = ( smlnum / eps ) / eps
204 bignum = one / smlnum
205
206
207
208 IF( rksel.EQ.1 ) THEN
209 rank = mn
210 ELSE IF( rksel.EQ.2 ) THEN
211 rank = ( 3*mn ) / 4
212 DO 10 j = rank + 1, mn
213 s( j ) = zero
214 10 CONTINUE
215 ELSE
216 CALL xerbla(
'CQRT15', 2 )
217 END IF
218
219 IF( rank.GT.0 ) THEN
220
221
222
223 s( 1 ) = one
224 DO 30 j = 2, rank
225 20 CONTINUE
227 IF( temp.GT.svmin ) THEN
228 s( j ) = abs( temp )
229 ELSE
230 GO TO 20
231 END IF
232 30 CONTINUE
233 CALL slaord(
'Decreasing', rank, s, 1 )
234
235
236
237 CALL clarnv( 2, iseed, m, work )
239 CALL claset(
'Full', m, rank, czero, cone, a, lda )
240 CALL clarf(
'Left', m, rank, work, 1, cmplx( two ), a, lda,
241 $ work( m+1 ) )
242
243
244
245
246
247 CALL clarnv( 2, iseed, rank*nrhs, work )
248 CALL cgemm(
'No transpose',
'No transpose', m, nrhs, rank,
249 $ cone, a, lda, work, rank, czero, b, ldb )
250
251
252
253
254
255 DO 40 j = 1, rank
256 CALL csscal( m, s( j ), a( 1, j ), 1 )
257 40 CONTINUE
258 IF( rank.LT.n )
259 $
CALL claset(
'Full', m, n-rank, czero, czero,
260 $ a( 1, rank+1 ), lda )
261 CALL claror(
'Right',
'No initialization', m, n, a, lda, iseed,
262 $ work, info )
263
264 ELSE
265
266
267
268
269
270 DO 50 j = 1, mn
271 s( j ) = zero
272 50 CONTINUE
273 CALL claset(
'Full', m, n, czero, czero, a, lda )
274 CALL claset(
'Full', m, nrhs, czero, czero, b, ldb )
275
276 END IF
277
278
279
280 IF( scale.NE.1 ) THEN
281 norma =
clange(
'Max', m, n, a, lda, dummy )
282 IF( norma.NE.zero ) THEN
283 IF( scale.EQ.2 ) THEN
284
285
286
287 CALL clascl(
'General', 0, 0, norma, bignum, m, n, a,
288 $ lda, info )
289 CALL slascl(
'General', 0, 0, norma, bignum, mn, 1, s,
290 $ mn, info )
291 CALL clascl(
'General', 0, 0, norma, bignum, m, nrhs, b,
292 $ ldb, info )
293 ELSE IF( scale.EQ.3 ) THEN
294
295
296
297 CALL clascl(
'General', 0, 0, norma, smlnum, m, n, a,
298 $ lda, info )
299 CALL slascl(
'General', 0, 0, norma, smlnum, mn, 1, s,
300 $ mn, info )
301 CALL clascl(
'General', 0, 0, norma, smlnum, m, nrhs, b,
302 $ ldb, info )
303 ELSE
304 CALL xerbla(
'CQRT15', 1 )
305 RETURN
306 END IF
307 END IF
308 END IF
309
310 norma =
sasum( mn, s, 1 )
311 normb =
clange(
'One-norm', m, nrhs, b, ldb, dummy )
312
313 RETURN
314
315
316
subroutine xerbla(srname, info)
subroutine claror(side, init, m, n, a, lda, iseed, x, info)
CLAROR
real function sasum(n, sx, incx)
SASUM
subroutine cgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
CGEMM
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 clarf(side, m, n, v, incv, tau, c, ldc, work)
CLARF applies an elementary reflector to a general rectangular matrix.
subroutine clarnv(idist, iseed, n, x)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
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.
subroutine slascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine claset(uplo, m, n, alpha, beta, a, lda)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
real(wp) function scnrm2(n, x, incx)
SCNRM2
subroutine csscal(n, sa, cx, incx)
CSSCAL
subroutine slaord(job, n, x, incx)
SLAORD
real function slarnd(idist, iseed)
SLARND