LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ spptrf()

subroutine spptrf ( character uplo,
integer n,
real, dimension( * ) ap,
integer info )

SPPTRF

Download SPPTRF + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> SPPTRF computes the Cholesky factorization of a real symmetric
!> positive definite matrix A stored in packed format.
!>
!> The factorization has the form
!>    A = U**T * U,  if UPLO = 'U', or
!>    A = L  * L**T,  if UPLO = 'L',
!> where U is an upper triangular matrix and L is lower triangular.
!> 
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>          = 'U':  Upper triangle of A is stored;
!>          = 'L':  Lower triangle of A is stored.
!> 
[in]N
!>          N is INTEGER
!>          The order of the matrix A.  N >= 0.
!> 
[in,out]AP
!>          AP is REAL array, dimension (N*(N+1)/2)
!>          On entry, the upper or lower triangle of the symmetric matrix
!>          A, packed columnwise in a linear array.  The j-th column of A
!>          is stored in the array AP as follows:
!>          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
!>          if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n.
!>          See below for further details.
!>
!>          On exit, if INFO = 0, the triangular factor U or L from the
!>          Cholesky factorization A = U**T*U or A = L*L**T, in the same
!>          storage format as A.
!> 
[out]INFO
!>          INFO is INTEGER
!>          = 0:  successful exit
!>          < 0:  if INFO = -i, the i-th argument had an illegal value
!>          > 0:  if INFO = i, the leading principal minor of order i
!>                is not positive, and the factorization could not be
!>                completed.
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  The packed storage scheme is illustrated by the following example
!>  when N = 4, UPLO = 'U':
!>
!>  Two-dimensional storage of the symmetric matrix A:
!>
!>     a11 a12 a13 a14
!>         a22 a23 a24
!>             a33 a34     (aij = aji)
!>                 a44
!>
!>  Packed storage of the upper triangle of A:
!>
!>  AP = [ a11, a12, a22, a13, a23, a33, a14, a24, a34, a44 ]
!> 

Definition at line 116 of file spptrf.f.

117*
118* -- LAPACK computational routine --
119* -- LAPACK is a software package provided by Univ. of Tennessee, --
120* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
121*
122* .. Scalar Arguments ..
123 CHARACTER UPLO
124 INTEGER INFO, N
125* ..
126* .. Array Arguments ..
127 REAL AP( * )
128* ..
129*
130* =====================================================================
131*
132* .. Parameters ..
133 REAL ONE, ZERO
134 parameter( one = 1.0e+0, zero = 0.0e+0 )
135* ..
136* .. Local Scalars ..
137 LOGICAL UPPER
138 INTEGER J, JC, JJ
139 REAL AJJ
140* ..
141* .. External Functions ..
142 LOGICAL LSAME
143 REAL SDOT
144 EXTERNAL lsame, sdot
145* ..
146* .. External Subroutines ..
147 EXTERNAL sscal, sspr, stpsv, xerbla
148* ..
149* .. Intrinsic Functions ..
150 INTRINSIC sqrt
151* ..
152* .. Executable Statements ..
153*
154* Test the input parameters.
155*
156 info = 0
157 upper = lsame( uplo, 'U' )
158 IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
159 info = -1
160 ELSE IF( n.LT.0 ) THEN
161 info = -2
162 END IF
163 IF( info.NE.0 ) THEN
164 CALL xerbla( 'SPPTRF', -info )
165 RETURN
166 END IF
167*
168* Quick return if possible
169*
170 IF( n.EQ.0 )
171 $ RETURN
172*
173 IF( upper ) THEN
174*
175* Compute the Cholesky factorization A = U**T*U.
176*
177 jj = 0
178 DO 10 j = 1, n
179 jc = jj + 1
180 jj = jj + j
181*
182* Compute elements 1:J-1 of column J.
183*
184 IF( j.GT.1 )
185 $ CALL stpsv( 'Upper', 'Transpose', 'Non-unit', j-1, ap,
186 $ ap( jc ), 1 )
187*
188* Compute U(J,J) and test for non-positive-definiteness.
189*
190 ajj = ap( jj ) - sdot( j-1, ap( jc ), 1, ap( jc ), 1 )
191 IF( ajj.LE.zero ) THEN
192 ap( jj ) = ajj
193 GO TO 30
194 END IF
195 ap( jj ) = sqrt( ajj )
196 10 CONTINUE
197 ELSE
198*
199* Compute the Cholesky factorization A = L*L**T.
200*
201 jj = 1
202 DO 20 j = 1, n
203*
204* Compute L(J,J) and test for non-positive-definiteness.
205*
206 ajj = ap( jj )
207 IF( ajj.LE.zero ) THEN
208 ap( jj ) = ajj
209 GO TO 30
210 END IF
211 ajj = sqrt( ajj )
212 ap( jj ) = ajj
213*
214* Compute elements J+1:N of column J and update the trailing
215* submatrix.
216*
217 IF( j.LT.n ) THEN
218 CALL sscal( n-j, one / ajj, ap( jj+1 ), 1 )
219 CALL sspr( 'Lower', n-j, -one, ap( jj+1 ), 1,
220 $ ap( jj+n-j+1 ) )
221 jj = jj + n - j + 1
222 END IF
223 20 CONTINUE
224 END IF
225 GO TO 40
226*
227 30 CONTINUE
228 info = j
229*
230 40 CONTINUE
231 RETURN
232*
233* End of SPPTRF
234*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
real function sdot(n, sx, incx, sy, incy)
SDOT
Definition sdot.f:82
subroutine sspr(uplo, n, alpha, x, incx, ap)
SSPR
Definition sspr.f:127
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine sscal(n, sa, sx, incx)
SSCAL
Definition sscal.f:79
subroutine stpsv(uplo, trans, diag, n, ap, x, incx)
STPSV
Definition stpsv.f:144
Here is the call graph for this function:
Here is the caller graph for this function: