Return to zpptri.f CVS log

Up to [local] / rpl / lapack / lapack

Cohérence Lapack.

    1: *> \brief \b ZPPTRI
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at 
    6: *            http://www.netlib.org/lapack/explore-html/ 
    7: *
    8: *> \htmlonly
    9: *> Download ZPPTRI + dependencies 
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zpptri.f"> 
   11: *> [TGZ]</a> 
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zpptri.f"> 
   13: *> [ZIP]</a> 
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zpptri.f"> 
   15: *> [TXT]</a>
   16: *> \endhtmlonly 
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE ZPPTRI( UPLO, N, AP, INFO )
   22: * 
   23: *       .. Scalar Arguments ..
   24: *       CHARACTER          UPLO
   25: *       INTEGER            INFO, N
   26: *       ..
   27: *       .. Array Arguments ..
   28: *       COMPLEX*16         AP( * )
   29: *       ..
   30: *  
   31: *
   32: *> \par Purpose:
   33: *  =============
   34: *>
   35: *> \verbatim
   36: *>
   37: *> ZPPTRI computes the inverse of a complex Hermitian positive definite
   38: *> matrix A using the Cholesky factorization A = U**H*U or A = L*L**H
   39: *> computed by ZPPTRF.
   40: *> \endverbatim
   41: *
   42: *  Arguments:
   43: *  ==========
   44: *
   45: *> \param[in] UPLO
   46: *> \verbatim
   47: *>          UPLO is CHARACTER*1
   48: *>          = 'U':  Upper triangular factor is stored in AP;
   49: *>          = 'L':  Lower triangular factor is stored in AP.
   50: *> \endverbatim
   51: *>
   52: *> \param[in] N
   53: *> \verbatim
   54: *>          N is INTEGER
   55: *>          The order of the matrix A.  N >= 0.
   56: *> \endverbatim
   57: *>
   58: *> \param[in,out] AP
   59: *> \verbatim
   60: *>          AP is COMPLEX*16 array, dimension (N*(N+1)/2)
   61: *>          On entry, the triangular factor U or L from the Cholesky
   62: *>          factorization A = U**H*U or A = L*L**H, packed columnwise as
   63: *>          a linear array.  The j-th column of U or L is stored in the
   64: *>          array AP as follows:
   65: *>          if UPLO = 'U', AP(i + (j-1)*j/2) = U(i,j) for 1<=i<=j;
   66: *>          if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = L(i,j) for j<=i<=n.
   67: *>
   68: *>          On exit, the upper or lower triangle of the (Hermitian)
   69: *>          inverse of A, overwriting the input factor U or L.
   70: *> \endverbatim
   71: *>
   72: *> \param[out] INFO
   73: *> \verbatim
   74: *>          INFO is INTEGER
   75: *>          = 0:  successful exit
   76: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
   77: *>          > 0:  if INFO = i, the (i,i) element of the factor U or L is
   78: *>                zero, and the inverse could not be computed.
   79: *> \endverbatim
   80: *
   81: *  Authors:
   82: *  ========
   83: *
   84: *> \author Univ. of Tennessee 
   85: *> \author Univ. of California Berkeley 
   86: *> \author Univ. of Colorado Denver 
   87: *> \author NAG Ltd. 
   88: *
   89: *> \date November 2011
   90: *
   91: *> \ingroup complex16OTHERcomputational
   92: *
   93: *  =====================================================================
   94:       SUBROUTINE ZPPTRI( UPLO, N, AP, INFO )
   95: *
   96: *  -- LAPACK computational routine (version 3.4.0) --
   97: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
   98: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
   99: *     November 2011
  100: *
  101: *     .. Scalar Arguments ..
  102:       CHARACTER          UPLO
  103:       INTEGER            INFO, N
  104: *     ..
  105: *     .. Array Arguments ..
  106:       COMPLEX*16         AP( * )
  107: *     ..
  108: *
  109: *  =====================================================================
  110: *
  111: *     .. Parameters ..
  112:       DOUBLE PRECISION   ONE
  113:       PARAMETER          ( ONE = 1.0D+0 )
  114: *     ..
  115: *     .. Local Scalars ..
  116:       LOGICAL            UPPER
  117:       INTEGER            J, JC, JJ, JJN
  118:       DOUBLE PRECISION   AJJ
  119: *     ..
  120: *     .. External Functions ..
  121:       LOGICAL            LSAME
  122:       COMPLEX*16         ZDOTC
  123:       EXTERNAL           LSAME, ZDOTC
  124: *     ..
  125: *     .. External Subroutines ..
  126:       EXTERNAL           XERBLA, ZDSCAL, ZHPR, ZTPMV, ZTPTRI
  127: *     ..
  128: *     .. Intrinsic Functions ..
  129:       INTRINSIC          DBLE
  130: *     ..
  131: *     .. Executable Statements ..
  132: *
  133: *     Test the input parameters.
  134: *
  135:       INFO = 0
  136:       UPPER = LSAME( UPLO, 'U' )
  137:       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
  138:          INFO = -1
  139:       ELSE IF( N.LT.0 ) THEN
  140:          INFO = -2
  141:       END IF
  142:       IF( INFO.NE.0 ) THEN
  143:          CALL XERBLA( 'ZPPTRI', -INFO )
  144:          RETURN
  145:       END IF
  146: *
  147: *     Quick return if possible
  148: *
  149:       IF( N.EQ.0 )
  150:      $   RETURN
  151: *
  152: *     Invert the triangular Cholesky factor U or L.
  153: *
  154:       CALL ZTPTRI( UPLO, 'Non-unit', N, AP, INFO )
  155:       IF( INFO.GT.0 )
  156:      $   RETURN
  157:       IF( UPPER ) THEN
  158: *
  159: *        Compute the product inv(U) * inv(U)**H.
  160: *
  161:          JJ = 0
  162:          DO 10 J = 1, N
  163:             JC = JJ + 1
  164:             JJ = JJ + J
  165:             IF( J.GT.1 )
  166:      $         CALL ZHPR( 'Upper', J-1, ONE, AP( JC ), 1, AP )
  167:             AJJ = AP( JJ )
  168:             CALL ZDSCAL( J, AJJ, AP( JC ), 1 )
  169:    10    CONTINUE
  170: *
  171:       ELSE
  172: *
  173: *        Compute the product inv(L)**H * inv(L).
  174: *
  175:          JJ = 1
  176:          DO 20 J = 1, N
  177:             JJN = JJ + N - J + 1
  178:             AP( JJ ) = DBLE( ZDOTC( N-J+1, AP( JJ ), 1, AP( JJ ), 1 ) )
  179:             IF( J.LT.N )
  180:      $         CALL ZTPMV( 'Lower', 'Conjugate transpose', 'Non-unit',
  181:      $                     N-J, AP( JJN ), AP( JJ+1 ), 1 )
  182:             JJ = JJN
  183:    20    CONTINUE
  184:       END IF
  185: *
  186:       RETURN
  187: *
  188: *     End of ZPPTRI
  189: *
  190:       END