1*> \brief \b SPOTRF2
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 SPOTRF2( UPLO, N, A, LDA, INFO )
12*
13*       .. Scalar Arguments ..
14*       CHARACTER          UPLO
15*       INTEGER            INFO, LDA, N
16*       ..
17*       .. Array Arguments ..
18*       REAL               A( LDA, * )
19*       ..
20*
21*
22*> \par Purpose:
23*  =============
24*>
25*> \verbatim
26*>
27*> SPOTRF2 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**T * U,  if UPLO = 'U', or
32*>    A = L  * L**T,  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 REAL 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**T*U or A = L*L**T.
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 realPOcomputational
105*
106*  =====================================================================
107      RECURSIVE SUBROUTINE SPOTRF2( 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      REAL               A( LDA, * )
120*     ..
121*
122*  =====================================================================
123*
124*     .. Parameters ..
125      REAL               ONE, ZERO
126      PARAMETER          ( ONE = 1.0E+0, ZERO=0.0E+0 )
127*     ..
128*     .. Local Scalars ..
129      LOGICAL            UPPER
130      INTEGER            N1, N2, IINFO
131*     ..
132*     .. External Functions ..
133      LOGICAL            LSAME, SISNAN
134      EXTERNAL           LSAME, SISNAN
135*     ..
136*     .. External Subroutines ..
137      EXTERNAL           SGEMM, SSYRK, XERBLA
138*     ..
139*     .. Intrinsic Functions ..
140      INTRINSIC          MAX, SQRT
141*     ..
142*     .. Executable Statements ..
143*
144*     Test the input parameters
145*
146      INFO = 0
147      UPPER = LSAME( UPLO, 'U' )
148      IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
149         INFO = -1
150      ELSE IF( N.LT.0 ) THEN
151         INFO = -2
152      ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
153         INFO = -4
154      END IF
155      IF( INFO.NE.0 ) THEN
156         CALL XERBLA( 'SPOTRF2', -INFO )
157         RETURN
158      END IF
159*
160*     Quick return if possible
161*
162      IF( N.EQ.0 )
163     $   RETURN
164*
165*     N=1 case
166*
167      IF( N.EQ.1 ) THEN
168*
169*        Test for non-positive-definiteness
170*
171         IF( A( 1, 1 ).LE.ZERO.OR.SISNAN( A( 1, 1 ) ) ) THEN
172            INFO = 1
173            RETURN
174         END IF
175*
176*        Factor
177*
178         A( 1, 1 ) = SQRT( A( 1, 1 ) )
179*
180*     Use recursive code
181*
182      ELSE
183         N1 = N/2
184         N2 = N-N1
185*
186*        Factor A11
187*
188         CALL SPOTRF2( UPLO, N1, A( 1, 1 ), LDA, IINFO )
189         IF ( IINFO.NE.0 ) THEN
190            INFO = IINFO
191            RETURN
192         END IF
193*
194*        Compute the Cholesky factorization A = U**T*U
195*
196         IF( UPPER ) THEN
197*
198*           Update and scale A12
199*
200            CALL STRSM( 'L', 'U', 'T', 'N', N1, N2, ONE,
201     $                  A( 1, 1 ), LDA, A( 1, N1+1 ), LDA )
202*
203*           Update and factor A22
204*
205            CALL SSYRK( UPLO, 'T', N2, N1, -ONE, A( 1, N1+1 ), LDA,
206     $                  ONE, A( N1+1, N1+1 ), LDA )
207            CALL SPOTRF2( UPLO, N2, A( N1+1, N1+1 ), LDA, IINFO )
208            IF ( IINFO.NE.0 ) THEN
209               INFO = IINFO + N1
210               RETURN
211            END IF
212*
213*        Compute the Cholesky factorization A = L*L**T
214*
215         ELSE
216*
217*           Update and scale A21
218*
219            CALL STRSM( 'R', 'L', 'T', 'N', N2, N1, ONE,
220     $                  A( 1, 1 ), LDA, A( N1+1, 1 ), LDA )
221*
222*           Update and factor A22
223*
224            CALL SSYRK( UPLO, 'N', N2, N1, -ONE, A( N1+1, 1 ), LDA,
225     $                  ONE, A( N1+1, N1+1 ), LDA )
226            CALL SPOTRF2( UPLO, N2, A( N1+1, N1+1 ), LDA, IINFO )
227            IF ( IINFO.NE.0 ) THEN
228               INFO = IINFO + N1
229               RETURN
230            END IF
231         END IF
232      END IF
233      RETURN
234*
235*     End of SPOTRF2
236*
237      END
238