3 * =========== DOCUMENTATION ===========
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
9 *> Download CPPTRF + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cpptrf.f">
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cpptrf.f">
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cpptrf.f">
21 * SUBROUTINE CPPTRF( UPLO, N, AP, INFO )
23 * .. Scalar Arguments ..
27 * .. Array Arguments ..
37 *> CPPTRF computes the Cholesky factorization of a complex Hermitian
38 *> positive definite matrix A stored in packed format.
40 *> The factorization has the form
41 *> A = U**H * U, if UPLO = 'U', or
42 *> A = L * L**H, if UPLO = 'L',
43 *> where U is an upper triangular matrix and L is lower triangular.
51 *> UPLO is CHARACTER*1
52 *> = 'U': Upper triangle of A is stored;
53 *> = 'L': Lower triangle of A is stored.
59 *> The order of the matrix A. N >= 0.
64 *> AP is COMPLEX array, dimension (N*(N+1)/2)
65 *> On entry, the upper or lower triangle of the Hermitian matrix
66 *> A, packed columnwise in a linear array. The j-th column of A
67 *> is stored in the array AP as follows:
68 *> if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
69 *> if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n.
70 *> See below for further details.
72 *> On exit, if INFO = 0, the triangular factor U or L from the
73 *> Cholesky factorization A = U**H*U or A = L*L**H, in the same
74 *> storage format as A.
80 *> = 0: successful exit
81 *> < 0: if INFO = -i, the i-th argument had an illegal value
82 *> > 0: if INFO = i, the leading minor of order i is not
83 *> positive definite, and the factorization could not be
90 *> \author Univ. of Tennessee
91 *> \author Univ. of California Berkeley
92 *> \author Univ. of Colorado Denver
95 *> \date November 2011
97 *> \ingroup complexOTHERcomputational
99 *> \par Further Details:
100 * =====================
104 *> The packed storage scheme is illustrated by the following example
105 *> when N = 4, UPLO = 'U':
107 *> Two-dimensional storage of the Hermitian matrix A:
111 *> a33 a34 (aij = conjg(aji))
114 *> Packed storage of the upper triangle of A:
116 *> AP = [ a11, a12, a22, a13, a23, a33, a14, a24, a34, a44 ]
119 * =====================================================================
120 SUBROUTINE CPPTRF( UPLO, N, AP, INFO )
122 * -- LAPACK computational routine (version 3.4.0) --
123 * -- LAPACK is a software package provided by Univ. of Tennessee, --
124 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
127 * .. Scalar Arguments ..
131 * .. Array Arguments ..
135 * =====================================================================
139 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 )
141 * .. Local Scalars ..
146 * .. External Functions ..
149 EXTERNAL LSAME, CDOTC
151 * .. External Subroutines ..
152 EXTERNAL CHPR, CSSCAL, CTPSV, XERBLA
154 * .. Intrinsic Functions ..
157 * .. Executable Statements ..
159 * Test the input parameters.
162 UPPER = LSAME( UPLO, 'U' )
163 IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
165 ELSE IF( N.LT.0 ) THEN
169 CALL XERBLA( 'CPPTRF', -INFO )
173 * Quick return if possible
180 * Compute the Cholesky factorization A = U**H * U.
187 * Compute elements 1:J-1 of column J.
190 $ CALL CTPSV( 'Upper', 'Conjugate transpose', 'Non-unit',
191 $ J-1, AP, AP( JC ), 1 )
193 * Compute U(J,J) and test for non-positive-definiteness.
195 AJJ = REAL( AP( JJ ) ) - CDOTC( J-1, AP( JC ), 1, AP( JC ),
197 IF( AJJ.LE.ZERO ) THEN
201 AP( JJ ) = SQRT( AJJ )
205 * Compute the Cholesky factorization A = L * L**H.
210 * Compute L(J,J) and test for non-positive-definiteness.
212 AJJ = REAL( AP( JJ ) )
213 IF( AJJ.LE.ZERO ) THEN
220 * Compute elements J+1:N of column J and update the trailing
224 CALL CSSCAL( N-J, ONE / AJJ, AP( JJ+1 ), 1 )
225 CALL CHPR( 'Lower', N-J, -ONE, AP( JJ+1 ), 1,