LAPACK  3.8.0
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.
Date
June 2017
Contributors:
Ming Gu and Ren-Cang Li, Computer Science Division, University of California at Berkeley, USA
Osni Marques, LBNL/NERSC, USA

Definition at line 189 of file zlalsd.f.

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