LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine sgelsd ( integer  M,
integer  N,
integer  NRHS,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( ldb, * )  B,
integer  LDB,
real, dimension( * )  S,
real  RCOND,
integer  RANK,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  INFO 
)

SGELSD computes the minimum-norm solution to a linear least squares problem for GE matrices

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

Purpose:
 SGELSD computes the minimum-norm solution to a real linear least
 squares problem:
     minimize 2-norm(| b - A*x |)
 using the singular value decomposition (SVD) of A. A is an M-by-N
 matrix which may be rank-deficient.

 Several right hand side vectors b and solution vectors x can be
 handled in a single call; they are stored as the columns of the
 M-by-NRHS right hand side matrix B and the N-by-NRHS solution
 matrix X.

 The problem is solved in three steps:
 (1) Reduce the coefficient matrix A to bidiagonal form with
     Householder transformations, reducing the original problem
     into a "bidiagonal least squares problem" (BLS)
 (2) Solve the BLS using a divide and conquer approach.
 (3) Apply back all the Householder tranformations to solve
     the original least squares problem.

 The effective rank of A is determined by treating as zero those
 singular values which are less than RCOND times the largest singular
 value.

 The divide and conquer algorithm 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 X-MP, Cray Y-MP, 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]M
          M is INTEGER
          The number of rows of A. M >= 0.
[in]N
          N is INTEGER
          The number of columns of A. N >= 0.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrices B and X. NRHS >= 0.
[in]A
          A is REAL array, dimension (LDA,N)
          On entry, the M-by-N matrix A.
          On exit, A has been destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).
[in,out]B
          B is REAL array, dimension (LDB,NRHS)
          On entry, the M-by-NRHS right hand side matrix B.
          On exit, B is overwritten by the N-by-NRHS solution
          matrix X.  If m >= n and RANK = n, the residual
          sum-of-squares for the solution in the i-th column is given
          by the sum of squares of elements n+1:m in that column.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B. LDB >= max(1,max(M,N)).
[out]S
          S is REAL array, dimension (min(M,N))
          The singular values of A in decreasing order.
          The condition number of A in the 2-norm = S(1)/S(min(m,n)).
[in]RCOND
          RCOND is REAL
          RCOND is used to determine the effective rank of A.
          Singular values S(i) <= RCOND*S(1) are treated as zero.
          If RCOND < 0, machine precision is used instead.
[out]RANK
          RANK is INTEGER
          The effective rank of A, i.e., the number of singular values
          which are greater than RCOND*S(1).
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK. LWORK must be at least 1.
          The exact minimum amount of workspace needed depends on M,
          N and NRHS. As long as LWORK is at least
              12*N + 2*N*SMLSIZ + 8*N*NLVL + N*NRHS + (SMLSIZ+1)**2,
          if M is greater than or equal to N or
              12*M + 2*M*SMLSIZ + 8*M*NLVL + M*NRHS + (SMLSIZ+1)**2,
          if M is less than N, the code will execute correctly.
          SMLSIZ is returned by ILAENV and is equal to the maximum
          size of the subproblems at the bottom of the computation
          tree (usually about 25), and
             NLVL = MAX( 0, INT( LOG_2( MIN( M,N )/(SMLSIZ+1) ) ) + 1 )
          For good performance, LWORK should generally be larger.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the array WORK and the
          minimum size of the array IWORK, and returns these values as
          the first entries of the WORK and IWORK arrays, and no error
          message related to LWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
          LIWORK >= max(1, 3*MINMN*NLVL + 11*MINMN),
          where MINMN = MIN( M,N ).
          On exit, if INFO = 0, IWORK(1) returns the minimum LIWORK.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value.
          > 0:  the algorithm for computing the SVD failed to converge;
                if INFO = i, i off-diagonal elements of an intermediate
                bidiagonal form did not converge to zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Contributors:
Ming Gu and Ren-Cang Li, Computer Science Division, University of California at Berkeley, USA
Osni Marques, LBNL/NERSC, USA

Definition at line 212 of file sgelsd.f.

212 *
213 * -- LAPACK driver routine (version 3.4.0) --
214 * -- LAPACK is a software package provided by Univ. of Tennessee, --
215 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
216 * November 2011
217 *
218 * .. Scalar Arguments ..
219  INTEGER info, lda, ldb, lwork, m, n, nrhs, rank
220  REAL rcond
221 * ..
222 * .. Array Arguments ..
223  INTEGER iwork( * )
224  REAL a( lda, * ), b( ldb, * ), s( * ), work( * )
225 * ..
226 *
227 * =====================================================================
228 *
229 * .. Parameters ..
230  REAL zero, one, two
231  parameter ( zero = 0.0e0, one = 1.0e0, two = 2.0e0 )
232 * ..
233 * .. Local Scalars ..
234  LOGICAL lquery
235  INTEGER iascl, ibscl, ie, il, itau, itaup, itauq,
236  $ ldwork, liwork, maxmn, maxwrk, minmn, minwrk,
237  $ mm, mnthr, nlvl, nwork, smlsiz, wlalsd
238  REAL anrm, bignum, bnrm, eps, sfmin, smlnum
239 * ..
240 * .. External Subroutines ..
241  EXTERNAL sgebrd, sgelqf, sgeqrf, slabad, slacpy, slalsd,
243 * ..
244 * .. External Functions ..
245  INTEGER ilaenv
246  REAL slamch, slange
247  EXTERNAL slamch, slange, ilaenv
248 * ..
249 * .. Intrinsic Functions ..
250  INTRINSIC int, log, max, min, real
251 * ..
252 * .. Executable Statements ..
253 *
254 * Test the input arguments.
255 *
256  info = 0
257  minmn = min( m, n )
258  maxmn = max( m, n )
259  lquery = ( lwork.EQ.-1 )
260  IF( m.LT.0 ) THEN
261  info = -1
262  ELSE IF( n.LT.0 ) THEN
263  info = -2
264  ELSE IF( nrhs.LT.0 ) THEN
265  info = -3
266  ELSE IF( lda.LT.max( 1, m ) ) THEN
267  info = -5
268  ELSE IF( ldb.LT.max( 1, maxmn ) ) THEN
269  info = -7
270  END IF
271 *
272 * Compute workspace.
273 * (Note: Comments in the code beginning "Workspace:" describe the
274 * minimal amount of workspace needed at that point in the code,
275 * as well as the preferred amount for good performance.
276 * NB refers to the optimal block size for the immediately
277 * following subroutine, as returned by ILAENV.)
278 *
279  IF( info.EQ.0 ) THEN
280  minwrk = 1
281  maxwrk = 1
282  liwork = 1
283  IF( minmn.GT.0 ) THEN
284  smlsiz = ilaenv( 9, 'SGELSD', ' ', 0, 0, 0, 0 )
285  mnthr = ilaenv( 6, 'SGELSD', ' ', m, n, nrhs, -1 )
286  nlvl = max( int( log( REAL( MINMN ) / REAL( SMLSIZ + 1 ) ) /
287  $ log( two ) ) + 1, 0 )
288  liwork = 3*minmn*nlvl + 11*minmn
289  mm = m
290  IF( m.GE.n .AND. m.GE.mnthr ) THEN
291 *
292 * Path 1a - overdetermined, with many more rows than
293 * columns.
294 *
295  mm = n
296  maxwrk = max( maxwrk, n + n*ilaenv( 1, 'SGEQRF', ' ', m,
297  $ n, -1, -1 ) )
298  maxwrk = max( maxwrk, n + nrhs*ilaenv( 1, 'SORMQR', 'LT',
299  $ m, nrhs, n, -1 ) )
300  END IF
301  IF( m.GE.n ) THEN
302 *
303 * Path 1 - overdetermined or exactly determined.
304 *
305  maxwrk = max( maxwrk, 3*n + ( mm + n )*ilaenv( 1,
306  $ 'SGEBRD', ' ', mm, n, -1, -1 ) )
307  maxwrk = max( maxwrk, 3*n + nrhs*ilaenv( 1, 'SORMBR',
308  $ 'QLT', mm, nrhs, n, -1 ) )
309  maxwrk = max( maxwrk, 3*n + ( n - 1 )*ilaenv( 1,
310  $ 'SORMBR', 'PLN', n, nrhs, n, -1 ) )
311  wlalsd = 9*n + 2*n*smlsiz + 8*n*nlvl + n*nrhs +
312  $ ( smlsiz + 1 )**2
313  maxwrk = max( maxwrk, 3*n + wlalsd )
314  minwrk = max( 3*n + mm, 3*n + nrhs, 3*n + wlalsd )
315  END IF
316  IF( n.GT.m ) THEN
317  wlalsd = 9*m + 2*m*smlsiz + 8*m*nlvl + m*nrhs +
318  $ ( smlsiz + 1 )**2
319  IF( n.GE.mnthr ) THEN
320 *
321 * Path 2a - underdetermined, with many more columns
322 * than rows.
323 *
324  maxwrk = m + m*ilaenv( 1, 'SGELQF', ' ', m, n, -1,
325  $ -1 )
326  maxwrk = max( maxwrk, m*m + 4*m + 2*m*ilaenv( 1,
327  $ 'SGEBRD', ' ', m, m, -1, -1 ) )
328  maxwrk = max( maxwrk, m*m + 4*m + nrhs*ilaenv( 1,
329  $ 'SORMBR', 'QLT', m, nrhs, m, -1 ) )
330  maxwrk = max( maxwrk, m*m + 4*m + ( m - 1 )*ilaenv( 1,
331  $ 'SORMBR', 'PLN', m, nrhs, m, -1 ) )
332  IF( nrhs.GT.1 ) THEN
333  maxwrk = max( maxwrk, m*m + m + m*nrhs )
334  ELSE
335  maxwrk = max( maxwrk, m*m + 2*m )
336  END IF
337  maxwrk = max( maxwrk, m + nrhs*ilaenv( 1, 'SORMLQ',
338  $ 'LT', n, nrhs, m, -1 ) )
339  maxwrk = max( maxwrk, m*m + 4*m + wlalsd )
340 ! XXX: Ensure the Path 2a case below is triggered. The workspace
341 ! calculation should use queries for all routines eventually.
342  maxwrk = max( maxwrk,
343  $ 4*m+m*m+max( m, 2*m-4, nrhs, n-3*m ) )
344  ELSE
345 *
346 * Path 2 - remaining underdetermined cases.
347 *
348  maxwrk = 3*m + ( n + m )*ilaenv( 1, 'SGEBRD', ' ', m,
349  $ n, -1, -1 )
350  maxwrk = max( maxwrk, 3*m + nrhs*ilaenv( 1, 'SORMBR',
351  $ 'QLT', m, nrhs, n, -1 ) )
352  maxwrk = max( maxwrk, 3*m + m*ilaenv( 1, 'SORMBR',
353  $ 'PLN', n, nrhs, m, -1 ) )
354  maxwrk = max( maxwrk, 3*m + wlalsd )
355  END IF
356  minwrk = max( 3*m + nrhs, 3*m + m, 3*m + wlalsd )
357  END IF
358  END IF
359  minwrk = min( minwrk, maxwrk )
360  work( 1 ) = maxwrk
361  iwork( 1 ) = liwork
362 *
363  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
364  info = -12
365  END IF
366  END IF
367 *
368  IF( info.NE.0 ) THEN
369  CALL xerbla( 'SGELSD', -info )
370  RETURN
371  ELSE IF( lquery ) THEN
372  RETURN
373  END IF
374 *
375 * Quick return if possible.
376 *
377  IF( m.EQ.0 .OR. n.EQ.0 ) THEN
378  rank = 0
379  RETURN
380  END IF
381 *
382 * Get machine parameters.
383 *
384  eps = slamch( 'P' )
385  sfmin = slamch( 'S' )
386  smlnum = sfmin / eps
387  bignum = one / smlnum
388  CALL slabad( smlnum, bignum )
389 *
390 * Scale A if max entry outside range [SMLNUM,BIGNUM].
391 *
392  anrm = slange( 'M', m, n, a, lda, work )
393  iascl = 0
394  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
395 *
396 * Scale matrix norm up to SMLNUM.
397 *
398  CALL slascl( 'G', 0, 0, anrm, smlnum, m, n, a, lda, info )
399  iascl = 1
400  ELSE IF( anrm.GT.bignum ) THEN
401 *
402 * Scale matrix norm down to BIGNUM.
403 *
404  CALL slascl( 'G', 0, 0, anrm, bignum, m, n, a, lda, info )
405  iascl = 2
406  ELSE IF( anrm.EQ.zero ) THEN
407 *
408 * Matrix all zero. Return zero solution.
409 *
410  CALL slaset( 'F', max( m, n ), nrhs, zero, zero, b, ldb )
411  CALL slaset( 'F', minmn, 1, zero, zero, s, 1 )
412  rank = 0
413  GO TO 10
414  END IF
415 *
416 * Scale B if max entry outside range [SMLNUM,BIGNUM].
417 *
418  bnrm = slange( 'M', m, nrhs, b, ldb, work )
419  ibscl = 0
420  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
421 *
422 * Scale matrix norm up to SMLNUM.
423 *
424  CALL slascl( 'G', 0, 0, bnrm, smlnum, m, nrhs, b, ldb, info )
425  ibscl = 1
426  ELSE IF( bnrm.GT.bignum ) THEN
427 *
428 * Scale matrix norm down to BIGNUM.
429 *
430  CALL slascl( 'G', 0, 0, bnrm, bignum, m, nrhs, b, ldb, info )
431  ibscl = 2
432  END IF
433 *
434 * If M < N make sure certain entries of B are zero.
435 *
436  IF( m.LT.n )
437  $ CALL slaset( 'F', n-m, nrhs, zero, zero, b( m+1, 1 ), ldb )
438 *
439 * Overdetermined case.
440 *
441  IF( m.GE.n ) THEN
442 *
443 * Path 1 - overdetermined or exactly determined.
444 *
445  mm = m
446  IF( m.GE.mnthr ) THEN
447 *
448 * Path 1a - overdetermined, with many more rows than columns.
449 *
450  mm = n
451  itau = 1
452  nwork = itau + n
453 *
454 * Compute A=Q*R.
455 * (Workspace: need 2*N, prefer N+N*NB)
456 *
457  CALL sgeqrf( m, n, a, lda, work( itau ), work( nwork ),
458  $ lwork-nwork+1, info )
459 *
460 * Multiply B by transpose(Q).
461 * (Workspace: need N+NRHS, prefer N+NRHS*NB)
462 *
463  CALL sormqr( 'L', 'T', m, nrhs, n, a, lda, work( itau ), b,
464  $ ldb, work( nwork ), lwork-nwork+1, info )
465 *
466 * Zero out below R.
467 *
468  IF( n.GT.1 ) THEN
469  CALL slaset( 'L', n-1, n-1, zero, zero, a( 2, 1 ), lda )
470  END IF
471  END IF
472 *
473  ie = 1
474  itauq = ie + n
475  itaup = itauq + n
476  nwork = itaup + n
477 *
478 * Bidiagonalize R in A.
479 * (Workspace: need 3*N+MM, prefer 3*N+(MM+N)*NB)
480 *
481  CALL sgebrd( mm, n, a, lda, s, work( ie ), work( itauq ),
482  $ work( itaup ), work( nwork ), lwork-nwork+1,
483  $ info )
484 *
485 * Multiply B by transpose of left bidiagonalizing vectors of R.
486 * (Workspace: need 3*N+NRHS, prefer 3*N+NRHS*NB)
487 *
488  CALL sormbr( 'Q', 'L', 'T', mm, nrhs, n, a, lda, work( itauq ),
489  $ b, ldb, work( nwork ), lwork-nwork+1, info )
490 *
491 * Solve the bidiagonal least squares problem.
492 *
493  CALL slalsd( 'U', smlsiz, n, nrhs, s, work( ie ), b, ldb,
494  $ rcond, rank, work( nwork ), iwork, info )
495  IF( info.NE.0 ) THEN
496  GO TO 10
497  END IF
498 *
499 * Multiply B by right bidiagonalizing vectors of R.
500 *
501  CALL sormbr( 'P', 'L', 'N', n, nrhs, n, a, lda, work( itaup ),
502  $ b, ldb, work( nwork ), lwork-nwork+1, info )
503 *
504  ELSE IF( n.GE.mnthr .AND. lwork.GE.4*m+m*m+
505  $ max( m, 2*m-4, nrhs, n-3*m, wlalsd ) ) THEN
506 *
507 * Path 2a - underdetermined, with many more columns than rows
508 * and sufficient workspace for an efficient algorithm.
509 *
510  ldwork = m
511  IF( lwork.GE.max( 4*m+m*lda+max( m, 2*m-4, nrhs, n-3*m ),
512  $ m*lda+m+m*nrhs, 4*m+m*lda+wlalsd ) )ldwork = lda
513  itau = 1
514  nwork = m + 1
515 *
516 * Compute A=L*Q.
517 * (Workspace: need 2*M, prefer M+M*NB)
518 *
519  CALL sgelqf( m, n, a, lda, work( itau ), work( nwork ),
520  $ lwork-nwork+1, info )
521  il = nwork
522 *
523 * Copy L to WORK(IL), zeroing out above its diagonal.
524 *
525  CALL slacpy( 'L', m, m, a, lda, work( il ), ldwork )
526  CALL slaset( 'U', m-1, m-1, zero, zero, work( il+ldwork ),
527  $ ldwork )
528  ie = il + ldwork*m
529  itauq = ie + m
530  itaup = itauq + m
531  nwork = itaup + m
532 *
533 * Bidiagonalize L in WORK(IL).
534 * (Workspace: need M*M+5*M, prefer M*M+4*M+2*M*NB)
535 *
536  CALL sgebrd( m, m, work( il ), ldwork, s, work( ie ),
537  $ work( itauq ), work( itaup ), work( nwork ),
538  $ lwork-nwork+1, info )
539 *
540 * Multiply B by transpose of left bidiagonalizing vectors of L.
541 * (Workspace: need M*M+4*M+NRHS, prefer M*M+4*M+NRHS*NB)
542 *
543  CALL sormbr( 'Q', 'L', 'T', m, nrhs, m, work( il ), ldwork,
544  $ work( itauq ), b, ldb, work( nwork ),
545  $ lwork-nwork+1, info )
546 *
547 * Solve the bidiagonal least squares problem.
548 *
549  CALL slalsd( 'U', smlsiz, m, nrhs, s, work( ie ), b, ldb,
550  $ rcond, rank, work( nwork ), iwork, info )
551  IF( info.NE.0 ) THEN
552  GO TO 10
553  END IF
554 *
555 * Multiply B by right bidiagonalizing vectors of L.
556 *
557  CALL sormbr( 'P', 'L', 'N', m, nrhs, m, work( il ), ldwork,
558  $ work( itaup ), b, ldb, work( nwork ),
559  $ lwork-nwork+1, info )
560 *
561 * Zero out below first M rows of B.
562 *
563  CALL slaset( 'F', n-m, nrhs, zero, zero, b( m+1, 1 ), ldb )
564  nwork = itau + m
565 *
566 * Multiply transpose(Q) by B.
567 * (Workspace: need M+NRHS, prefer M+NRHS*NB)
568 *
569  CALL sormlq( 'L', 'T', n, nrhs, m, a, lda, work( itau ), b,
570  $ ldb, work( nwork ), lwork-nwork+1, info )
571 *
572  ELSE
573 *
574 * Path 2 - remaining underdetermined cases.
575 *
576  ie = 1
577  itauq = ie + m
578  itaup = itauq + m
579  nwork = itaup + m
580 *
581 * Bidiagonalize A.
582 * (Workspace: need 3*M+N, prefer 3*M+(M+N)*NB)
583 *
584  CALL sgebrd( m, n, a, lda, s, work( ie ), work( itauq ),
585  $ work( itaup ), work( nwork ), lwork-nwork+1,
586  $ info )
587 *
588 * Multiply B by transpose of left bidiagonalizing vectors.
589 * (Workspace: need 3*M+NRHS, prefer 3*M+NRHS*NB)
590 *
591  CALL sormbr( 'Q', 'L', 'T', m, nrhs, n, a, lda, work( itauq ),
592  $ b, ldb, work( nwork ), lwork-nwork+1, info )
593 *
594 * Solve the bidiagonal least squares problem.
595 *
596  CALL slalsd( 'L', smlsiz, m, nrhs, s, work( ie ), b, ldb,
597  $ rcond, rank, work( nwork ), iwork, info )
598  IF( info.NE.0 ) THEN
599  GO TO 10
600  END IF
601 *
602 * Multiply B by right bidiagonalizing vectors of A.
603 *
604  CALL sormbr( 'P', 'L', 'N', n, nrhs, m, a, lda, work( itaup ),
605  $ b, ldb, work( nwork ), lwork-nwork+1, info )
606 *
607  END IF
608 *
609 * Undo scaling.
610 *
611  IF( iascl.EQ.1 ) THEN
612  CALL slascl( 'G', 0, 0, anrm, smlnum, n, nrhs, b, ldb, info )
613  CALL slascl( 'G', 0, 0, smlnum, anrm, minmn, 1, s, minmn,
614  $ info )
615  ELSE IF( iascl.EQ.2 ) THEN
616  CALL slascl( 'G', 0, 0, anrm, bignum, n, nrhs, b, ldb, info )
617  CALL slascl( 'G', 0, 0, bignum, anrm, minmn, 1, s, minmn,
618  $ info )
619  END IF
620  IF( ibscl.EQ.1 ) THEN
621  CALL slascl( 'G', 0, 0, smlnum, bnrm, n, nrhs, b, ldb, info )
622  ELSE IF( ibscl.EQ.2 ) THEN
623  CALL slascl( 'G', 0, 0, bignum, bnrm, n, nrhs, b, ldb, info )
624  END IF
625 *
626  10 CONTINUE
627  work( 1 ) = maxwrk
628  iwork( 1 ) = liwork
629  RETURN
630 *
631 * End of SGELSD
632 *
subroutine sormqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMQR
Definition: sormqr.f:170
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.
Definition: slascl.f:145
subroutine slabad(SMALL, LARGE)
SLABAD
Definition: slabad.f:76
subroutine slalsd(UPLO, SMLSIZ, N, NRHS, D, E, B, LDB, RCOND, RANK, WORK, IWORK, INFO)
SLALSD uses the singular value decomposition of A to solve the least squares problem.
Definition: slalsd.f:181
subroutine sgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
SGEQRF
Definition: sgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: slaset.f:112
subroutine sormlq(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMLQ
Definition: sormlq.f:170
subroutine sormbr(VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMBR
Definition: sormbr.f:198
real function slange(NORM, M, N, A, LDA, WORK)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: slange.f:116
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine sgebrd(M, N, A, LDA, D, E, TAUQ, TAUP, WORK, LWORK, INFO)
SGEBRD
Definition: sgebrd.f:207
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine sgelqf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
SGELQF
Definition: sgelqf.f:137

Here is the call graph for this function:

Here is the caller graph for this function: