#include "ESMF_LapackBlas.inc" !> \brief \b DTRSM ! ! =========== DOCUMENTATION =========== ! ! Online html documentation available at ! http://www.netlib.org/lapack/explore-html/ ! ! Definition: ! =========== ! ! SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB) ! ! .. Scalar Arguments .. ! DOUBLE PRECISION ALPHA ! INTEGER LDA,LDB,M,N ! CHARACTER DIAG,SIDE,TRANSA,UPLO ! .. ! .. Array Arguments .. ! DOUBLE PRECISION A(LDA,*),B(LDB,*) ! .. ! ! !> \par Purpose: ! ============= !> !> \verbatim !> !> DTRSM solves one of the matrix equations !> !> op( A )*X = alpha*B, or X*op( A ) = alpha*B, !> !> where alpha is a scalar, X and B are m by n matrices, A is a unit, or !> non-unit, upper or lower triangular matrix and op( A ) is one of !> !> op( A ) = A or op( A ) = A**T. !> !> The matrix X is overwritten on B. !> \endverbatim ! ! Arguments: ! ========== ! !> \param[in] SIDE !> \verbatim !> SIDE is CHARACTER*1 !> On entry, SIDE specifies whether op( A ) appears on the left !> or right of X as follows: !> !> SIDE = 'L' or 'l' op( A )*X = alpha*B. !> !> SIDE = 'R' or 'r' X*op( A ) = alpha*B. !> \endverbatim !> !> \param[in] UPLO !> \verbatim !> UPLO is CHARACTER*1 !> On entry, UPLO specifies whether the matrix A is an upper or !> lower triangular matrix as follows: !> !> UPLO = 'U' or 'u' A is an upper triangular matrix. !> !> UPLO = 'L' or 'l' A is a lower triangular matrix. !> \endverbatim !> !> \param[in] TRANSA !> \verbatim !> TRANSA is CHARACTER*1 !> On entry, TRANSA specifies the form of op( A ) to be used in !> the matrix multiplication as follows: !> !> TRANSA = 'N' or 'n' op( A ) = A. !> !> TRANSA = 'T' or 't' op( A ) = A**T. !> !> TRANSA = 'C' or 'c' op( A ) = A**T. !> \endverbatim !> !> \param[in] DIAG !> \verbatim !> DIAG is CHARACTER*1 !> On entry, DIAG specifies whether or not A is unit triangular !> as follows: !> !> DIAG = 'U' or 'u' A is assumed to be unit triangular. !> !> DIAG = 'N' or 'n' A is not assumed to be unit !> triangular. !> \endverbatim !> !> \param[in] M !> \verbatim !> M is INTEGER !> On entry, M specifies the number of rows of B. M must be at !> least zero. !> \endverbatim !> !> \param[in] N !> \verbatim !> N is INTEGER !> On entry, N specifies the number of columns of B. N must be !> at least zero. !> \endverbatim !> !> \param[in] ALPHA !> \verbatim !> ALPHA is DOUBLE PRECISION. !> On entry, ALPHA specifies the scalar alpha. When alpha is !> zero then A is not referenced and B need not be set before !> entry. !> \endverbatim !> !> \param[in] A !> \verbatim !> A is DOUBLE PRECISION array of DIMENSION ( LDA, k ), !> where k is m when SIDE = 'L' or 'l' !> and k is n when SIDE = 'R' or 'r'. !> Before entry with UPLO = 'U' or 'u', the leading k by k !> upper triangular part of the array A must contain the upper !> triangular matrix and the strictly lower triangular part of !> A is not referenced. !> Before entry with UPLO = 'L' or 'l', the leading k by k !> lower triangular part of the array A must contain the lower !> triangular matrix and the strictly upper triangular part of !> A is not referenced. !> Note that when DIAG = 'U' or 'u', the diagonal elements of !> A are not referenced either, but are assumed to be unity. !> \endverbatim !> !> \param[in] LDA !> \verbatim !> LDA is INTEGER !> On entry, LDA specifies the first dimension of A as declared !> in the calling (sub) program. When SIDE = 'L' or 'l' then !> LDA must be at least max( 1, m ), when SIDE = 'R' or 'r' !> then LDA must be at least max( 1, n ). !> \endverbatim !> !> \param[in,out] B !> \verbatim !> B is DOUBLE PRECISION array of DIMENSION ( LDB, n ). !> Before entry, the leading m by n part of the array B must !> contain the right-hand side matrix B, and on exit is !> overwritten by the solution matrix X. !> \endverbatim !> !> \param[in] LDB !> \verbatim !> LDB is INTEGER !> On entry, LDB specifies the first dimension of B as declared !> in the calling (sub) program. LDB must be at least !> max( 1, m ). !> \endverbatim ! ! Authors: ! ======== ! !> \author Univ. of Tennessee !> \author Univ. of California Berkeley !> \author Univ. of Colorado Denver !> \author NAG Ltd. ! !> \date November 2011 ! !> \ingroup double_blas_level3 ! !> \par Further Details: ! ===================== !> !> \verbatim !> !> Level 3 Blas routine. !> !> !> -- Written on 8-February-1989. !> Jack Dongarra, Argonne National Laboratory. !> Iain Duff, AERE Harwell. !> Jeremy Du Croz, Numerical Algorithms Group Ltd. !> Sven Hammarling, Numerical Algorithms Group Ltd. !> \endverbatim !> ! ===================================================================== SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB) ! ! -- Reference BLAS level3 routine (version 3.4.0) -- ! -- Reference BLAS is a software package provided by Univ. of Tennessee, -- ! -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- ! November 2011 ! ! .. Scalar Arguments .. DOUBLE PRECISION ALPHA INTEGER LDA,LDB,M,N CHARACTER DIAG,SIDE,TRANSA,UPLO ! .. ! .. Array Arguments .. DOUBLE PRECISION A(LDA,*),B(LDB,*) ! .. ! ! ===================================================================== ! ! .. External Functions .. LOGICAL LSAME EXTERNAL LSAME ! .. ! .. External Subroutines .. EXTERNAL XERBLA ! .. ! .. Intrinsic Functions .. INTRINSIC MAX ! .. ! .. Local Scalars .. DOUBLE PRECISION TEMP INTEGER I,INFO,J,K,NROWA LOGICAL LSIDE,NOUNIT,UPPER ! .. ! .. Parameters .. DOUBLE PRECISION ONE,ZERO PARAMETER (ONE=1.0D+0,ZERO=0.0D+0) ! .. ! ! Test the input parameters. ! LSIDE = LSAME(SIDE,'L') IF (LSIDE) THEN NROWA = M ELSE NROWA = N END IF NOUNIT = LSAME(DIAG,'N') UPPER = LSAME(UPLO,'U') ! INFO = 0 IF ((.NOT.LSIDE) .AND. (.NOT.LSAME(SIDE,'R'))) THEN INFO = 1 ELSE IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN INFO = 2 ELSE IF ((.NOT.LSAME(TRANSA,'N')) .AND. & & (.NOT.LSAME(TRANSA,'T')) .AND. & & (.NOT.LSAME(TRANSA,'C'))) THEN INFO = 3 ELSE IF ((.NOT.LSAME(DIAG,'U')) .AND. (.NOT.LSAME(DIAG,'N'))) THEN INFO = 4 ELSE IF (M.LT.0) THEN INFO = 5 ELSE IF (N.LT.0) THEN INFO = 6 ELSE IF (LDA.LT.MAX(1,NROWA)) THEN INFO = 9 ELSE IF (LDB.LT.MAX(1,M)) THEN INFO = 11 END IF IF (INFO.NE.0) THEN CALL XERBLA('DTRSM ',INFO) RETURN END IF ! ! Quick return if possible. ! IF (M.EQ.0 .OR. N.EQ.0) RETURN ! ! And when alpha.eq.zero. ! IF (ALPHA.EQ.ZERO) THEN DO 20 J = 1,N DO 10 I = 1,M B(I,J) = ZERO 10 CONTINUE 20 CONTINUE RETURN END IF ! ! Start the operations. ! IF (LSIDE) THEN IF (LSAME(TRANSA,'N')) THEN ! ! Form B := alpha*inv( A )*B. ! IF (UPPER) THEN DO 60 J = 1,N IF (ALPHA.NE.ONE) THEN DO 30 I = 1,M B(I,J) = ALPHA*B(I,J) 30 CONTINUE END IF DO 50 K = M,1,-1 IF (B(K,J).NE.ZERO) THEN IF (NOUNIT) B(K,J) = B(K,J)/A(K,K) DO 40 I = 1,K - 1 B(I,J) = B(I,J) - B(K,J)*A(I,K) 40 CONTINUE END IF 50 CONTINUE 60 CONTINUE ELSE DO 100 J = 1,N IF (ALPHA.NE.ONE) THEN DO 70 I = 1,M B(I,J) = ALPHA*B(I,J) 70 CONTINUE END IF DO 90 K = 1,M IF (B(K,J).NE.ZERO) THEN IF (NOUNIT) B(K,J) = B(K,J)/A(K,K) DO 80 I = K + 1,M B(I,J) = B(I,J) - B(K,J)*A(I,K) 80 CONTINUE END IF 90 CONTINUE 100 CONTINUE END IF ELSE ! ! Form B := alpha*inv( A**T )*B. ! IF (UPPER) THEN DO 130 J = 1,N DO 120 I = 1,M TEMP = ALPHA*B(I,J) DO 110 K = 1,I - 1 TEMP = TEMP - A(K,I)*B(K,J) 110 CONTINUE IF (NOUNIT) TEMP = TEMP/A(I,I) B(I,J) = TEMP 120 CONTINUE 130 CONTINUE ELSE DO 160 J = 1,N DO 150 I = M,1,-1 TEMP = ALPHA*B(I,J) DO 140 K = I + 1,M TEMP = TEMP - A(K,I)*B(K,J) 140 CONTINUE IF (NOUNIT) TEMP = TEMP/A(I,I) B(I,J) = TEMP 150 CONTINUE 160 CONTINUE END IF END IF ELSE IF (LSAME(TRANSA,'N')) THEN ! ! Form B := alpha*B*inv( A ). ! IF (UPPER) THEN DO 210 J = 1,N IF (ALPHA.NE.ONE) THEN DO 170 I = 1,M B(I,J) = ALPHA*B(I,J) 170 CONTINUE END IF DO 190 K = 1,J - 1 IF (A(K,J).NE.ZERO) THEN DO 180 I = 1,M B(I,J) = B(I,J) - A(K,J)*B(I,K) 180 CONTINUE END IF 190 CONTINUE IF (NOUNIT) THEN TEMP = ONE/A(J,J) DO 200 I = 1,M B(I,J) = TEMP*B(I,J) 200 CONTINUE END IF 210 CONTINUE ELSE DO 260 J = N,1,-1 IF (ALPHA.NE.ONE) THEN DO 220 I = 1,M B(I,J) = ALPHA*B(I,J) 220 CONTINUE END IF DO 240 K = J + 1,N IF (A(K,J).NE.ZERO) THEN DO 230 I = 1,M B(I,J) = B(I,J) - A(K,J)*B(I,K) 230 CONTINUE END IF 240 CONTINUE IF (NOUNIT) THEN TEMP = ONE/A(J,J) DO 250 I = 1,M B(I,J) = TEMP*B(I,J) 250 CONTINUE END IF 260 CONTINUE END IF ELSE ! ! Form B := alpha*B*inv( A**T ). ! IF (UPPER) THEN DO 310 K = N,1,-1 IF (NOUNIT) THEN TEMP = ONE/A(K,K) DO 270 I = 1,M B(I,K) = TEMP*B(I,K) 270 CONTINUE END IF DO 290 J = 1,K - 1 IF (A(J,K).NE.ZERO) THEN TEMP = A(J,K) DO 280 I = 1,M B(I,J) = B(I,J) - TEMP*B(I,K) 280 CONTINUE END IF 290 CONTINUE IF (ALPHA.NE.ONE) THEN DO 300 I = 1,M B(I,K) = ALPHA*B(I,K) 300 CONTINUE END IF 310 CONTINUE ELSE DO 360 K = 1,N IF (NOUNIT) THEN TEMP = ONE/A(K,K) DO 320 I = 1,M B(I,K) = TEMP*B(I,K) 320 CONTINUE END IF DO 340 J = K + 1,N IF (A(J,K).NE.ZERO) THEN TEMP = A(J,K) DO 330 I = 1,M B(I,J) = B(I,J) - TEMP*B(I,K) 330 CONTINUE END IF 340 CONTINUE IF (ALPHA.NE.ONE) THEN DO 350 I = 1,M B(I,K) = ALPHA*B(I,K) 350 CONTINUE END IF 360 CONTINUE END IF END IF END IF ! RETURN ! ! End of DTRSM . ! END