LAPACK  3.10.1
LAPACK: Linear Algebra PACKage

◆ zlalsd()

subroutine zlalsd ( character  UPLO,
integer  SMLSIZ,
integer  N,
integer  NRHS,
double precision, dimension( * )  D,
double precision, dimension( * )  E,
complex*16, dimension( ldb, * )  B,
integer  LDB,
double precision  RCOND,
integer  RANK,
complex*16, dimension( * )  WORK,
double precision, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer  INFO 
)

ZLALSD uses the singular value decomposition of A to solve the least squares problem.

Download ZLALSD + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZLALSD uses the singular value decomposition of A to solve the least
 squares problem of finding X to minimize the Euclidean norm of each
 column of A*X-B, where A is N-by-N upper bidiagonal, and X and B
 are N-by-NRHS. The solution X overwrites B.

 The singular values of A smaller than RCOND times the largest
 singular value are treated as zero in solving the least squares
 problem; in this case a minimum norm solution is returned.
 The actual singular values are returned in D in ascending order.

 This code makes very mild assumptions about floating point
 arithmetic. It will work on machines with a guard digit in
 add/subtract, or on those binary machines without guard digits
 which subtract like the Cray XMP, Cray YMP, Cray C 90, or Cray 2.
 It could conceivably fail on hexadecimal or decimal machines
 without guard digits, but we know of none.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
         = 'U': D and E define an upper bidiagonal matrix.
         = 'L': D and E define a  lower bidiagonal matrix.
[in]SMLSIZ
          SMLSIZ is INTEGER
         The maximum size of the subproblems at the bottom of the
         computation tree.
[in]N
          N is INTEGER
         The dimension of the  bidiagonal matrix.  N >= 0.
[in]NRHS
          NRHS is INTEGER
         The number of columns of B. NRHS must be at least 1.
[in,out]D
          D is DOUBLE PRECISION array, dimension (N)
         On entry D contains the main diagonal of the bidiagonal
         matrix. On exit, if INFO = 0, D contains its singular values.
[in,out]E
          E is DOUBLE PRECISION array, dimension (N-1)
         Contains the super-diagonal entries of the bidiagonal matrix.
         On exit, E has been destroyed.
[in,out]B
          B is COMPLEX*16 array, dimension (LDB,NRHS)
         On input, B contains the right hand sides of the least
         squares problem. On output, B contains the solution X.
[in]LDB
          LDB is INTEGER
         The leading dimension of B in the calling subprogram.
         LDB must be at least max(1,N).
[in]RCOND
          RCOND is DOUBLE PRECISION
         The singular values of A less than or equal to RCOND times
         the largest singular value are treated as zero in solving
         the least squares problem. If RCOND is negative,
         machine precision is used instead.
         For example, if diag(S)*X=B were the least squares problem,
         where diag(S) is a diagonal matrix of singular values, the
         solution would be X(i) = B(i) / S(i) if S(i) is greater than
         RCOND*max(S), and X(i) = 0 if S(i) is less than or equal to
         RCOND*max(S).
[out]RANK
          RANK is INTEGER
         The number of singular values of A greater than RCOND times
         the largest singular value.
[out]WORK
          WORK is COMPLEX*16 array, dimension (N * NRHS)
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension at least
         (9*N + 2*N*SMLSIZ + 8*N*NLVL + 3*SMLSIZ*NRHS +
         MAX( (SMLSIZ+1)**2, N*(1+NRHS) + 2*NRHS ),
         where
         NLVL = MAX( 0, INT( LOG_2( MIN( M,N )/(SMLSIZ+1) ) ) + 1 )
[out]IWORK
          IWORK is INTEGER array, dimension at least
         (3*N*NLVL + 11*N).
[out]INFO
          INFO is INTEGER
         = 0:  successful exit.
         < 0:  if INFO = -i, the i-th argument had an illegal value.
         > 0:  The algorithm failed to compute a singular value while
               working on the submatrix lying in rows and columns
               INFO/(N+1) through MOD(INFO,N+1).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Contributors:
Ming Gu and Ren-Cang Li, Computer Science Division, University of California at Berkeley, USA
Osni Marques, LBNL/NERSC, USA

Definition at line 185 of file zlalsd.f.

187 *
188 * -- LAPACK computational routine --
189 * -- LAPACK is a software package provided by Univ. of Tennessee, --
190 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
191 *
192 * .. Scalar Arguments ..
193  CHARACTER UPLO
194  INTEGER INFO, LDB, N, NRHS, RANK, SMLSIZ
195  DOUBLE PRECISION RCOND
196 * ..
197 * .. Array Arguments ..
198  INTEGER IWORK( * )
199  DOUBLE PRECISION D( * ), E( * ), RWORK( * )
200  COMPLEX*16 B( LDB, * ), WORK( * )
201 * ..
202 *
203 * =====================================================================
204 *
205 * .. Parameters ..
206  DOUBLE PRECISION ZERO, ONE, TWO
207  parameter( zero = 0.0d0, one = 1.0d0, two = 2.0d0 )
208  COMPLEX*16 CZERO
209  parameter( czero = ( 0.0d0, 0.0d0 ) )
210 * ..
211 * .. Local Scalars ..
212  INTEGER BX, BXST, C, DIFL, DIFR, GIVCOL, GIVNUM,
213  $ GIVPTR, I, ICMPQ1, ICMPQ2, IRWB, IRWIB, IRWRB,
214  $ IRWU, IRWVT, IRWWRK, IWK, J, JCOL, JIMAG,
215  $ JREAL, JROW, K, NLVL, NM1, NRWORK, NSIZE, NSUB,
216  $ PERM, POLES, S, SIZEI, SMLSZP, SQRE, ST, ST1,
217  $ U, VT, Z
218  DOUBLE PRECISION CS, EPS, ORGNRM, RCND, R, SN, TOL
219 * ..
220 * .. External Functions ..
221  INTEGER IDAMAX
222  DOUBLE PRECISION DLAMCH, DLANST
223  EXTERNAL idamax, dlamch, dlanst
224 * ..
225 * .. External Subroutines ..
226  EXTERNAL dgemm, dlartg, dlascl, dlasda, dlasdq, dlaset,
228  $ zlascl, zlaset
229 * ..
230 * .. Intrinsic Functions ..
231  INTRINSIC abs, dble, dcmplx, dimag, int, log, sign
232 * ..
233 * .. Executable Statements ..
234 *
235 * Test the input parameters.
236 *
237  info = 0
238 *
239  IF( n.LT.0 ) THEN
240  info = -3
241  ELSE IF( nrhs.LT.1 ) THEN
242  info = -4
243  ELSE IF( ( ldb.LT.1 ) .OR. ( ldb.LT.n ) ) THEN
244  info = -8
245  END IF
246  IF( info.NE.0 ) THEN
247  CALL xerbla( 'ZLALSD', -info )
248  RETURN
249  END IF
250 *
251  eps = dlamch( 'Epsilon' )
252 *
253 * Set up the tolerance.
254 *
255  IF( ( rcond.LE.zero ) .OR. ( rcond.GE.one ) ) THEN
256  rcnd = eps
257  ELSE
258  rcnd = rcond
259  END IF
260 *
261  rank = 0
262 *
263 * Quick return if possible.
264 *
265  IF( n.EQ.0 ) THEN
266  RETURN
267  ELSE IF( n.EQ.1 ) THEN
268  IF( d( 1 ).EQ.zero ) THEN
269  CALL zlaset( 'A', 1, nrhs, czero, czero, b, ldb )
270  ELSE
271  rank = 1
272  CALL zlascl( 'G', 0, 0, d( 1 ), one, 1, nrhs, b, ldb, info )
273  d( 1 ) = abs( d( 1 ) )
274  END IF
275  RETURN
276  END IF
277 *
278 * Rotate the matrix if it is lower bidiagonal.
279 *
280  IF( uplo.EQ.'L' ) THEN
281  DO 10 i = 1, n - 1
282  CALL dlartg( d( i ), e( i ), cs, sn, r )
283  d( i ) = r
284  e( i ) = sn*d( i+1 )
285  d( i+1 ) = cs*d( i+1 )
286  IF( nrhs.EQ.1 ) THEN
287  CALL zdrot( 1, b( i, 1 ), 1, b( i+1, 1 ), 1, cs, sn )
288  ELSE
289  rwork( i*2-1 ) = cs
290  rwork( i*2 ) = sn
291  END IF
292  10 CONTINUE
293  IF( nrhs.GT.1 ) THEN
294  DO 30 i = 1, nrhs
295  DO 20 j = 1, n - 1
296  cs = rwork( j*2-1 )
297  sn = rwork( j*2 )
298  CALL zdrot( 1, b( j, i ), 1, b( j+1, i ), 1, cs, sn )
299  20 CONTINUE
300  30 CONTINUE
301  END IF
302  END IF
303 *
304 * Scale.
305 *
306  nm1 = n - 1
307  orgnrm = dlanst( 'M', n, d, e )
308  IF( orgnrm.EQ.zero ) THEN
309  CALL zlaset( 'A', n, nrhs, czero, czero, b, ldb )
310  RETURN
311  END IF
312 *
313  CALL dlascl( 'G', 0, 0, orgnrm, one, n, 1, d, n, info )
314  CALL dlascl( 'G', 0, 0, orgnrm, one, nm1, 1, e, nm1, info )
315 *
316 * If N is smaller than the minimum divide size SMLSIZ, then solve
317 * the problem with another solver.
318 *
319  IF( n.LE.smlsiz ) THEN
320  irwu = 1
321  irwvt = irwu + n*n
322  irwwrk = irwvt + n*n
323  irwrb = irwwrk
324  irwib = irwrb + n*nrhs
325  irwb = irwib + n*nrhs
326  CALL dlaset( 'A', n, n, zero, one, rwork( irwu ), n )
327  CALL dlaset( 'A', n, n, zero, one, rwork( irwvt ), n )
328  CALL dlasdq( 'U', 0, n, n, n, 0, d, e, rwork( irwvt ), n,
329  $ rwork( irwu ), n, rwork( irwwrk ), 1,
330  $ rwork( irwwrk ), info )
331  IF( info.NE.0 ) THEN
332  RETURN
333  END IF
334 *
335 * In the real version, B is passed to DLASDQ and multiplied
336 * internally by Q**H. Here B is complex and that product is
337 * computed below in two steps (real and imaginary parts).
338 *
339  j = irwb - 1
340  DO 50 jcol = 1, nrhs
341  DO 40 jrow = 1, n
342  j = j + 1
343  rwork( j ) = dble( b( jrow, jcol ) )
344  40 CONTINUE
345  50 CONTINUE
346  CALL dgemm( 'T', 'N', n, nrhs, n, one, rwork( irwu ), n,
347  $ rwork( irwb ), n, zero, rwork( irwrb ), n )
348  j = irwb - 1
349  DO 70 jcol = 1, nrhs
350  DO 60 jrow = 1, n
351  j = j + 1
352  rwork( j ) = dimag( b( jrow, jcol ) )
353  60 CONTINUE
354  70 CONTINUE
355  CALL dgemm( 'T', 'N', n, nrhs, n, one, rwork( irwu ), n,
356  $ rwork( irwb ), n, zero, rwork( irwib ), n )
357  jreal = irwrb - 1
358  jimag = irwib - 1
359  DO 90 jcol = 1, nrhs
360  DO 80 jrow = 1, n
361  jreal = jreal + 1
362  jimag = jimag + 1
363  b( jrow, jcol ) = dcmplx( rwork( jreal ),
364  $ rwork( jimag ) )
365  80 CONTINUE
366  90 CONTINUE
367 *
368  tol = rcnd*abs( d( idamax( n, d, 1 ) ) )
369  DO 100 i = 1, n
370  IF( d( i ).LE.tol ) THEN
371  CALL zlaset( 'A', 1, nrhs, czero, czero, b( i, 1 ), ldb )
372  ELSE
373  CALL zlascl( 'G', 0, 0, d( i ), one, 1, nrhs, b( i, 1 ),
374  $ ldb, info )
375  rank = rank + 1
376  END IF
377  100 CONTINUE
378 *
379 * Since B is complex, the following call to DGEMM is performed
380 * in two steps (real and imaginary parts). That is for V * B
381 * (in the real version of the code V**H is stored in WORK).
382 *
383 * CALL DGEMM( 'T', 'N', N, NRHS, N, ONE, WORK, N, B, LDB, ZERO,
384 * $ WORK( NWORK ), N )
385 *
386  j = irwb - 1
387  DO 120 jcol = 1, nrhs
388  DO 110 jrow = 1, n
389  j = j + 1
390  rwork( j ) = dble( b( jrow, jcol ) )
391  110 CONTINUE
392  120 CONTINUE
393  CALL dgemm( 'T', 'N', n, nrhs, n, one, rwork( irwvt ), n,
394  $ rwork( irwb ), n, zero, rwork( irwrb ), n )
395  j = irwb - 1
396  DO 140 jcol = 1, nrhs
397  DO 130 jrow = 1, n
398  j = j + 1
399  rwork( j ) = dimag( b( jrow, jcol ) )
400  130 CONTINUE
401  140 CONTINUE
402  CALL dgemm( 'T', 'N', n, nrhs, n, one, rwork( irwvt ), n,
403  $ rwork( irwb ), n, zero, rwork( irwib ), n )
404  jreal = irwrb - 1
405  jimag = irwib - 1
406  DO 160 jcol = 1, nrhs
407  DO 150 jrow = 1, n
408  jreal = jreal + 1
409  jimag = jimag + 1
410  b( jrow, jcol ) = dcmplx( rwork( jreal ),
411  $ rwork( jimag ) )
412  150 CONTINUE
413  160 CONTINUE
414 *
415 * Unscale.
416 *
417  CALL dlascl( 'G', 0, 0, one, orgnrm, n, 1, d, n, info )
418  CALL dlasrt( 'D', n, d, info )
419  CALL zlascl( 'G', 0, 0, orgnrm, one, n, nrhs, b, ldb, info )
420 *
421  RETURN
422  END IF
423 *
424 * Book-keeping and setting up some constants.
425 *
426  nlvl = int( log( dble( n ) / dble( smlsiz+1 ) ) / log( two ) ) + 1
427 *
428  smlszp = smlsiz + 1
429 *
430  u = 1
431  vt = 1 + smlsiz*n
432  difl = vt + smlszp*n
433  difr = difl + nlvl*n
434  z = difr + nlvl*n*2
435  c = z + nlvl*n
436  s = c + n
437  poles = s + n
438  givnum = poles + 2*nlvl*n
439  nrwork = givnum + 2*nlvl*n
440  bx = 1
441 *
442  irwrb = nrwork
443  irwib = irwrb + smlsiz*nrhs
444  irwb = irwib + smlsiz*nrhs
445 *
446  sizei = 1 + n
447  k = sizei + n
448  givptr = k + n
449  perm = givptr + n
450  givcol = perm + nlvl*n
451  iwk = givcol + nlvl*n*2
452 *
453  st = 1
454  sqre = 0
455  icmpq1 = 1
456  icmpq2 = 0
457  nsub = 0
458 *
459  DO 170 i = 1, n
460  IF( abs( d( i ) ).LT.eps ) THEN
461  d( i ) = sign( eps, d( i ) )
462  END IF
463  170 CONTINUE
464 *
465  DO 240 i = 1, nm1
466  IF( ( abs( e( i ) ).LT.eps ) .OR. ( i.EQ.nm1 ) ) THEN
467  nsub = nsub + 1
468  iwork( nsub ) = st
469 *
470 * Subproblem found. First determine its size and then
471 * apply divide and conquer on it.
472 *
473  IF( i.LT.nm1 ) THEN
474 *
475 * A subproblem with E(I) small for I < NM1.
476 *
477  nsize = i - st + 1
478  iwork( sizei+nsub-1 ) = nsize
479  ELSE IF( abs( e( i ) ).GE.eps ) THEN
480 *
481 * A subproblem with E(NM1) not too small but I = NM1.
482 *
483  nsize = n - st + 1
484  iwork( sizei+nsub-1 ) = nsize
485  ELSE
486 *
487 * A subproblem with E(NM1) small. This implies an
488 * 1-by-1 subproblem at D(N), which is not solved
489 * explicitly.
490 *
491  nsize = i - st + 1
492  iwork( sizei+nsub-1 ) = nsize
493  nsub = nsub + 1
494  iwork( nsub ) = n
495  iwork( sizei+nsub-1 ) = 1
496  CALL zcopy( nrhs, b( n, 1 ), ldb, work( bx+nm1 ), n )
497  END IF
498  st1 = st - 1
499  IF( nsize.EQ.1 ) THEN
500 *
501 * This is a 1-by-1 subproblem and is not solved
502 * explicitly.
503 *
504  CALL zcopy( nrhs, b( st, 1 ), ldb, work( bx+st1 ), n )
505  ELSE IF( nsize.LE.smlsiz ) THEN
506 *
507 * This is a small subproblem and is solved by DLASDQ.
508 *
509  CALL dlaset( 'A', nsize, nsize, zero, one,
510  $ rwork( vt+st1 ), n )
511  CALL dlaset( 'A', nsize, nsize, zero, one,
512  $ rwork( u+st1 ), n )
513  CALL dlasdq( 'U', 0, nsize, nsize, nsize, 0, d( st ),
514  $ e( st ), rwork( vt+st1 ), n, rwork( u+st1 ),
515  $ n, rwork( nrwork ), 1, rwork( nrwork ),
516  $ info )
517  IF( info.NE.0 ) THEN
518  RETURN
519  END IF
520 *
521 * In the real version, B is passed to DLASDQ and multiplied
522 * internally by Q**H. Here B is complex and that product is
523 * computed below in two steps (real and imaginary parts).
524 *
525  j = irwb - 1
526  DO 190 jcol = 1, nrhs
527  DO 180 jrow = st, st + nsize - 1
528  j = j + 1
529  rwork( j ) = dble( b( jrow, jcol ) )
530  180 CONTINUE
531  190 CONTINUE
532  CALL dgemm( 'T', 'N', nsize, nrhs, nsize, one,
533  $ rwork( u+st1 ), n, rwork( irwb ), nsize,
534  $ zero, rwork( irwrb ), nsize )
535  j = irwb - 1
536  DO 210 jcol = 1, nrhs
537  DO 200 jrow = st, st + nsize - 1
538  j = j + 1
539  rwork( j ) = dimag( b( jrow, jcol ) )
540  200 CONTINUE
541  210 CONTINUE
542  CALL dgemm( 'T', 'N', nsize, nrhs, nsize, one,
543  $ rwork( u+st1 ), n, rwork( irwb ), nsize,
544  $ zero, rwork( irwib ), nsize )
545  jreal = irwrb - 1
546  jimag = irwib - 1
547  DO 230 jcol = 1, nrhs
548  DO 220 jrow = st, st + nsize - 1
549  jreal = jreal + 1
550  jimag = jimag + 1
551  b( jrow, jcol ) = dcmplx( rwork( jreal ),
552  $ rwork( jimag ) )
553  220 CONTINUE
554  230 CONTINUE
555 *
556  CALL zlacpy( 'A', nsize, nrhs, b( st, 1 ), ldb,
557  $ work( bx+st1 ), n )
558  ELSE
559 *
560 * A large problem. Solve it using divide and conquer.
561 *
562  CALL dlasda( icmpq1, smlsiz, nsize, sqre, d( st ),
563  $ e( st ), rwork( u+st1 ), n, rwork( vt+st1 ),
564  $ iwork( k+st1 ), rwork( difl+st1 ),
565  $ rwork( difr+st1 ), rwork( z+st1 ),
566  $ rwork( poles+st1 ), iwork( givptr+st1 ),
567  $ iwork( givcol+st1 ), n, iwork( perm+st1 ),
568  $ rwork( givnum+st1 ), rwork( c+st1 ),
569  $ rwork( s+st1 ), rwork( nrwork ),
570  $ iwork( iwk ), info )
571  IF( info.NE.0 ) THEN
572  RETURN
573  END IF
574  bxst = bx + st1
575  CALL zlalsa( icmpq2, smlsiz, nsize, nrhs, b( st, 1 ),
576  $ ldb, work( bxst ), n, rwork( u+st1 ), n,
577  $ rwork( vt+st1 ), iwork( k+st1 ),
578  $ rwork( difl+st1 ), rwork( difr+st1 ),
579  $ rwork( z+st1 ), rwork( poles+st1 ),
580  $ iwork( givptr+st1 ), iwork( givcol+st1 ), n,
581  $ iwork( perm+st1 ), rwork( givnum+st1 ),
582  $ rwork( c+st1 ), rwork( s+st1 ),
583  $ rwork( nrwork ), iwork( iwk ), info )
584  IF( info.NE.0 ) THEN
585  RETURN
586  END IF
587  END IF
588  st = i + 1
589  END IF
590  240 CONTINUE
591 *
592 * Apply the singular values and treat the tiny ones as zero.
593 *
594  tol = rcnd*abs( d( idamax( n, d, 1 ) ) )
595 *
596  DO 250 i = 1, n
597 *
598 * Some of the elements in D can be negative because 1-by-1
599 * subproblems were not solved explicitly.
600 *
601  IF( abs( d( i ) ).LE.tol ) THEN
602  CALL zlaset( 'A', 1, nrhs, czero, czero, work( bx+i-1 ), n )
603  ELSE
604  rank = rank + 1
605  CALL zlascl( 'G', 0, 0, d( i ), one, 1, nrhs,
606  $ work( bx+i-1 ), n, info )
607  END IF
608  d( i ) = abs( d( i ) )
609  250 CONTINUE
610 *
611 * Now apply back the right singular vectors.
612 *
613  icmpq2 = 1
614  DO 320 i = 1, nsub
615  st = iwork( i )
616  st1 = st - 1
617  nsize = iwork( sizei+i-1 )
618  bxst = bx + st1
619  IF( nsize.EQ.1 ) THEN
620  CALL zcopy( nrhs, work( bxst ), n, b( st, 1 ), ldb )
621  ELSE IF( nsize.LE.smlsiz ) THEN
622 *
623 * Since B and BX are complex, the following call to DGEMM
624 * is performed in two steps (real and imaginary parts).
625 *
626 * CALL DGEMM( 'T', 'N', NSIZE, NRHS, NSIZE, ONE,
627 * $ RWORK( VT+ST1 ), N, RWORK( BXST ), N, ZERO,
628 * $ B( ST, 1 ), LDB )
629 *
630  j = bxst - n - 1
631  jreal = irwb - 1
632  DO 270 jcol = 1, nrhs
633  j = j + n
634  DO 260 jrow = 1, nsize
635  jreal = jreal + 1
636  rwork( jreal ) = dble( work( j+jrow ) )
637  260 CONTINUE
638  270 CONTINUE
639  CALL dgemm( 'T', 'N', nsize, nrhs, nsize, one,
640  $ rwork( vt+st1 ), n, rwork( irwb ), nsize, zero,
641  $ rwork( irwrb ), nsize )
642  j = bxst - n - 1
643  jimag = irwb - 1
644  DO 290 jcol = 1, nrhs
645  j = j + n
646  DO 280 jrow = 1, nsize
647  jimag = jimag + 1
648  rwork( jimag ) = dimag( work( j+jrow ) )
649  280 CONTINUE
650  290 CONTINUE
651  CALL dgemm( 'T', 'N', nsize, nrhs, nsize, one,
652  $ rwork( vt+st1 ), n, rwork( irwb ), nsize, zero,
653  $ rwork( irwib ), nsize )
654  jreal = irwrb - 1
655  jimag = irwib - 1
656  DO 310 jcol = 1, nrhs
657  DO 300 jrow = st, st + nsize - 1
658  jreal = jreal + 1
659  jimag = jimag + 1
660  b( jrow, jcol ) = dcmplx( rwork( jreal ),
661  $ rwork( jimag ) )
662  300 CONTINUE
663  310 CONTINUE
664  ELSE
665  CALL zlalsa( icmpq2, smlsiz, nsize, nrhs, work( bxst ), n,
666  $ b( st, 1 ), ldb, rwork( u+st1 ), n,
667  $ rwork( vt+st1 ), iwork( k+st1 ),
668  $ rwork( difl+st1 ), rwork( difr+st1 ),
669  $ rwork( z+st1 ), rwork( poles+st1 ),
670  $ iwork( givptr+st1 ), iwork( givcol+st1 ), n,
671  $ iwork( perm+st1 ), rwork( givnum+st1 ),
672  $ rwork( c+st1 ), rwork( s+st1 ),
673  $ rwork( nrwork ), iwork( iwk ), info )
674  IF( info.NE.0 ) THEN
675  RETURN
676  END IF
677  END IF
678  320 CONTINUE
679 *
680 * Unscale and sort the singular values.
681 *
682  CALL dlascl( 'G', 0, 0, one, orgnrm, n, 1, d, n, info )
683  CALL dlasrt( 'D', n, d, info )
684  CALL zlascl( 'G', 0, 0, orgnrm, one, n, nrhs, b, ldb, info )
685 *
686  RETURN
687 *
688 * End of ZLALSD
689 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
double precision function dlanst(NORM, N, D, E)
DLANST returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: dlanst.f:100
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: dlascl.f:143
subroutine dlartg(f, g, c, s, r)
DLARTG generates a plane rotation with real cosine and real sine.
Definition: dlartg.f90:113
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: dlaset.f:110
subroutine dlasda(ICOMPQ, SMLSIZ, N, SQRE, D, E, U, LDU, VT, K, DIFL, DIFR, Z, POLES, GIVPTR, GIVCOL, LDGCOL, PERM, GIVNUM, C, S, WORK, IWORK, INFO)
DLASDA computes the singular value decomposition (SVD) of a real upper bidiagonal matrix with diagona...
Definition: dlasda.f:273
subroutine dlasdq(UPLO, SQRE, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, LDU, C, LDC, WORK, INFO)
DLASDQ computes the SVD of a real bidiagonal matrix with diagonal d and off-diagonal e....
Definition: dlasdq.f:211
integer function idamax(N, DX, INCX)
IDAMAX
Definition: idamax.f:71
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dlasrt(ID, N, D, INFO)
DLASRT sorts numbers in increasing or decreasing order.
Definition: dlasrt.f:88
subroutine zdrot(N, ZX, INCX, ZY, INCY, C, S)
ZDROT
Definition: zdrot.f:98
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:81
subroutine zlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: zlascl.f:143
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: zlaset.f:106
subroutine zlalsa(ICOMPQ, SMLSIZ, N, NRHS, B, LDB, BX, LDBX, U, LDU, VT, K, DIFL, DIFR, Z, POLES, GIVPTR, GIVCOL, LDGCOL, PERM, GIVNUM, C, S, RWORK, IWORK, INFO)
ZLALSA computes the SVD of the coefficient matrix in compact form. Used by sgelsd.
Definition: zlalsa.f:267
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:187
Here is the call graph for this function:
Here is the caller graph for this function: