*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download CHESV_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chesv_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chesv_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chesv_aasen.f">
+*> Download CHESV_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chesv_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chesv_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chesv_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE CHESV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
* LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER INFO, LDA, LDB, LWORK, N, NRHS
* INTEGER IPIV( * )
* COMPLEX A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> A = U * T * U**H, if UPLO = 'U', or
*> A = L * T * L**H, if UPLO = 'L',
*> where U (or L) is a product of permutation and unit upper (lower)
-*> triangular matrices, and T is Hermitian and tridiagonal. The factored form
+*> triangular matrices, and T is Hermitian and tridiagonal. The factored form
*> of A is then used to solve the system of equations A * X = B.
*> \endverbatim
*
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup complexHEsolve
*
-* @generated from zhesv_aasen.f, fortran z -> c, Mon Oct 3 01:48:05 2016
-*
* =====================================================================
SUBROUTINE CHESV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
$ LWORK, INFO )
*
-* -- LAPACK driver routine (version 3.4.0) --
+* -- LAPACK driver routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download CHETRF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chetrf_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chetrf_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chetrf_aasen.f">
+*> Download CHETRF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chetrf_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chetrf_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chetrf_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE CHETRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, LDA, LWORK, INFO
* ..
* .. Array Arguments ..
* INTEGER IPIV( * )
-* COMPLEX A( LDA, * ), WORK( * )
+* COMPLEX A( LDA, * ), WORK( * )
* ..
*
*> \par Purpose:
*>
*> \verbatim
*>
-*> CHETRF_AASEN computes the factorization of a real hermitian matrix A
+*> CHETRF_AASEN computes the factorization of a complex hermitian matrix A
*> using the Aasen's algorithm. The form of the factorization is
*>
*> A = U*T*U**T or A = L*T*L**T
*> triangular part of A is not referenced.
*>
*> On exit, the tridiagonal matrix is stored in the diagonals
-*> and the subdiagonals of A just below (or above) the diagonals,
+*> and the subdiagonals of A just below (or above) the diagonals,
*> and L is stored below (or above) the subdiaonals, when UPLO
*> is 'L' (or 'U').
*> \endverbatim
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
-*> \ingroup complexSYcomputational
-*
-* @generated from zhetrf_aasen.f, fortran z -> c, Sun Oct 2 22:29:10 2016
+*> \ingroup complexHEcomputational
*
* =====================================================================
SUBROUTINE CHETRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
* ..
* .. Array Arguments ..
INTEGER IPIV( * )
- COMPLEX A( LDA, * ), WORK( * )
+ COMPLEX A( LDA, * ), WORK( * )
* ..
*
* =====================================================================
* .. Parameters ..
- COMPLEX ZERO, ONE
+ COMPLEX ZERO, ONE
PARAMETER ( ZERO = (0.0E+0, 0.0E+0), ONE = (1.0E+0, 0.0E+0) )
*
* .. Local Scalars ..
LOGICAL LQUERY, UPPER
INTEGER J, LWKOPT, IINFO
INTEGER NB, MJ, NJ, K1, K2, J1, J2, J3, JB
- COMPLEX ALPHA
+ COMPLEX ALPHA
* ..
* .. External Functions ..
LOGICAL LSAME
*
J = 0
10 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J + 1
*
* Panel factorization
*
- CALL CLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL CLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( MAX(1, J), J+1 ), LDA,
- $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
+ $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
$ IINFO )
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL CSWAP( J1-K1-2, A( 1, J2 ), 1,
+ CALL CSWAP( J1-K1-2, A( 1, J2 ), 1,
$ A( 1, IPIV(J2) ), 1 )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
+* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
* WORK stores the current block of the auxiriarly matrix H
*
IF( J.LT.N ) THEN
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block of J2-th block row with CGEMM
*
- CALL CGEMM( 'Conjugate transpose', 'Transpose',
+ CALL CGEMM( 'Conjugate transpose', 'Transpose',
$ NJ, N-J3+1, JB+1,
$ -ONE, A( J1-K2, J2 ), LDA,
$ WORK( (J3-J1+1)+K1*N ), N,
* Factorize A as L*D*L**T using the lower triangle of A
* .....................................................
*
-* copy first column A(1:N, 1) into H(1:N, 1)
+* copy first column A(1:N, 1) into H(1:N, 1)
* (stored in WORK(1:N))
*
CALL CCOPY( N, A( 1, 1 ), 1, WORK( 1 ), 1 )
*
J = 0
11 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J+1
*
* Panel factorization
*
- CALL CLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL CLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( J+1, MAX(1, J) ), LDA,
$ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ), IINFO)
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL CSWAP( J1-K1-2, A( J2, 1 ), LDA,
+ CALL CSWAP( J1-K1-2, A( J2, 1 ), LDA,
$ A( IPIV(J2), 1 ), LDA )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* A(J2+1, J1-1) stores L(J2+1, J1) and
+* A(J2+1, J1-1) stores L(J2+1, J1) and
* WORK(J2+1, 1) stores H(J2+1, 1)
*
IF( J.LT.N ) THEN
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block of J2-th block column with CGEMM
*
- CALL CGEMM( 'No transpose', 'Conjugate transpose',
+ CALL CGEMM( 'No transpose', 'Conjugate transpose',
$ N-J3+1, NJ, JB+1,
$ -ONE, WORK( (J3-J1+1)+K1*N ), N,
$ A( J2, J1-K2 ), LDA,
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download CHETRS_AASEN + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chetrs_aasen.f">
-*> [TGZ]</a>
+*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chetrs_aasen.f">
-*> [ZIP]</a>
+*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chetrs_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE CHETRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
* WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, NRHS, LDA, LDB, LWORK, INFO
* INTEGER IPIV( * )
* COMPLEX A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
-*> \ingroup complexSYcomputational
-*
-* @generated from zhetrs_aasen.f, fortran z -> c, Fri Sep 23 00:09:52 2016
+*> \ingroup complexHEcomputational
*
* =====================================================================
SUBROUTINE CHETRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
$ WORK, LWORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
PARAMETER ( ONE = 1.0E+0 )
* ..
* .. Local Scalars ..
- LOGICAL UPPER
- INTEGER K, KP
+ LOGICAL LQUERY, UPPER
+ INTEGER K, KP, LWKOPT
* ..
* .. External Functions ..
LOGICAL LSAME
*
INFO = 0
UPPER = LSAME( UPLO, 'U' )
+ LQUERY = ( LWORK.EQ.-1 )
IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -5
ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
INFO = -8
- ELSE IF( LWORK.LT.(3*N-2) ) THEN
+ ELSE IF( LWORK.LT.(3*N-2) .AND. .NOT.LQUERY ) THEN
INFO = -10
END IF
IF( INFO.NE.0 ) THEN
CALL XERBLA( 'CHETRS_AASEN', -INFO )
RETURN
+ ELSE IF( LQUERY ) THEN
+ LWKOPT = (3*N-2)
+ WORK( 1 ) = LWKOPT
+ RETURN
END IF
*
* Quick return if possible
*
* Compute (L \P**T * B) -> B [ (L \P**T * B) ]
*
- CALL CTRSM( 'L', 'L', 'N', 'U', N-1, NRHS, ONE, A( 2, 1), LDA,
+ CALL CTRSM( 'L', 'L', 'N', 'U', N-1, NRHS, ONE, A( 2, 1), LDA,
$ B(2, 1), LDB)
*
* Compute T \ B -> B [ T \ (L \P**T * B) ]
$ INFO)
*
* Compute (L**T \ B) -> B [ L**T \ (T \ (L \P**T * B) ) ]
-*
+*
CALL CTRSM( 'L', 'L', 'C', 'U', N-1, NRHS, ONE, A( 2, 1 ), LDA,
$ B( 2, 1 ), LDB)
*
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download CLAHEF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/clahef_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/clahef_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/clahef_aasen.f">
+*> Download CLAHEF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/clahef_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/clahef_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/clahef_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
-* SUBROUTINE CLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+* SUBROUTINE CLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
* H, LDH, WORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER J1, M, NB, LDA, LDH, INFO
* INTEGER IPIV( * )
* COMPLEX A( LDA, * ), H( LDH, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> last row, or column, of the previous panel. The first row, or column,
*> of A is set to be the first row, or column, of an identity matrix,
*> which is used to factorize the first panel.
-*>
+*>
*> The resulting J-th row of U, or J-th column of L, is stored in the
-*> (J-1)-th row, or column, of A (without the unit diatonals), while
+*> (J-1)-th row, or column, of A (without the unit diatonals), while
*> the diagonal and subdiagonal of A are overwritten by those of T.
*>
*> \endverbatim
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup complexSYcomputational
*
-* @generated from zlahef_aasen.f, fortran z -> c, Sun Oct 2 22:41:33 2016
-*
* =====================================================================
- SUBROUTINE CLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+ SUBROUTINE CLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
$ H, LDH, WORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* .. Local Scalars ..
INTEGER J, K, K1, I1, I2
- COMPLEX PIV, ALPHA
+ COMPLEX PIV, ALPHA
* ..
* .. External Functions ..
LOGICAL LSAME
*
A( K, J ) = REAL( WORK( 1 ) )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L(J, (J+1):N)
* where A(J, J) stores T(J, J) and A(J-1, (J+1):N) stores U(J, (J+1):N)
*
IF( (J1+J-1).GT.1 ) THEN
- ALPHA = -A( K, J )
- CALL CAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
+ ALPHA = -A( K, J )
+ CALL CAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL CSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
+ CALL CSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
$ A( J1+I1, I2 ), 1 )
CALL CLACGV( I2-I1, A( J1+I1-1, I1+1 ), LDA )
CALL CLACGV( I2-I1-1, A( J1+I1, I2 ), 1 )
*
* Swap A(I1, I2+1:N) with A(I2, I2+1:N)
*
- CALL CSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
+ CALL CSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
$ A( J1+I2-1, I2+1 ), LDA )
*
* Swap A(I1, I1) with A(I2,I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL CSWAP( I1-K1+1, A( 1, I1 ), 1,
+ CALL CSWAP( I1-K1+1, A( 1, I1 ), 1,
$ A( 1, I2 ), 1 )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J, J+1) = T(J, J+1)
*
A( K, J+1 ) = WORK( 2 )
- IF( (A( K, J ).EQ.ZERO ) .AND.
+ IF( (A( K, J ).EQ.ZERO ) .AND.
$ ( (J.EQ.M) .OR. (A( K, J+1 ).EQ.ZERO))) THEN
IF(INFO .EQ. 0) THEN
INFO = J
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J:N, J),
+* Copy A(J+1:N, J+1) into H(J:N, J),
*
- CALL CCOPY( M-J, A( K+1, J+1 ), LDA,
+ CALL CCOPY( M-J, A( K+1, J+1 ), LDA,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL CCOPY( M-J-1, WORK( 3 ), 1, A( K, J+2 ), LDA )
CALL CSCAL( M-J-1, ALPHA, A( K, J+2 ), LDA )
ELSE
- CALL CLASET( 'Full', 1, M-J-1, ZERO, ZERO,
+ CALL CLASET( 'Full', 1, M-J-1, ZERO, ZERO,
$ A( K, J+2 ), LDA)
END IF
ELSE
*
A( J, K ) = REAL( WORK( 1 ) )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L((J+1):N, J)
* where A(J, J) = T(J, J) and A((J+1):N, J-1) = L((J+1):N, J)
*
IF( (J1+J-1).GT.1 ) THEN
ALPHA = -A( J, K )
- CALL CAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
+ CALL CAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL CSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
+ CALL CSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
$ A( I2, J1+I1 ), LDA )
CALL CLACGV( I2-I1, A( I1+1, J1+I1-1 ), 1 )
CALL CLACGV( I2-I1-1, A( I2, J1+I1 ), LDA )
*
* Swap A(I2+1:N, I1) with A(I2+1:N, I2)
*
- CALL CSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
+ CALL CSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
$ A( I2+1, J1+I2-1 ), 1 )
*
* Swap A(I1, I1) with A(I2, I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL CSWAP( I1-K1+1, A( I1, 1 ), LDA,
+ CALL CSWAP( I1-K1+1, A( I1, 1 ), LDA,
$ A( I2, 1 ), LDA )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J+1, J) = T(J+1, J)
*
A( J+1, K ) = WORK( 2 )
- IF( (A( J, K ).EQ.ZERO) .AND.
+ IF( (A( J, K ).EQ.ZERO) .AND.
$ ( (J.EQ.M) .OR. (A( J+1, K ).EQ.ZERO)) ) THEN
- IF (INFO .EQ. 0)
+ IF (INFO .EQ. 0)
$ INFO = J
END IF
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J+1:N, J),
+* Copy A(J+1:N, J+1) into H(J+1:N, J),
*
- CALL CCOPY( M-J, A( J+1, K+1 ), 1,
+ CALL CCOPY( M-J, A( J+1, K+1 ), 1,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL CCOPY( M-J-1, WORK( 3 ), 1, A( J+2, K ), 1 )
CALL CSCAL( M-J-1, ALPHA, A( J+2, K ), 1 )
ELSE
- CALL CLASET( 'Full', M-J-1, 1, ZERO, ZERO,
+ CALL CLASET( 'Full', M-J-1, 1, ZERO, ZERO,
$ A( J+2, K ), LDA )
END IF
ELSE
- IF( (A( J, K ).EQ.ZERO) .AND. (J.EQ.M)
+ IF( (A( J, K ).EQ.ZERO) .AND. (J.EQ.M)
$ .AND. (INFO.EQ.0) ) INFO = J
END IF
J = J + 1
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download DLASYF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasyf_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasyf_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasyf_aasen.f">
+*> Download DLASYF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasyf_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasyf_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasyf_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
-* SUBROUTINE DLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+* SUBROUTINE DLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
* H, LDH, WORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER J1, M, NB, LDA, LDH, INFO
* INTEGER IPIV( * )
* DOUBLE PRECISION A( LDA, * ), H( LDH, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> last row, or column, of the previous panel. The first row, or column,
*> of A is set to be the first row, or column, of an identity matrix,
*> which is used to factorize the first panel.
-*>
+*>
*> The resulting J-th row of U, or J-th column of L, is stored in the
-*> (J-1)-th row, or column, of A (without the unit diatonals), while
+*> (J-1)-th row, or column, of A (without the unit diatonals), while
*> the diagonal and subdiagonal of A are overwritten by those of T.
*>
*> \endverbatim
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup doubleSYcomputational
*
-* @precisions fortran d -> s
-*
* =====================================================================
- SUBROUTINE DLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+ SUBROUTINE DLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
$ H, LDH, WORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* .. Local Scalars ..
INTEGER J, K, K1, I1, I2
- DOUBLE PRECISION PIV, ALPHA
+ DOUBLE PRECISION PIV, ALPHA
* ..
* .. External Functions ..
LOGICAL LSAME
*
A( K, J ) = WORK( 1 )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L(J, (J+1):N)
* where A(J, J) stores T(J, J) and A(J-1, (J+1):N) stores U(J, (J+1):N)
*
IF( (J1+J-1).GT.1 ) THEN
- ALPHA = -A( K, J )
- CALL DAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
+ ALPHA = -A( K, J )
+ CALL DAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL DSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
+ CALL DSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
$ A( J1+I1, I2 ), 1 )
*
* Swap A(I1, I2+1:N) with A(I2, I2+1:N)
*
- CALL DSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
+ CALL DSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
$ A( J1+I2-1, I2+1 ), LDA )
*
* Swap A(I1, I1) with A(I2,I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL DSWAP( I1-K1+1, A( 1, I1 ), 1,
+ CALL DSWAP( I1-K1+1, A( 1, I1 ), 1,
$ A( 1, I2 ), 1 )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J, J+1) = T(J, J+1)
*
A( K, J+1 ) = WORK( 2 )
- IF( (A( K, J ).EQ.ZERO ) .AND.
+ IF( (A( K, J ).EQ.ZERO ) .AND.
$ ( (J.EQ.M) .OR. (A( K, J+1 ).EQ.ZERO))) THEN
IF(INFO .EQ. 0) THEN
INFO = J
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J:N, J),
+* Copy A(J+1:N, J+1) into H(J:N, J),
*
- CALL DCOPY( M-J, A( K+1, J+1 ), LDA,
+ CALL DCOPY( M-J, A( K+1, J+1 ), LDA,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL DCOPY( M-J-1, WORK( 3 ), 1, A( K, J+2 ), LDA )
CALL DSCAL( M-J-1, ALPHA, A( K, J+2 ), LDA )
ELSE
- CALL DLASET( 'Full', 1, M-J-1, ZERO, ZERO,
+ CALL DLASET( 'Full', 1, M-J-1, ZERO, ZERO,
$ A( K, J+2 ), LDA)
END IF
ELSE
*
A( J, K ) = WORK( 1 )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L((J+1):N, J)
* where A(J, J) = T(J, J) and A((J+1):N, J-1) = L((J+1):N, J)
*
IF( (J1+J-1).GT.1 ) THEN
- ALPHA = -A( J, K )
- CALL DAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
+ ALPHA = -A( J, K )
+ CALL DAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL DSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
+ CALL DSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
$ A( I2, J1+I1 ), LDA )
*
* Swap A(I2+1:N, I1) with A(I2+1:N, I2)
*
- CALL DSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
+ CALL DSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
$ A( I2+1, J1+I2-1 ), 1 )
*
* Swap A(I1, I1) with A(I2, I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL DSWAP( I1-K1+1, A( I1, 1 ), LDA,
+ CALL DSWAP( I1-K1+1, A( I1, 1 ), LDA,
$ A( I2, 1 ), LDA )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J+1, J) = T(J+1, J)
*
A( J+1, K ) = WORK( 2 )
- IF( (A( J, K ).EQ.ZERO) .AND.
+ IF( (A( J, K ).EQ.ZERO) .AND.
$ ( (J.EQ.M) .OR. (A( J+1, K ).EQ.ZERO)) ) THEN
- IF (INFO .EQ. 0)
+ IF (INFO .EQ. 0)
$ INFO = J
END IF
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J+1:N, J),
+* Copy A(J+1:N, J+1) into H(J+1:N, J),
*
- CALL DCOPY( M-J, A( J+1, K+1 ), 1,
+ CALL DCOPY( M-J, A( J+1, K+1 ), 1,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL DCOPY( M-J-1, WORK( 3 ), 1, A( J+2, K ), 1 )
CALL DSCAL( M-J-1, ALPHA, A( J+2, K ), 1 )
ELSE
- CALL DLASET( 'Full', M-J-1, 1, ZERO, ZERO,
+ CALL DLASET( 'Full', M-J-1, 1, ZERO, ZERO,
$ A( J+2, K ), LDA )
END IF
ELSE
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download DSYSV_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsysv_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsysv_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsysv_aasen.f">
+*> Download DSYSV_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsysv_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsysv_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsysv_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE DSYSV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
* LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, NRHS, LDA, LDB, LWORK, INFO
* INTEGER IPIV( * )
* DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> A = U * T * U**T, if UPLO = 'U', or
*> A = L * T * L**T, if UPLO = 'L',
*> where U (or L) is a product of permutation and unit upper (lower)
-*> triangular matrices, and T is symmetric tridiagonal. The factored
+*> triangular matrices, and T is symmetric tridiagonal. The factored
*> form of A is then used to solve the system of equations A * X = B.
*> \endverbatim
*
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
-*> The length of WORK. LWORK >= MAX(2*N, 3*N-2), and for
-*> the best performance, LWORK >= max(1,N*NB), where NB is
+*> The length of WORK. LWORK >= MAX(2*N, 3*N-2), and for
+*> the best performance, LWORK >= max(1,N*NB), where NB is
*> the optimal blocksize for DSYTRF_AASEN.
*>
*> If LWORK = -1, then a workspace query is assumed; the routine
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup doubleSYsolve
*
-* @precisions fortran d -> s
-*
* =====================================================================
SUBROUTINE DSYSV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
$ LWORK, INFO )
*
-* -- LAPACK driver routine (version 3.4.0) --
+* -- LAPACK driver routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download DSYTRF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsytrf_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsytrf_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsytrf_aasen.f">
+*> Download DSYTRF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsytrf_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsytrf_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsytrf_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE DSYTRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, LDA, LWORK, INFO
*> triangular part of A is not referenced.
*>
*> On exit, the tridiagonal matrix is stored in the diagonals
-*> and the subdiagonals of A just below (or above) the diagonals,
+*> and the subdiagonals of A just below (or above) the diagonals,
*> and L is stored below (or above) the subdiaonals, when UPLO
*> is 'L' (or 'U').
*> \endverbatim
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup doubleSYcomputational
*
-* @precisions fortran d -> s
-*
* =====================================================================
SUBROUTINE DSYTRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
J = 0
10 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J + 1
*
* Panel factorization
*
- CALL DLASYF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL DLASYF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( MAX(1, J), J+1 ), LDA,
- $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
+ $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
$ IINFO )
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL DSWAP( J1-K1-2, A( 1, J2 ), 1,
+ CALL DSWAP( J1-K1-2, A( 1, J2 ), 1,
$ A( 1, IPIV(J2) ), 1 )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
+* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
* WORK stores the current block of the auxiriarly matrix H
*
IF( J.LT.N ) THEN
*
ALPHA = A( J, J+1 )
A( J, J+1 ) = ONE
- CALL DCOPY( N-J, A( J-1, J+1 ), LDA,
+ CALL DCOPY( N-J, A( J-1, J+1 ), LDA,
$ WORK( (J+1-J1+1)+JB*N ), 1 )
CALL DSCAL( N-J, ALPHA, WORK( (J+1-J1+1)+JB*N ), 1 )
*
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
+* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
* while K1=0 and K2=1 for the rest
*
IF( J1.GT.1 ) THEN
* Not first panel
*
K2 = 1
- ELSE
+ ELSE
*
* First panel
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block of J2-th block row with DGEMM
*
- CALL DGEMM( 'Transpose', 'Transpose',
+ CALL DGEMM( 'Transpose', 'Transpose',
$ NJ, N-J3+1, JB+1,
$ -ONE, A( J1-K2, J2 ), LDA,
$ WORK( J3-J1+1+K1*N ), N,
* Factorize A as L*D*L**T using the lower triangle of A
* .....................................................
*
-* copy first column A(1:N, 1) into H(1:N, 1)
+* copy first column A(1:N, 1) into H(1:N, 1)
* (stored in WORK(1:N))
*
CALL DCOPY( N, A( 1, 1 ), 1, WORK( 1 ), 1 )
*
J = 0
11 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J+1
*
* Panel factorization
*
- CALL DLASYF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL DLASYF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( J+1, MAX(1, J) ), LDA,
$ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ), IINFO)
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL DSWAP( J1-K1-2, A( J2, 1 ), LDA,
+ CALL DSWAP( J1-K1-2, A( J2, 1 ), LDA,
$ A( IPIV(J2), 1 ), LDA )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* A(J2+1, J1-1) stores L(J2+1, J1) and
+* A(J2+1, J1-1) stores L(J2+1, J1) and
* WORK(J2+1, 1) stores H(J2+1, 1)
*
IF( J.LT.N ) THEN
*
ALPHA = A( J+1, J )
A( J+1, J ) = ONE
- CALL DCOPY( N-J, A( J+1, J-1 ), 1,
+ CALL DCOPY( N-J, A( J+1, J-1 ), 1,
$ WORK( (J+1-J1+1)+JB*N ), 1 )
CALL DSCAL( N-J, ALPHA, WORK( (J+1-J1+1)+JB*N ), 1 )
*
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
+* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
* while K1=0 and K2=1 for the rest
*
IF( J1.GT.1 ) THEN
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block in J2-th block column with DGEMM
*
- CALL DGEMM( 'No transpose', 'Transpose',
+ CALL DGEMM( 'No transpose', 'Transpose',
$ N-J3+1, NJ, JB+1,
$ -ONE, WORK( J3-J1+1+K1*N ), N,
$ A( J2, J1-K2 ), LDA,
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download DSYTRS_AASEN + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsytrs_aasen.f">
-*> [TGZ]</a>
+*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsytrs_aasen.f">
-*> [ZIP]</a>
+*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsytrs_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE DSYTRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
* WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, NRHS, LDA, LDB, LWORK, INFO
* INTEGER IPIV( * )
* DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup doubleSYcomputational
*
-* @precisions fortran d -> s
-*
* =====================================================================
SUBROUTINE DSYTRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
$ WORK, LWORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
PARAMETER ( ONE = 1.0D+0 )
* ..
* .. Local Scalars ..
- LOGICAL UPPER
- INTEGER K, KP
+ LOGICAL LQUERY, UPPER
+ INTEGER K, KP, LWKOPT
* ..
* .. External Functions ..
LOGICAL LSAME
*
INFO = 0
UPPER = LSAME( UPLO, 'U' )
+ LQUERY = ( LWORK.EQ.-1 )
IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -5
ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
INFO = -8
- ELSE IF( LWORK.LT.(3*N-2) ) THEN
+ ELSE IF( LWORK.LT.(3*N-2) .AND. .NOT.LQUERY ) THEN
INFO = -10
END IF
IF( INFO.NE.0 ) THEN
CALL XERBLA( 'DSYTRS_AASEN', -INFO )
RETURN
+ ELSE IF( LQUERY ) THEN
+ LWKOPT = (3*N-2)
+ WORK( 1 ) = LWKOPT
+ RETURN
END IF
*
* Quick return if possible
$ INFO)
*
* Compute (L**T \ B) -> B [ L**T \ (T \ (L \P**T * B) ) ]
-*
+*
CALL DTRSM( 'L', 'L', 'T', 'U', N-1, NRHS, ONE, A( 2, 1 ), LDA,
$ B( 2, 1 ), LDB)
*
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download SLASYF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slasyf_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slasyf_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slasyf_aasen.f">
+*> Download SLASYF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slasyf_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slasyf_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slasyf_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
-* SUBROUTINE SLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+* SUBROUTINE SLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
* H, LDH, WORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER J1, M, NB, LDA, LDH, INFO
* INTEGER IPIV( * )
* REAL A( LDA, * ), H( LDH, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> last row, or column, of the previous panel. The first row, or column,
*> of A is set to be the first row, or column, of an identity matrix,
*> which is used to factorize the first panel.
-*>
+*>
*> The resulting J-th row of U, or J-th column of L, is stored in the
-*> (J-1)-th row, or column, of A (without the unit diatonals), while
+*> (J-1)-th row, or column, of A (without the unit diatonals), while
*> the diagonal and subdiagonal of A are overwritten by those of T.
*>
*> \endverbatim
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup realSYcomputational
*
-* @generated from dlasyf_aasen.f, fortran d -> s, Sun Oct 2 22:57:56 2016
-*
* =====================================================================
- SUBROUTINE SLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+ SUBROUTINE SLASYF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
$ H, LDH, WORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* .. Local Scalars ..
INTEGER J, K, K1, I1, I2
- REAL PIV, ALPHA
+ REAL PIV, ALPHA
* ..
* .. External Functions ..
LOGICAL LSAME
*
A( K, J ) = WORK( 1 )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L(J, (J+1):N)
* where A(J, J) stores T(J, J) and A(J-1, (J+1):N) stores U(J, (J+1):N)
*
IF( (J1+J-1).GT.1 ) THEN
- ALPHA = -A( K, J )
- CALL SAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
+ ALPHA = -A( K, J )
+ CALL SAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL SSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
+ CALL SSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
$ A( J1+I1, I2 ), 1 )
*
* Swap A(I1, I2+1:N) with A(I2, I2+1:N)
*
- CALL SSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
+ CALL SSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
$ A( J1+I2-1, I2+1 ), LDA )
*
* Swap A(I1, I1) with A(I2,I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL SSWAP( I1-K1+1, A( 1, I1 ), 1,
+ CALL SSWAP( I1-K1+1, A( 1, I1 ), 1,
$ A( 1, I2 ), 1 )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J, J+1) = T(J, J+1)
*
A( K, J+1 ) = WORK( 2 )
- IF( (A( K, J ).EQ.ZERO ) .AND.
+ IF( (A( K, J ).EQ.ZERO ) .AND.
$ ( (J.EQ.M) .OR. (A( K, J+1 ).EQ.ZERO))) THEN
IF(INFO .EQ. 0) THEN
INFO = J
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J:N, J),
+* Copy A(J+1:N, J+1) into H(J:N, J),
*
- CALL SCOPY( M-J, A( K+1, J+1 ), LDA,
+ CALL SCOPY( M-J, A( K+1, J+1 ), LDA,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL SCOPY( M-J-1, WORK( 3 ), 1, A( K, J+2 ), LDA )
CALL SSCAL( M-J-1, ALPHA, A( K, J+2 ), LDA )
ELSE
- CALL SLASET( 'Full', 1, M-J-1, ZERO, ZERO,
+ CALL SLASET( 'Full', 1, M-J-1, ZERO, ZERO,
$ A( K, J+2 ), LDA)
END IF
ELSE
*
A( J, K ) = WORK( 1 )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L((J+1):N, J)
* where A(J, J) = T(J, J) and A((J+1):N, J-1) = L((J+1):N, J)
*
IF( (J1+J-1).GT.1 ) THEN
- ALPHA = -A( J, K )
- CALL SAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
+ ALPHA = -A( J, K )
+ CALL SAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL SSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
+ CALL SSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
$ A( I2, J1+I1 ), LDA )
*
* Swap A(I2+1:N, I1) with A(I2+1:N, I2)
*
- CALL SSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
+ CALL SSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
$ A( I2+1, J1+I2-1 ), 1 )
*
* Swap A(I1, I1) with A(I2, I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL SSWAP( I1-K1+1, A( I1, 1 ), LDA,
+ CALL SSWAP( I1-K1+1, A( I1, 1 ), LDA,
$ A( I2, 1 ), LDA )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J+1, J) = T(J+1, J)
*
A( J+1, K ) = WORK( 2 )
- IF( (A( J, K ).EQ.ZERO) .AND.
+ IF( (A( J, K ).EQ.ZERO) .AND.
$ ( (J.EQ.M) .OR. (A( J+1, K ).EQ.ZERO)) ) THEN
- IF (INFO .EQ. 0)
+ IF (INFO .EQ. 0)
$ INFO = J
END IF
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J+1:N, J),
+* Copy A(J+1:N, J+1) into H(J+1:N, J),
*
- CALL SCOPY( M-J, A( J+1, K+1 ), 1,
+ CALL SCOPY( M-J, A( J+1, K+1 ), 1,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL SCOPY( M-J-1, WORK( 3 ), 1, A( J+2, K ), 1 )
CALL SSCAL( M-J-1, ALPHA, A( J+2, K ), 1 )
ELSE
- CALL SLASET( 'Full', M-J-1, 1, ZERO, ZERO,
+ CALL SLASET( 'Full', M-J-1, 1, ZERO, ZERO,
$ A( J+2, K ), LDA )
END IF
ELSE
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download SSYSV_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssysv_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssysv_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssysv_aasen.f">
+*> Download SSYSV_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssysv_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssysv_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssysv_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE SSYSV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
* LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, NRHS, LDA, LDB, LWORK, INFO
* INTEGER IPIV( * )
* REAL A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> A = U * T * U**T, if UPLO = 'U', or
*> A = L * T * L**T, if UPLO = 'L',
*> where U (or L) is a product of permutation and unit upper (lower)
-*> triangular matrices, and T is symmetric tridiagonal. The factored
+*> triangular matrices, and T is symmetric tridiagonal. The factored
*> form of A is then used to solve the system of equations A * X = B.
*> \endverbatim
*
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
-*> The length of WORK. LWORK >= MAX(2*N, 3*N-2), and for
-*> the best performance, LWORK >= max(1,N*NB), where NB is
+*> The length of WORK. LWORK >= MAX(2*N, 3*N-2), and for
+*> the best performance, LWORK >= max(1,N*NB), where NB is
*> the optimal blocksize for SSYTRF_AASEN.
*>
*> If LWORK = -1, then a workspace query is assumed; the routine
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup realSYsolve
*
-* @generated from dsysv_aasen.f, fortran d -> s, Mon Oct 3 01:04:05 2016
-*
* =====================================================================
SUBROUTINE SSYSV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
$ LWORK, INFO )
*
-* -- LAPACK driver routine (version 3.4.0) --
+* -- LAPACK driver routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download SSYTRF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssytrf_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssytrf_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssytrf_aasen.f">
+*> Download SSYTRF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssytrf_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssytrf_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssytrf_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE SSYTRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, LDA, LWORK, INFO
*> triangular part of A is not referenced.
*>
*> On exit, the tridiagonal matrix is stored in the diagonals
-*> and the subdiagonals of A just below (or above) the diagonals,
+*> and the subdiagonals of A just below (or above) the diagonals,
*> and L is stored below (or above) the subdiaonals, when UPLO
*> is 'L' (or 'U').
*> \endverbatim
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup realSYcomputational
*
-* @generated from dsytrf_aasen.f, fortran d -> s, Sun Oct 2 22:27:17 2016
-*
* =====================================================================
SUBROUTINE SSYTRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
J = 0
10 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J + 1
*
* Panel factorization
*
- CALL SLASYF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL SLASYF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( MAX(1, J), J+1 ), LDA,
- $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
+ $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
$ IINFO )
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL SSWAP( J1-K1-2, A( 1, J2 ), 1,
+ CALL SSWAP( J1-K1-2, A( 1, J2 ), 1,
$ A( 1, IPIV(J2) ), 1 )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
+* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
* WORK stores the current block of the auxiriarly matrix H
*
IF( J.LT.N ) THEN
*
ALPHA = A( J, J+1 )
A( J, J+1 ) = ONE
- CALL SCOPY( N-J, A( J-1, J+1 ), LDA,
+ CALL SCOPY( N-J, A( J-1, J+1 ), LDA,
$ WORK( (J+1-J1+1)+JB*N ), 1 )
CALL SSCAL( N-J, ALPHA, WORK( (J+1-J1+1)+JB*N ), 1 )
*
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
+* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
* while K1=0 and K2=1 for the rest
*
IF( J1.GT.1 ) THEN
* Not first panel
*
K2 = 1
- ELSE
+ ELSE
*
* First panel
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block of J2-th block row with SGEMM
*
- CALL SGEMM( 'Transpose', 'Transpose',
+ CALL SGEMM( 'Transpose', 'Transpose',
$ NJ, N-J3+1, JB+1,
$ -ONE, A( J1-K2, J2 ), LDA,
$ WORK( J3-J1+1+K1*N ), N,
* Factorize A as L*D*L**T using the lower triangle of A
* .....................................................
*
-* copy first column A(1:N, 1) into H(1:N, 1)
+* copy first column A(1:N, 1) into H(1:N, 1)
* (stored in WORK(1:N))
*
CALL SCOPY( N, A( 1, 1 ), 1, WORK( 1 ), 1 )
*
J = 0
11 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J+1
*
* Panel factorization
*
- CALL SLASYF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL SLASYF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( J+1, MAX(1, J) ), LDA,
$ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ), IINFO)
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL SSWAP( J1-K1-2, A( J2, 1 ), LDA,
+ CALL SSWAP( J1-K1-2, A( J2, 1 ), LDA,
$ A( IPIV(J2), 1 ), LDA )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* A(J2+1, J1-1) stores L(J2+1, J1) and
+* A(J2+1, J1-1) stores L(J2+1, J1) and
* WORK(J2+1, 1) stores H(J2+1, 1)
*
IF( J.LT.N ) THEN
*
ALPHA = A( J+1, J )
A( J+1, J ) = ONE
- CALL SCOPY( N-J, A( J+1, J-1 ), 1,
+ CALL SCOPY( N-J, A( J+1, J-1 ), 1,
$ WORK( (J+1-J1+1)+JB*N ), 1 )
CALL SSCAL( N-J, ALPHA, WORK( (J+1-J1+1)+JB*N ), 1 )
*
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
+* explicitly stored, e.g., K1=1 and K2= 0 for the first panel,
* while K1=0 and K2=1 for the rest
*
IF( J1.GT.1 ) THEN
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block in J2-th block column with SGEMM
*
- CALL SGEMM( 'No transpose', 'Transpose',
+ CALL SGEMM( 'No transpose', 'Transpose',
$ N-J3+1, NJ, JB+1,
$ -ONE, WORK( J3-J1+1+K1*N ), N,
$ A( J2, J1-K2 ), LDA,
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download SSYTRS_AASEN + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssytrs_aasen.f">
-*> [TGZ]</a>
+*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssytrs_aasen.f">
-*> [ZIP]</a>
+*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssytrs_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE SSYTRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
* WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, NRHS, LDA, LDB, LWORK, INFO
* INTEGER IPIV( * )
* REAL A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup realSYcomputational
*
-* @generated from dsytrs_aasen.f, fortran d -> s, Wed Sep 21 16:39:24 2016
-*
* =====================================================================
SUBROUTINE SSYTRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
$ WORK, LWORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* =====================================================================
*
- REAL ONE
+ REAL ONE
PARAMETER ( ONE = 1.0E+0 )
* ..
* .. Local Scalars ..
- LOGICAL UPPER
- INTEGER K, KP
+ LOGICAL LQUERY, UPPER
+ INTEGER K, KP, LWKOPT
* ..
* .. External Functions ..
LOGICAL LSAME
*
INFO = 0
UPPER = LSAME( UPLO, 'U' )
+ LQUERY = ( LWORK.EQ.-1 )
IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -5
ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
INFO = -8
- ELSE IF( LWORK.LT.(3*N-2) ) THEN
+ ELSE IF( LWORK.LT.(3*N-2) .AND. .NOT.LQUERY ) THEN
INFO = -10
END IF
IF( INFO.NE.0 ) THEN
CALL XERBLA( 'SSYTRS_AASEN', -INFO )
RETURN
+ ELSE IF( LQUERY ) THEN
+ LWKOPT = (3*N-2)
+ WORK( 1 ) = LWKOPT
+ RETURN
END IF
*
* Quick return if possible
END IF
CALL SGTSV(N, NRHS, WORK(1), WORK(N), WORK(2*N), B, LDB,
$ INFO)
-*
+*
*
* Compute (U**T \ B) -> B [ U**T \ (T \ (U \P**T * B) ) ]
*
*
* Compute (L \P**T * B) -> B [ (L \P**T * B) ]
*
- CALL STRSM( 'L', 'L', 'N', 'U', N-1, NRHS, ONE, A( 2, 1), LDA,
+ CALL STRSM( 'L', 'L', 'N', 'U', N-1, NRHS, ONE, A( 2, 1), LDA,
$ B(2, 1), LDB)
*
* Compute T \ B -> B [ T \ (L \P**T * B) ]
$ INFO)
*
* Compute (L**T \ B) -> B [ L**T \ (T \ (L \P**T * B) ) ]
-*
+*
CALL STRSM( 'L', 'L', 'T', 'U', N-1, NRHS, ONE, A( 2, 1 ), LDA,
$ B( 2, 1 ), LDB)
*
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download ZHESV_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhesv_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhesv_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhesv_aasen.f">
+*> Download ZHESV_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhesv_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhesv_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhesv_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE ZHESV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
* LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER INFO, LDA, LDB, LWORK, N, NRHS
* INTEGER IPIV( * )
* COMPLEX*16 A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> A = U * T * U**H, if UPLO = 'U', or
*> A = L * T * L**H, if UPLO = 'L',
*> where U (or L) is a product of permutation and unit upper (lower)
-*> triangular matrices, and T is Hermitian and tridiagonal. The factored form
+*> triangular matrices, and T is Hermitian and tridiagonal. The factored form
*> of A is then used to solve the system of equations A * X = B.
*> \endverbatim
*
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
*> \ingroup complex16HEsolve
*
-* @precisions fortran z -> c
-*
* =====================================================================
SUBROUTINE ZHESV_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
$ LWORK, INFO )
*
-* -- LAPACK driver routine (version 3.4.0) --
+* -- LAPACK driver routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download ZHETRF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhetrf_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhetrf_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhetrf_aasen.f">
+*> Download ZHETRF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhetrf_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhetrf_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhetrf_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE ZHETRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, LDA, LWORK, INFO
*>
*> \verbatim
*>
-*> ZHETRF_AASEN computes the factorization of a real hermitian matrix A
+*> ZHETRF_AASEN computes the factorization of a complex hermitian matrix A
*> using the Aasen's algorithm. The form of the factorization is
*>
*> A = U*T*U**T or A = L*T*L**T
*> triangular part of A is not referenced.
*>
*> On exit, the tridiagonal matrix is stored in the diagonals
-*> and the subdiagonals of A just below (or above) the diagonals,
+*> and the subdiagonals of A just below (or above) the diagonals,
*> and L is stored below (or above) the subdiaonals, when UPLO
*> is 'L' (or 'U').
*> \endverbatim
*> \param[out] IPIV
*> \verbatim
*> 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
+*> 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).
*> \endverbatim
*>
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
-*> \ingroup complex16SYcomputational
-*
-* @precisions fortran z -> c
+*> \ingroup complex16HEcomputational
*
* =====================================================================
SUBROUTINE ZHETRF_AASEN( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
J = 0
10 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J + 1
*
* Panel factorization
*
- CALL ZLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL ZLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( MAX(1, J), J+1 ), LDA,
- $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
+ $ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ),
$ IINFO )
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL ZSWAP( J1-K1-2, A( 1, J2 ), 1,
+ CALL ZSWAP( J1-K1-2, A( 1, J2 ), 1,
$ A( 1, IPIV(J2) ), 1 )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
+* the row A(J1-1, J2-1:N) stores U(J1, J2+1:N) and
* WORK stores the current block of the auxiriarly matrix H
*
IF( J.LT.N ) THEN
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block of J2-th block row with ZGEMM
*
- CALL ZGEMM( 'Conjugate transpose', 'Transpose',
+ CALL ZGEMM( 'Conjugate transpose', 'Transpose',
$ NJ, N-J3+1, JB+1,
$ -ONE, A( J1-K2, J2 ), LDA,
$ WORK( (J3-J1+1)+K1*N ), N,
* Factorize A as L*D*L**T using the lower triangle of A
* .....................................................
*
-* copy first column A(1:N, 1) into H(1:N, 1)
+* copy first column A(1:N, 1) into H(1:N, 1)
* (stored in WORK(1:N))
*
CALL ZCOPY( N, A( 1, 1 ), 1, WORK( 1 ), 1 )
*
J = 0
11 CONTINUE
- IF( J.GE.N )
+ IF( J.GE.N )
$ GO TO 20
*
* each step of the main loop
* J is the last column of the previous panel
* J1 is the first column of the current panel
* K1 identifies if the previous column of the panel has been
-* explicitly stored, e.g., K1=1 for the first panel, and
+* explicitly stored, e.g., K1=1 for the first panel, and
* K1=0 for the rest
*
J1 = J+1
*
* Panel factorization
*
- CALL ZLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
+ CALL ZLAHEF_AASEN( UPLO, 2-K1, N-J, JB,
$ A( J+1, MAX(1, J) ), LDA,
$ IPIV( J+1 ), WORK, N, WORK( N*NB+1 ), IINFO)
IF( (IINFO.GT.0) .AND. (INFO.EQ.0) ) THEN
INFO = IINFO+J
- ENDIF
+ ENDIF
*
* Ajust IPIV and apply it back (J-th step picks (J+1)-th pivot)
*
DO J2 = J+2, MIN(N, J+JB+1)
IPIV( J2 ) = IPIV( J2 ) + J
IF( (J2.NE.IPIV(J2)) .AND. ((J1-K1).GT.2) ) THEN
- CALL ZSWAP( J1-K1-2, A( J2, 1 ), LDA,
+ CALL ZSWAP( J1-K1-2, A( J2, 1 ), LDA,
$ A( IPIV(J2), 1 ), LDA )
END IF
END DO
J = J + JB
*
* Trailing submatrix update, where
-* A(J2+1, J1-1) stores L(J2+1, J1) and
+* A(J2+1, J1-1) stores L(J2+1, J1) and
* WORK(J2+1, 1) stores H(J2+1, 1)
*
IF( J.LT.N ) THEN
*
K2 = 0
*
-* First update skips the first column
+* First update skips the first column
*
JB = JB - 1
END IF
*
* Update off-diagonal block of J2-th block column with ZGEMM
*
- CALL ZGEMM( 'No transpose', 'Conjugate transpose',
+ CALL ZGEMM( 'No transpose', 'Conjugate transpose',
$ N-J3+1, NJ, JB+1,
$ -ONE, WORK( (J3-J1+1)+K1*N ), N,
$ A( J2, J1-K2 ), LDA,
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download ZHETRS_AASEN + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhetrs_aasen.f">
-*> [TGZ]</a>
+*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhetrs_aasen.f">
-*> [ZIP]</a>
+*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhetrs_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE ZHETRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
* WORK, LWORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER N, NRHS, LDA, LDB, LWORK, INFO
* INTEGER IPIV( * )
* COMPLEX*16 A( LDA, * ), B( LDB, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
-*> \ingroup complex16SYcomputational
-*
-* @precisions fortran z -> c
+*> \ingroup complex16HEcomputational
*
* =====================================================================
SUBROUTINE ZHETRS_AASEN( UPLO, N, NRHS, A, LDA, IPIV, B, LDB,
$ WORK, LWORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
PARAMETER ( ONE = 1.0D+0 )
* ..
* .. Local Scalars ..
- LOGICAL UPPER
- INTEGER K, KP
+ LOGICAL LQUERY, UPPER
+ INTEGER K, KP, LWKOPT
* ..
* .. External Functions ..
LOGICAL LSAME
*
INFO = 0
UPPER = LSAME( UPLO, 'U' )
+ LQUERY = ( LWORK.EQ.-1 )
IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -5
ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
INFO = -8
- ELSE IF( LWORK.LT.(3*N-2) ) THEN
+ ELSE IF( LWORK.LT.(3*N-2) .AND. .NOT.LQUERY ) THEN
INFO = -10
END IF
IF( INFO.NE.0 ) THEN
CALL XERBLA( 'ZHETRS_AASEN', -INFO )
RETURN
+ ELSE IF( LQUERY ) THEN
+ LWKOPT = (3*N-2)
+ WORK( 1 ) = LWKOPT
+ RETURN
END IF
*
* Quick return if possible
$ B( 2, 1 ), LDB)
*
* Compute T \ B -> B [ T \ (U \P**T * B) ]
-*
+*
CALL ZLACPY( 'F', 1, N, A(1, 1), LDA+1, WORK(N), 1)
IF( N.GT.1 ) THEN
CALL ZLACPY( 'F', 1, N-1, A( 1, 2 ), LDA+1, WORK( 2*N ), 1)
END IF
CALL ZGTSV(N, NRHS, WORK(1), WORK(N), WORK(2*N), B, LDB,
$ INFO)
-*
+*
* Compute (U**T \ B) -> B [ U**T \ (T \ (U \P**T * B) ) ]
*
CALL ZTRSM( 'L', 'U', 'N', 'U', N-1, NRHS, ONE, A( 1, 2 ), LDA,
END IF
CALL ZGTSV(N, NRHS, WORK(1), WORK(N), WORK(2*N), B, LDB,
$ INFO)
-*
+*
* Compute (L**T \ B) -> B [ L**T \ (T \ (L \P**T * B) ) ]
-*
+*
CALL ZTRSM( 'L', 'L', 'C', 'U', N-1, NRHS, ONE, A( 2, 1 ), LDA,
$ B( 2, 1 ), LDB)
*
*
* =========== DOCUMENTATION ===========
*
-* Online html documentation available at
-* http://www.netlib.org/lapack/explore-html/
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
-*> Download ZLAHEF_AASEN + dependencies
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlahef_aasen.f">
-*> [TGZ]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlahef_aasen.f">
-*> [ZIP]</a>
-*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlahef_aasen.f">
+*> Download ZLAHEF_AASEN + dependencies
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlahef_aasen.f">
+*> [TGZ]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlahef_aasen.f">
+*> [ZIP]</a>
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlahef_aasen.f">
*> [TXT]</a>
-*> \endhtmlonly
+*> \endhtmlonly
*
* Definition:
* ===========
*
-* SUBROUTINE ZLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+* SUBROUTINE ZLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
* H, LDH, WORK, INFO )
-*
+*
* .. Scalar Arguments ..
* CHARACTER UPLO
* INTEGER J1, M, NB, LDA, LDH, INFO
* INTEGER IPIV( * )
* COMPLEX*16 A( LDA, * ), H( LDH, * ), WORK( * )
* ..
-*
+*
*
*> \par Purpose:
* =============
*> last row, or column, of the previous panel. The first row, or column,
*> of A is set to be the first row, or column, of an identity matrix,
*> which is used to factorize the first panel.
-*>
+*>
*> The resulting J-th row of U, or J-th column of L, is stored in the
-*> (J-1)-th row, or column, of A (without the unit diatonals), while
+*> (J-1)-th row, or column, of A (without the unit diatonals), while
*> the diagonal and subdiagonal of A are overwritten by those of T.
*>
*> \endverbatim
* Authors:
* ========
*
-*> \author Univ. of Tennessee
-*> \author Univ. of California Berkeley
-*> \author Univ. of Colorado Denver
-*> \author NAG Ltd.
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
*> \date November 2016
*
-*> \ingroup complex16SYcomputational
-*
-* @precisions fortran z -> c
+*> \ingroup complex16HEcomputational
*
* =====================================================================
- SUBROUTINE ZLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
+ SUBROUTINE ZLAHEF_AASEN( UPLO, J1, M, NB, A, LDA, IPIV,
$ H, LDH, WORK, INFO )
*
-* -- LAPACK computational routine (version 3.4.0) --
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* November 2016
*
* .. Local Scalars ..
INTEGER J, K, K1, I1, I2
- COMPLEX*16 PIV, ALPHA
+ COMPLEX*16 PIV, ALPHA
* ..
* .. External Functions ..
LOGICAL LSAME
*
A( K, J ) = DBLE( WORK( 1 ) )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L(J, (J+1):N)
* where A(J, J) stores T(J, J) and A(J-1, (J+1):N) stores U(J, (J+1):N)
*
IF( (J1+J-1).GT.1 ) THEN
- ALPHA = -A( K, J )
- CALL ZAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
+ ALPHA = -A( K, J )
+ CALL ZAXPY( M-J, ALPHA, A( K-1, J+1 ), LDA,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL ZSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
+ CALL ZSWAP( I2-I1-1, A( J1+I1-1, I1+1 ), LDA,
$ A( J1+I1, I2 ), 1 )
CALL ZLACGV( I2-I1, A( J1+I1-1, I1+1 ), LDA )
CALL ZLACGV( I2-I1-1, A( J1+I1, I2 ), 1 )
*
* Swap A(I1, I2+1:N) with A(I2, I2+1:N)
*
- CALL ZSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
+ CALL ZSWAP( M-I2, A( J1+I1-1, I2+1 ), LDA,
$ A( J1+I2-1, I2+1 ), LDA )
*
* Swap A(I1, I1) with A(I2,I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL ZSWAP( I1-K1+1, A( 1, I1 ), 1,
+ CALL ZSWAP( I1-K1+1, A( 1, I1 ), 1,
$ A( 1, I2 ), 1 )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J, J+1) = T(J, J+1)
*
A( K, J+1 ) = WORK( 2 )
- IF( (A( K, J ).EQ.ZERO ) .AND.
+ IF( (A( K, J ).EQ.ZERO ) .AND.
$ ( (J.EQ.M) .OR. (A( K, J+1 ).EQ.ZERO))) THEN
IF(INFO .EQ. 0) THEN
INFO = J
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J:N, J),
+* Copy A(J+1:N, J+1) into H(J:N, J),
*
- CALL ZCOPY( M-J, A( K+1, J+1 ), LDA,
+ CALL ZCOPY( M-J, A( K+1, J+1 ), LDA,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL ZCOPY( M-J-1, WORK( 3 ), 1, A( K, J+2 ), LDA )
CALL ZSCAL( M-J-1, ALPHA, A( K, J+2 ), LDA )
ELSE
- CALL ZLASET( 'Full', 1, M-J-1, ZERO, ZERO,
+ CALL ZLASET( 'Full', 1, M-J-1, ZERO, ZERO,
$ A( K, J+2 ), LDA)
END IF
ELSE
*
A( J, K ) = DBLE( WORK( 1 ) )
*
- IF( J.LT.M ) THEN
+ IF( J.LT.M ) THEN
*
* Compute WORK(2:N) = T(J, J) L((J+1):N, J)
* where A(J, J) = T(J, J) and A((J+1):N, J-1) = L((J+1):N, J)
*
IF( (J1+J-1).GT.1 ) THEN
ALPHA = -A( J, K )
- CALL ZAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
+ CALL ZAXPY( M-J, ALPHA, A( J+1, K-1 ), 1,
$ WORK( 2 ), 1 )
ENDIF
*
*
I1 = I1+J-1
I2 = I2+J-1
- CALL ZSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
+ CALL ZSWAP( I2-I1-1, A( I1+1, J1+I1-1 ), 1,
$ A( I2, J1+I1 ), LDA )
CALL ZLACGV( I2-I1, A( I1+1, J1+I1-1 ), 1 )
CALL ZLACGV( I2-I1-1, A( I2, J1+I1 ), LDA )
*
* Swap A(I2+1:N, I1) with A(I2+1:N, I2)
*
- CALL ZSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
+ CALL ZSWAP( M-I2, A( I2+1, J1+I1-1 ), 1,
$ A( I2+1, J1+I2-1 ), 1 )
*
* Swap A(I1, I1) with A(I2, I2)
* Swap L(1:I1-1, I1) with L(1:I1-1, I2),
* skipping the first column
*
- CALL ZSWAP( I1-K1+1, A( I1, 1 ), LDA,
+ CALL ZSWAP( I1-K1+1, A( I1, 1 ), LDA,
$ A( I2, 1 ), LDA )
END IF
- ELSE
+ ELSE
IPIV( J+1 ) = J+1
ENDIF
*
* Set A(J+1, J) = T(J+1, J)
*
A( J+1, K ) = WORK( 2 )
- IF( (A( J, K ).EQ.ZERO) .AND.
+ IF( (A( J, K ).EQ.ZERO) .AND.
$ ( (J.EQ.M) .OR. (A( J+1, K ).EQ.ZERO)) ) THEN
- IF (INFO .EQ. 0)
+ IF (INFO .EQ. 0)
$ INFO = J
END IF
*
IF( J.LT.NB ) THEN
*
-* Copy A(J+1:N, J+1) into H(J+1:N, J),
+* Copy A(J+1:N, J+1) into H(J+1:N, J),
*
- CALL ZCOPY( M-J, A( J+1, K+1 ), 1,
+ CALL ZCOPY( M-J, A( J+1, K+1 ), 1,
$ H( J+1, J+1 ), 1 )
END IF
*
CALL ZCOPY( M-J-1, WORK( 3 ), 1, A( J+2, K ), 1 )
CALL ZSCAL( M-J-1, ALPHA, A( J+2, K ), 1 )
ELSE
- CALL ZLASET( 'Full', M-J-1, 1, ZERO, ZERO,
+ CALL ZLASET( 'Full', M-J-1, 1, ZERO, ZERO,
$ A( J+2, K ), LDA )
END IF
ELSE
- IF( (A( J, K ).EQ.ZERO) .AND. (J.EQ.M)
+ IF( (A( J, K ).EQ.ZERO) .AND. (J.EQ.M)
$ .AND. (INFO.EQ.0) ) INFO = J
END IF
J = J + 1
projects / platform / upstream / lapack.git / commitdiff