LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ ssytrf_aa_2stage()

subroutine ssytrf_aa_2stage ( character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  TB,
integer  LTB,
integer, dimension( * )  IPIV,
integer, dimension( * )  IPIV2,
real, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

SSYTRF_AA_2STAGE

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

Purpose:
 SSYTRF_AA_2STAGE computes the factorization of a real symmetric matrix A
 using the Aasen's algorithm.  The form of the factorization is

    A = U*T*U**T  or  A = L*T*L**T

 where U (or L) is a product of permutation and unit upper (lower)
 triangular matrices, and T is a symmetric band matrix with the
 bandwidth of NB (NB is internally selected and stored in TB( 1 ), and T is 
 LU factorized with partial pivoting).

 This is the blocked version of the algorithm, calling Level 3 BLAS.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA,N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
          N-by-N upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading N-by-N lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.

          On exit, L is stored below (or above) the subdiaonal blocks,
          when UPLO  is 'L' (or 'U').
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]TB
          TB is REAL array, dimension (LTB)
          On exit, details of the LU factorization of the band matrix.
[in]LTB
          The size of the array TB. LTB >= 4*N, internally
          used to select NB such that LTB >= (3*NB+1)*N.

          If LTB = -1, then a workspace query is assumed; the
          routine only calculates the optimal size of LTB, 
          returns this value as the first entry of TB, and
          no error message related to LTB is issued by XERBLA.
[out]IPIV
          IPIV is INTEGER array, dimension (N)
          On exit, it contains the details of the interchanges, i.e.,
          the row and column k of A were interchanged with the
          row and column IPIV(k).
[out]IPIV2
          IPIV is INTEGER array, dimension (N)
          On exit, it contains the details of the interchanges, i.e.,
          the row and column k of T were interchanged with the
          row and column IPIV(k).
[out]WORK
          WORK is REAL workspace of size LWORK
[in]LWORK
          The size of WORK. LWORK >= N, internally used to select NB
          such that LWORK >= N*NB.

          If LWORK = -1, then a workspace query is assumed; the
          routine only calculates the optimal size of the WORK array,
          returns this value as the first entry of the WORK array, and
          no error message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value.
          > 0:  if INFO = i, band LU factorization failed on i-th column
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2017

Definition at line 160 of file ssytrf_aa_2stage.f.

160 *
161 * -- LAPACK computational routine (version 3.8.0) --
162 * -- LAPACK is a software package provided by Univ. of Tennessee, --
163 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
164 * November 2017
165 *
166  IMPLICIT NONE
167 *
168 * .. Scalar Arguments ..
169  CHARACTER uplo
170  INTEGER n, lda, ltb, lwork, info
171 * ..
172 * .. Array Arguments ..
173  INTEGER ipiv( * ), ipiv2( * )
174  REAL a( lda, * ), tb( * ), work( * )
175 * ..
176 *
177 * =====================================================================
178 * .. Parameters ..
179  REAL zero, one
180  parameter( zero = 0.0e+0, one = 1.0e+0 )
181 *
182 * .. Local Scalars ..
183  LOGICAL upper, tquery, wquery
184  INTEGER i, j, k, i1, i2, td
185  INTEGER ldtb, nb, kb, jb, nt, iinfo
186  REAL piv
187 * ..
188 * .. External Functions ..
189  LOGICAL lsame
190  INTEGER ilaenv
191  EXTERNAL lsame, ilaenv
192 * ..
193 * .. External Subroutines ..
194  EXTERNAL xerbla, scopy, slacgv, slacpy,
195  $ slaset, sgbtrf, sgemm, sgetrf,
196  $ ssygst, sswap, strsm
197 * ..
198 * .. Intrinsic Functions ..
199  INTRINSIC min, max
200 * ..
201 * .. Executable Statements ..
202 *
203 * Test the input parameters.
204 *
205  info = 0
206  upper = lsame( uplo, 'U' )
207  wquery = ( lwork.EQ.-1 )
208  tquery = ( ltb.EQ.-1 )
209  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
210  info = -1
211  ELSE IF( n.LT.0 ) THEN
212  info = -2
213  ELSE IF( lda.LT.max( 1, n ) ) THEN
214  info = -4
215  ELSE IF ( ltb .LT. 4*n .AND. .NOT.tquery ) THEN
216  info = -6
217  ELSE IF ( lwork .LT. n .AND. .NOT.wquery ) THEN
218  info = -10
219  END IF
220 *
221  IF( info.NE.0 ) THEN
222  CALL xerbla( 'SSYTRF_AA_2STAGE', -info )
223  RETURN
224  END IF
225 *
226 * Answer the query
227 *
228  nb = ilaenv( 1, 'SSYTRF_AA_2STAGE', uplo, n, -1, -1, -1 )
229  IF( info.EQ.0 ) THEN
230  IF( tquery ) THEN
231  tb( 1 ) = (3*nb+1)*n
232  END IF
233  IF( wquery ) THEN
234  work( 1 ) = n*nb
235  END IF
236  END IF
237  IF( tquery .OR. wquery ) THEN
238  RETURN
239  END IF
240 *
241 * Quick return
242 *
243  IF ( n.EQ.0 ) THEN
244  RETURN
245  ENDIF
246 *
247 * Determine the number of the block size
248 *
249  ldtb = ltb/n
250  IF( ldtb .LT. 3*nb+1 ) THEN
251  nb = (ldtb-1)/3
252  END IF
253  IF( lwork .LT. nb*n ) THEN
254  nb = lwork/n
255  END IF
256 *
257 * Determine the number of the block columns
258 *
259  nt = (n+nb-1)/nb
260  td = 2*nb
261  kb = min(nb, n)
262 *
263 * Initialize vectors/matrices
264 *
265  DO j = 1, kb
266  ipiv( j ) = j
267  END DO
268 *
269 * Save NB
270 *
271  tb( 1 ) = nb
272 *
273  IF( upper ) THEN
274 *
275 * .....................................................
276 * Factorize A as L*D*L**T using the upper triangle of A
277 * .....................................................
278 *
279  DO j = 0, nt-1
280 *
281 * Generate Jth column of W and H
282 *
283  kb = min(nb, n-j*nb)
284  DO i = 1, j-1
285  IF( i.EQ.1 ) THEN
286 * H(I,J) = T(I,I)*U(I,J) + T(I+1,I)*U(I+1,J)
287  IF( i .EQ. (j-1) ) THEN
288  jb = nb+kb
289  ELSE
290  jb = 2*nb
291  END IF
292  CALL sgemm( 'NoTranspose', 'NoTranspose',
293  $ nb, kb, jb,
294  $ one, tb( td+1 + (i*nb)*ldtb ), ldtb-1,
295  $ a( (i-1)*nb+1, j*nb+1 ), lda,
296  $ zero, work( i*nb+1 ), n )
297  ELSE
298 * H(I,J) = T(I,I-1)*U(I-1,J) + T(I,I)*U(I,J) + T(I,I+1)*U(I+1,J)
299  IF( i .EQ. j-1) THEN
300  jb = 2*nb+kb
301  ELSE
302  jb = 3*nb
303  END IF
304  CALL sgemm( 'NoTranspose', 'NoTranspose',
305  $ nb, kb, jb,
306  $ one, tb( td+nb+1 + ((i-1)*nb)*ldtb ),
307  $ ldtb-1,
308  $ a( (i-2)*nb+1, j*nb+1 ), lda,
309  $ zero, work( i*nb+1 ), n )
310  END IF
311  END DO
312 *
313 * Compute T(J,J)
314 *
315  CALL slacpy( 'Upper', kb, kb, a( j*nb+1, j*nb+1 ), lda,
316  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
317  IF( j.GT.1 ) THEN
318 * T(J,J) = U(1:J,J)'*H(1:J)
319  CALL sgemm( 'Transpose', 'NoTranspose',
320  $ kb, kb, (j-1)*nb,
321  $ -one, a( 1, j*nb+1 ), lda,
322  $ work( nb+1 ), n,
323  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
324 * T(J,J) += U(J,J)'*T(J,J-1)*U(J-1,J)
325  CALL sgemm( 'Transpose', 'NoTranspose',
326  $ kb, nb, kb,
327  $ one, a( (j-1)*nb+1, j*nb+1 ), lda,
328  $ tb( td+nb+1 + ((j-1)*nb)*ldtb ), ldtb-1,
329  $ zero, work( 1 ), n )
330  CALL sgemm( 'NoTranspose', 'NoTranspose',
331  $ kb, kb, nb,
332  $ -one, work( 1 ), n,
333  $ a( (j-2)*nb+1, j*nb+1 ), lda,
334  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
335  END IF
336  IF( j.GT.0 ) THEN
337  CALL ssygst( 1, 'Upper', kb,
338  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1,
339  $ a( (j-1)*nb+1, j*nb+1 ), lda, iinfo )
340  END IF
341 *
342 * Expand T(J,J) into full format
343 *
344  DO i = 1, kb
345  DO k = i+1, kb
346  tb( td+(k-i)+1 + (j*nb+i-1)*ldtb )
347  $ = tb( td-(k-(i+1)) + (j*nb+k-1)*ldtb )
348  END DO
349  END DO
350 *
351  IF( j.LT.nt-1 ) THEN
352  IF( j.GT.0 ) THEN
353 *
354 * Compute H(J,J)
355 *
356  IF( j.EQ.1 ) THEN
357  CALL sgemm( 'NoTranspose', 'NoTranspose',
358  $ kb, kb, kb,
359  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1,
360  $ a( (j-1)*nb+1, j*nb+1 ), lda,
361  $ zero, work( j*nb+1 ), n )
362  ELSE
363  CALL sgemm( 'NoTranspose', 'NoTranspose',
364  $ kb, kb, nb+kb,
365  $ one, tb( td+nb+1 + ((j-1)*nb)*ldtb ),
366  $ ldtb-1,
367  $ a( (j-2)*nb+1, j*nb+1 ), lda,
368  $ zero, work( j*nb+1 ), n )
369  END IF
370 *
371 * Update with the previous column
372 *
373  CALL sgemm( 'Transpose', 'NoTranspose',
374  $ nb, n-(j+1)*nb, j*nb,
375  $ -one, work( nb+1 ), n,
376  $ a( 1, (j+1)*nb+1 ), lda,
377  $ one, a( j*nb+1, (j+1)*nb+1 ), lda )
378  END IF
379 *
380 * Copy panel to workspace to call SGETRF
381 *
382  DO k = 1, nb
383  CALL scopy( n-(j+1)*nb,
384  $ a( j*nb+k, (j+1)*nb+1 ), lda,
385  $ work( 1+(k-1)*n ), 1 )
386  END DO
387 *
388 * Factorize panel
389 *
390  CALL sgetrf( n-(j+1)*nb, nb,
391  $ work, n,
392  $ ipiv( (j+1)*nb+1 ), iinfo )
393 c IF (IINFO.NE.0 .AND. INFO.EQ.0) THEN
394 c INFO = IINFO+(J+1)*NB
395 c END IF
396 *
397 * Copy panel back
398 *
399  DO k = 1, nb
400  CALL scopy( n-(j+1)*nb,
401  $ work( 1+(k-1)*n ), 1,
402  $ a( j*nb+k, (j+1)*nb+1 ), lda )
403  END DO
404 *
405 * Compute T(J+1, J), zero out for GEMM update
406 *
407  kb = min(nb, n-(j+1)*nb)
408  CALL slaset( 'Full', kb, nb, zero, zero,
409  $ tb( td+nb+1 + (j*nb)*ldtb), ldtb-1 )
410  CALL slacpy( 'Upper', kb, nb,
411  $ work, n,
412  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
413  IF( j.GT.0 ) THEN
414  CALL strsm( 'R', 'U', 'N', 'U', kb, nb, one,
415  $ a( (j-1)*nb+1, j*nb+1 ), lda,
416  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
417  END IF
418 *
419 * Copy T(J,J+1) into T(J+1, J), both upper/lower for GEMM
420 * updates
421 *
422  DO k = 1, nb
423  DO i = 1, kb
424  tb( td-nb+k-i+1 + (j*nb+nb+i-1)*ldtb )
425  $ = tb( td+nb+i-k+1 + (j*nb+k-1)*ldtb )
426  END DO
427  END DO
428  CALL slaset( 'Lower', kb, nb, zero, one,
429  $ a( j*nb+1, (j+1)*nb+1), lda )
430 *
431 * Apply pivots to trailing submatrix of A
432 *
433  DO k = 1, kb
434 * > Adjust ipiv
435  ipiv( (j+1)*nb+k ) = ipiv( (j+1)*nb+k ) + (j+1)*nb
436 *
437  i1 = (j+1)*nb+k
438  i2 = ipiv( (j+1)*nb+k )
439  IF( i1.NE.i2 ) THEN
440 * > Apply pivots to previous columns of L
441  CALL sswap( k-1, a( (j+1)*nb+1, i1 ), 1,
442  $ a( (j+1)*nb+1, i2 ), 1 )
443 * > Swap A(I1+1:M, I1) with A(I2, I1+1:M)
444  CALL sswap( i2-i1-1, a( i1, i1+1 ), lda,
445  $ a( i1+1, i2 ), 1 )
446 * > Swap A(I2+1:M, I1) with A(I2+1:M, I2)
447  CALL sswap( n-i2, a( i1, i2+1 ), lda,
448  $ a( i2, i2+1 ), lda )
449 * > Swap A(I1, I1) with A(I2, I2)
450  piv = a( i1, i1 )
451  a( i1, i1 ) = a( i2, i2 )
452  a( i2, i2 ) = piv
453 * > Apply pivots to previous columns of L
454  IF( j.GT.0 ) THEN
455  CALL sswap( j*nb, a( 1, i1 ), 1,
456  $ a( 1, i2 ), 1 )
457  END IF
458  ENDIF
459  END DO
460  END IF
461  END DO
462  ELSE
463 *
464 * .....................................................
465 * Factorize A as L*D*L**T using the lower triangle of A
466 * .....................................................
467 *
468  DO j = 0, nt-1
469 *
470 * Generate Jth column of W and H
471 *
472  kb = min(nb, n-j*nb)
473  DO i = 1, j-1
474  IF( i.EQ.1 ) THEN
475 * H(I,J) = T(I,I)*L(J,I)' + T(I+1,I)'*L(J,I+1)'
476  IF( i .EQ. (j-1) ) THEN
477  jb = nb+kb
478  ELSE
479  jb = 2*nb
480  END IF
481  CALL sgemm( 'NoTranspose', 'Transpose',
482  $ nb, kb, jb,
483  $ one, tb( td+1 + (i*nb)*ldtb ), ldtb-1,
484  $ a( j*nb+1, (i-1)*nb+1 ), lda,
485  $ zero, work( i*nb+1 ), n )
486  ELSE
487 * H(I,J) = T(I,I-1)*L(J,I-1)' + T(I,I)*L(J,I)' + T(I,I+1)*L(J,I+1)'
488  IF( i .EQ. j-1) THEN
489  jb = 2*nb+kb
490  ELSE
491  jb = 3*nb
492  END IF
493  CALL sgemm( 'NoTranspose', 'Transpose',
494  $ nb, kb, jb,
495  $ one, tb( td+nb+1 + ((i-1)*nb)*ldtb ),
496  $ ldtb-1,
497  $ a( j*nb+1, (i-2)*nb+1 ), lda,
498  $ zero, work( i*nb+1 ), n )
499  END IF
500  END DO
501 *
502 * Compute T(J,J)
503 *
504  CALL slacpy( 'Lower', kb, kb, a( j*nb+1, j*nb+1 ), lda,
505  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
506  IF( j.GT.1 ) THEN
507 * T(J,J) = L(J,1:J)*H(1:J)
508  CALL sgemm( 'NoTranspose', 'NoTranspose',
509  $ kb, kb, (j-1)*nb,
510  $ -one, a( j*nb+1, 1 ), lda,
511  $ work( nb+1 ), n,
512  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
513 * T(J,J) += L(J,J)*T(J,J-1)*L(J,J-1)'
514  CALL sgemm( 'NoTranspose', 'NoTranspose',
515  $ kb, nb, kb,
516  $ one, a( j*nb+1, (j-1)*nb+1 ), lda,
517  $ tb( td+nb+1 + ((j-1)*nb)*ldtb ), ldtb-1,
518  $ zero, work( 1 ), n )
519  CALL sgemm( 'NoTranspose', 'Transpose',
520  $ kb, kb, nb,
521  $ -one, work( 1 ), n,
522  $ a( j*nb+1, (j-2)*nb+1 ), lda,
523  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
524  END IF
525  IF( j.GT.0 ) THEN
526  CALL ssygst( 1, 'Lower', kb,
527  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1,
528  $ a( j*nb+1, (j-1)*nb+1 ), lda, iinfo )
529  END IF
530 *
531 * Expand T(J,J) into full format
532 *
533  DO i = 1, kb
534  DO k = i+1, kb
535  tb( td-(k-(i+1)) + (j*nb+k-1)*ldtb )
536  $ = tb( td+(k-i)+1 + (j*nb+i-1)*ldtb )
537  END DO
538  END DO
539 *
540  IF( j.LT.nt-1 ) THEN
541  IF( j.GT.0 ) THEN
542 *
543 * Compute H(J,J)
544 *
545  IF( j.EQ.1 ) THEN
546  CALL sgemm( 'NoTranspose', 'Transpose',
547  $ kb, kb, kb,
548  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1,
549  $ a( j*nb+1, (j-1)*nb+1 ), lda,
550  $ zero, work( j*nb+1 ), n )
551  ELSE
552  CALL sgemm( 'NoTranspose', 'Transpose',
553  $ kb, kb, nb+kb,
554  $ one, tb( td+nb+1 + ((j-1)*nb)*ldtb ),
555  $ ldtb-1,
556  $ a( j*nb+1, (j-2)*nb+1 ), lda,
557  $ zero, work( j*nb+1 ), n )
558  END IF
559 *
560 * Update with the previous column
561 *
562  CALL sgemm( 'NoTranspose', 'NoTranspose',
563  $ n-(j+1)*nb, nb, j*nb,
564  $ -one, a( (j+1)*nb+1, 1 ), lda,
565  $ work( nb+1 ), n,
566  $ one, a( (j+1)*nb+1, j*nb+1 ), lda )
567  END IF
568 *
569 * Factorize panel
570 *
571  CALL sgetrf( n-(j+1)*nb, nb,
572  $ a( (j+1)*nb+1, j*nb+1 ), lda,
573  $ ipiv( (j+1)*nb+1 ), iinfo )
574 c IF (IINFO.NE.0 .AND. INFO.EQ.0) THEN
575 c INFO = IINFO+(J+1)*NB
576 c END IF
577 *
578 * Compute T(J+1, J), zero out for GEMM update
579 *
580  kb = min(nb, n-(j+1)*nb)
581  CALL slaset( 'Full', kb, nb, zero, zero,
582  $ tb( td+nb+1 + (j*nb)*ldtb), ldtb-1 )
583  CALL slacpy( 'Upper', kb, nb,
584  $ a( (j+1)*nb+1, j*nb+1 ), lda,
585  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
586  IF( j.GT.0 ) THEN
587  CALL strsm( 'R', 'L', 'T', 'U', kb, nb, one,
588  $ a( j*nb+1, (j-1)*nb+1 ), lda,
589  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
590  END IF
591 *
592 * Copy T(J+1,J) into T(J, J+1), both upper/lower for GEMM
593 * updates
594 *
595  DO k = 1, nb
596  DO i = 1, kb
597  tb( td-nb+k-i+1 + (j*nb+nb+i-1)*ldtb ) =
598  $ tb( td+nb+i-k+1 + (j*nb+k-1)*ldtb )
599  END DO
600  END DO
601  CALL slaset( 'Upper', kb, nb, zero, one,
602  $ a( (j+1)*nb+1, j*nb+1), lda )
603 *
604 * Apply pivots to trailing submatrix of A
605 *
606  DO k = 1, kb
607 * > Adjust ipiv
608  ipiv( (j+1)*nb+k ) = ipiv( (j+1)*nb+k ) + (j+1)*nb
609 *
610  i1 = (j+1)*nb+k
611  i2 = ipiv( (j+1)*nb+k )
612  IF( i1.NE.i2 ) THEN
613 * > Apply pivots to previous columns of L
614  CALL sswap( k-1, a( i1, (j+1)*nb+1 ), lda,
615  $ a( i2, (j+1)*nb+1 ), lda )
616 * > Swap A(I1+1:M, I1) with A(I2, I1+1:M)
617  CALL sswap( i2-i1-1, a( i1+1, i1 ), 1,
618  $ a( i2, i1+1 ), lda )
619 * > Swap A(I2+1:M, I1) with A(I2+1:M, I2)
620  CALL sswap( n-i2, a( i2+1, i1 ), 1,
621  $ a( i2+1, i2 ), 1 )
622 * > Swap A(I1, I1) with A(I2, I2)
623  piv = a( i1, i1 )
624  a( i1, i1 ) = a( i2, i2 )
625  a( i2, i2 ) = piv
626 * > Apply pivots to previous columns of L
627  IF( j.GT.0 ) THEN
628  CALL sswap( j*nb, a( i1, 1 ), lda,
629  $ a( i2, 1 ), lda )
630  END IF
631  ENDIF
632  END DO
633 *
634 * Apply pivots to previous columns of L
635 *
636 c CALL SLASWP( J*NB, A( 1, 1 ), LDA,
637 c $ (J+1)*NB+1, (J+1)*NB+KB, IPIV, 1 )
638  END IF
639  END DO
640  END IF
641 *
642 * Factor the band matrix
643  CALL sgbtrf( n, n, nb, nb, tb, ldtb, ipiv2, info )
644 *
645 * End of SSYTRF_AA_2STAGE
646 *
subroutine strsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRSM
Definition: strsm.f:183
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:189
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: tstiee.f:83
subroutine sgetrf(M, N, A, LDA, IPIV, INFO)
SGETRF
Definition: sgetrf.f:110
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
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
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine ssygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
SSYGST
Definition: ssygst.f:129
subroutine sswap(N, SX, INCX, SY, INCY)
SSWAP
Definition: sswap.f:84
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 scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:84
subroutine sgbtrf(M, N, KL, KU, AB, LDAB, IPIV, INFO)
SGBTRF
Definition: sgbtrf.f:146
Here is the call graph for this function:
Here is the caller graph for this function: