Return to zpotrf2.f CVS log

Up to [local] / rpl / lapack / lapack

Cohérence Lapack.

    1: *> \brief \b ZPOTRF2
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at 
    6: *            http://www.netlib.org/lapack/explore-html/ 
    7: *
    8: *  Definition:
    9: *  ===========
   10: *
   11: *       RECURSIVE SUBROUTINE ZPOTRF2( UPLO, N, A, LDA, INFO )
   12: * 
   13: *       .. Scalar Arguments ..
   14: *       CHARACTER          UPLO
   15: *       INTEGER            INFO, LDA, N
   16: *       ..
   17: *       .. Array Arguments ..
   18: *       COMPLEX*16         A( LDA, * )
   19: *       ..
   20: *  
   21: *
   22: *> \par Purpose:
   23: *  =============
   24: *>
   25: *> \verbatim
   26: *>
   27: *> ZPOTRF2 computes the Cholesky factorization of a real symmetric
   28: *> positive definite matrix A using the recursive algorithm.
   29: *>
   30: *> The factorization has the form
   31: *>    A = U**H * U,  if UPLO = 'U', or
   32: *>    A = L  * L**H,  if UPLO = 'L',
   33: *> where U is an upper triangular matrix and L is lower triangular.
   34: *>
   35: *> This is the recursive version of the algorithm. It divides
   36: *> the matrix into four submatrices:
   37: *>
   38: *>        [  A11 | A12  ]  where A11 is n1 by n1 and A22 is n2 by n2
   39: *>    A = [ -----|----- ]  with n1 = n/2
   40: *>        [  A21 | A22  ]       n2 = n-n1
   41: *>
   42: *> The subroutine calls itself to factor A11. Update and scale A21
   43: *> or A12, update A22 then call itself to factor A22.
   44: *> 
   45: *> \endverbatim
   46: *
   47: *  Arguments:
   48: *  ==========
   49: *
   50: *> \param[in] UPLO
   51: *> \verbatim
   52: *>          UPLO is CHARACTER*1
   53: *>          = 'U':  Upper triangle of A is stored;
   54: *>          = 'L':  Lower triangle of A is stored.
   55: *> \endverbatim
   56: *>
   57: *> \param[in] N
   58: *> \verbatim
   59: *>          N is INTEGER
   60: *>          The order of the matrix A.  N >= 0.
   61: *> \endverbatim
   62: *>
   63: *> \param[in,out] A
   64: *> \verbatim
   65: *>          A is COMPLEX*16 array, dimension (LDA,N)
   66: *>          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
   67: *>          N-by-N upper triangular part of A contains the upper
   68: *>          triangular part of the matrix A, and the strictly lower
   69: *>          triangular part of A is not referenced.  If UPLO = 'L', the
   70: *>          leading N-by-N lower triangular part of A contains the lower
   71: *>          triangular part of the matrix A, and the strictly upper
   72: *>          triangular part of A is not referenced.
   73: *>
   74: *>          On exit, if INFO = 0, the factor U or L from the Cholesky
   75: *>          factorization A = U**H*U or A = L*L**H.
   76: *> \endverbatim
   77: *>
   78: *> \param[in] LDA
   79: *> \verbatim
   80: *>          LDA is INTEGER
   81: *>          The leading dimension of the array A.  LDA >= max(1,N).
   82: *> \endverbatim
   83: *>
   84: *> \param[out] INFO
   85: *> \verbatim
   86: *>          INFO is INTEGER
   87: *>          = 0:  successful exit
   88: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
   89: *>          > 0:  if INFO = i, the leading minor of order i is not
   90: *>                positive definite, and the factorization could not be
   91: *>                completed.
   92: *> \endverbatim
   93: *
   94: *  Authors:
   95: *  ========
   96: *
   97: *> \author Univ. of Tennessee 
   98: *> \author Univ. of California Berkeley 
   99: *> \author Univ. of Colorado Denver 
  100: *> \author NAG Ltd. 
  101: *
  102: *> \date November 2015
  103: *
  104: *> \ingroup complex16POcomputational
  105: *
  106: *  =====================================================================
  107:       RECURSIVE SUBROUTINE ZPOTRF2( UPLO, N, A, LDA, INFO )
  108: *
  109: *  -- LAPACK computational routine (version 3.6.0) --
  110: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  111: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  112: *     November 2015
  113: *
  114: *     .. Scalar Arguments ..
  115:       CHARACTER          UPLO
  116:       INTEGER            INFO, LDA, N
  117: *     ..
  118: *     .. Array Arguments ..
  119:       COMPLEX*16         A( LDA, * )
  120: *     ..
  121: *
  122: *  =====================================================================
  123: *
  124: *     .. Parameters ..
  125:       DOUBLE PRECISION   ONE, ZERO
  126:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
  127:       COMPLEX*16         CONE
  128:       PARAMETER          ( CONE = (1.0D+0, 0.0D+0) )
  129: *     ..
  130: *     .. Local Scalars ..
  131:       LOGICAL            UPPER            
  132:       INTEGER            N1, N2, IINFO
  133:       DOUBLE PRECISION   AJJ
  134: *     ..
  135: *     .. External Functions ..
  136:       LOGICAL            LSAME, DISNAN
  137:       EXTERNAL           LSAME, DISNAN
  138: *     ..
  139: *     .. External Subroutines ..
  140:       EXTERNAL           ZHERK, ZTRSM, XERBLA
  141: *     ..
  142: *     .. Intrinsic Functions ..
  143:       INTRINSIC          MAX, DBLE, SQRT
  144: *     ..
  145: *     .. Executable Statements ..
  146: *
  147: *     Test the input parameters
  148: *
  149:       INFO = 0
  150:       UPPER = LSAME( UPLO, 'U' )
  151:       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
  152:          INFO = -1
  153:       ELSE IF( N.LT.0 ) THEN
  154:          INFO = -2
  155:       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
  156:          INFO = -4
  157:       END IF
  158:       IF( INFO.NE.0 ) THEN
  159:          CALL XERBLA( 'ZPOTRF2', -INFO )
  160:          RETURN
  161:       END IF
  162: *
  163: *     Quick return if possible
  164: *
  165:       IF( N.EQ.0 )
  166:      $   RETURN
  167: *
  168: *     N=1 case
  169: *
  170:       IF( N.EQ.1 ) THEN
  171: *
  172: *        Test for non-positive-definiteness
  173: *
  174:          AJJ = DBLE( A( 1, 1 ) )
  175:          IF( AJJ.LE.ZERO.OR.DISNAN( AJJ ) ) THEN
  176:             INFO = 1
  177:             RETURN
  178:          END IF
  179: *
  180: *        Factor
  181: *
  182:          A( 1, 1 ) = SQRT( AJJ )
  183: *
  184: *     Use recursive code
  185: *
  186:       ELSE
  187:          N1 = N/2
  188:          N2 = N-N1
  189: *
  190: *        Factor A11
  191: *
  192:          CALL ZPOTRF2( UPLO, N1, A( 1, 1 ), LDA, IINFO )
  193:          IF ( IINFO.NE.0 ) THEN
  194:             INFO = IINFO
  195:             RETURN
  196:          END IF    
  197: *
  198: *        Compute the Cholesky factorization A = U**H*U
  199: *
  200:          IF( UPPER ) THEN
  201: *
  202: *           Update and scale A12
  203: *
  204:             CALL ZTRSM( 'L', 'U', 'C', 'N', N1, N2, CONE,
  205:      $                  A( 1, 1 ), LDA, A( 1, N1+1 ), LDA )
  206: *
  207: *           Update and factor A22
  208: *          
  209:             CALL ZHERK( UPLO, 'C', N2, N1, -ONE, A( 1, N1+1 ), LDA,
  210:      $                  ONE, A( N1+1, N1+1 ), LDA )
  211:             CALL ZPOTRF2( UPLO, N2, A( N1+1, N1+1 ), LDA, IINFO )
  212:             IF ( IINFO.NE.0 ) THEN
  213:                INFO = IINFO + N1
  214:                RETURN
  215:             END IF
  216: *
  217: *        Compute the Cholesky factorization A = L*L**H
  218: *
  219:          ELSE
  220: *
  221: *           Update and scale A21
  222: *
  223:             CALL ZTRSM( 'R', 'L', 'C', 'N', N2, N1, CONE,
  224:      $                  A( 1, 1 ), LDA, A( N1+1, 1 ), LDA )
  225: *
  226: *           Update and factor A22
  227: *
  228:             CALL ZHERK( UPLO, 'N', N2, N1, -ONE, A( N1+1, 1 ), LDA,
  229:      $                  ONE, A( N1+1, N1+1 ), LDA )
  230:             CALL ZPOTRF2( UPLO, N2, A( N1+1, N1+1 ), LDA, IINFO )
  231:             IF ( IINFO.NE.0 ) THEN
  232:                INFO = IINFO + N1
  233:                RETURN
  234:             END IF
  235:          END IF
  236:       END IF
  237:       RETURN
  238: *
  239: *     End of ZPOTRF2
  240: *
  241:       END