LAPACK  3.9.0
LAPACK: Linear Algebra PACKage

◆ dsytrf_aa_2stage()

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

DSYTRF_AA_2STAGE

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

Purpose:
 DSYTRF_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*T*U  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 DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (LTB)
          On exit, details of the LU factorization of the band matrix.
[in]LTB
          LTB is INTEGER
          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
          IPIV2 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 IPIV2(k).
[out]WORK
          WORK is DOUBLE PRECISION workspace of size LWORK
[in]LWORK
          LWORK is INTEGER
          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 162 of file dsytrf_aa_2stage.f.

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