SUBROUTINE SLLTBI( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLLTBI computes an L*L**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked left-looking (i-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, IB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SLLTI, SSYRK, STRSM
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SLLTI( N, A, LDA, INFO )
      ELSE
         DO 10 I = 1, N, NB
            IB = MIN( NB, N-I+1 )
            CALL STRSM( 'Right', 'Lower', 'Transpose', 'Non-unit', IB,
     $                  I-1, ONE, A, LDA, A( I, 1 ), LDA )
            CALL SSYRK( 'Lower', 'No transpose', IB, I-1, -ONE,
     $                  A( I, 1 ), LDA, ONE, A( I, I ), LDA )
            CALL SLLTI( IB, A( I, I ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + I - 1
      RETURN
*
*     End of SLLTBI
*
      END
      SUBROUTINE SLLTI( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLLTI  computes an L*L**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked left-looking (i-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ZERO
      PARAMETER          ( ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           STRSV
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 I = 1, N
         CALL STRSV( 'Lower', 'No transpose', 'Non-unit', I-1, A, LDA,
     $               A( I, 1 ), LDA )
         A( I, I ) = A( I, I ) - SDOT( I-1, A( I, 1 ), LDA, A( I, 1 ),
     $               LDA )
         IF( A( I, I ).LE.ZERO )
     $      GO TO 20
         A( I, I ) = SQRT( A( I, I ) )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = I
      RETURN
*
*     End of SLLTI
*
      END
      SUBROUTINE SLLTBJ( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLLTBJ computes an L*L**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked top-looking (j-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J, JB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMM, SLLTJ, SSYRK, STRSM
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SLLTJ( N, A, LDA, INFO )
      ELSE
         DO 10 J = 1, N, NB
            JB = MIN( NB, N-J+1 )
            CALL SSYRK( 'Lower', 'No transpose', JB, J-1, -ONE,
     $                  A( J, 1 ), LDA, ONE, A( J, J ), LDA )
            CALL SLLTJ( JB, A( J, J ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
            IF( J+JB.LE.N ) THEN
               CALL SGEMM( 'No transpose', 'Transpose', N-J-JB+1, JB,
     $                     J-1, -ONE, A( J+JB, 1 ), LDA, A( J, 1 ), LDA,
     $                     ONE, A( J+JB, J ), LDA )
               CALL STRSM( 'Right', 'Lower', 'Transpose', 'Non-unit',
     $                     N-J-JB+1, JB, ONE, A( J, J ), LDA,
     $                     A( J+JB, J ), LDA )
            END IF
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + J - 1
      RETURN
*
*     End of SLLTBJ
*
      END
      SUBROUTINE SLLTJ( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLLTJ  computes an L*L**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked top-looking (j-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMV, SSCAL
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 J = 1, N
         A( J, J ) = A( J, J ) - SDOT( J-1, A( J, 1 ), LDA, A( J, 1 ),
     $               LDA )
         IF( A( J, J ).LE.ZERO )
     $      GO TO 20
         A( J, J ) = SQRT( A( J, J ) )
         IF( J.LT.N ) THEN
            CALL SGEMV( 'No transpose', N-J, J-1, -ONE, A( J+1, 1 ),
     $                  LDA, A( J, 1 ), LDA, ONE, A( J+1, J ), 1 )
            CALL SSCAL( N-J, ONE/A( J, J ), A( J+1, J ), 1 )
         END IF
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = J
      RETURN
*
*     End of SLLTJ
*
      END
      SUBROUTINE SLLTBK( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLLTBK computes an L*L**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked right-looking (k-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K, KB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SLLTK, SSYRK, STRSM
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SLLTK( N, A, LDA, INFO )
      ELSE
         DO 10 K = 1, N, NB
            KB = MIN( NB, N-K+1 )
            CALL SLLTK( KB, A( K, K ), LDA, INFO )
            IF( K+KB.LE.N ) THEN
               CALL STRSM( 'Right', 'Lower', 'Transpose', 'Non-unit',
     $                     N-K-KB+1, KB, ONE, A( K, K ), LDA,
     $                     A( K+KB, K ), LDA )
               CALL SSYRK( 'Lower', 'No transpose', N-K-KB+1, KB, -ONE,
     $                     A( K+KB, K ), LDA, ONE, A( K+KB, K+KB ),
     $                     LDA )
            END IF
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + K - 1
      RETURN
*
*     End of SLLTBK
*
      END
      SUBROUTINE SLLTK( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLLTK  computes an L*L**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked right-looking (k-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSCAL, SSYR
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 K = 1, N
         IF( A( K, K ).LE.ZERO )
     $      GO TO 20
         A( K, K ) = SQRT( A( K, K ) )
         IF( K.LT.N ) THEN
            CALL SSCAL( N-K, ONE/A( K, K ), A( K+1, K ), 1 )
            CALL SSYR( 'Lower', N-K, -ONE, A( K+1, K ), 1,
     $                 A( K+1, K+1 ), LDA )
         END IF
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = K
      RETURN
*
*     End of SLLTK
*
      END
      SUBROUTINE SLTLBI( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLTLBI computes an L**T*L factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked i-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, IB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMM, SLTLI, SSYRK, STRSM
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SLTLI( N, A, LDA, INFO )
      ELSE
         DO 10 I = ( ( N-1 )/NB )*NB + 1, 1, -NB
            IB = MIN( NB, N-I+1 )
            IF( I+IB.LE.N )
     $         CALL SSYRK( 'Lower', 'Transpose', IB, N-I-IB+1, -ONE,
     $                     A( I+IB, I ), LDA, ONE, A( I, I ), LDA )
            CALL SLTLI( IB, A( I, I ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
            IF( I+IB.LE.N )
     $         CALL SGEMM( 'Transpose', 'No transpose', IB, I-1,
     $                     N-I-IB+1, -ONE, A( I+IB, I ), LDA,
     $                     A( I+IB, 1 ), LDA, ONE, A( I, 1 ), LDA )
            CALL STRSM( 'Left', 'Lower', 'Transpose', 'Non-unit', IB,
     $                  I-1, ONE, A( I, I ), LDA, A( I, 1 ), LDA )
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + I - 1
      RETURN
*
*     End of SLTLBI
*
      END
      SUBROUTINE SLTLI( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLTLI  computes an L**T*L factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked i-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMV, SSCAL
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 I = N, 1, -1
         IF( I.LT.N )
     $      A( I, I ) = A( I, I ) - SDOT( N-I, A( I+1, I ), 1,
     $                  A( I+1, I ), 1 )
         IF( A( I, I ).LE.ZERO )
     $      GO TO 20
         A( I, I ) = SQRT( A( I, I ) )
         IF( I.LT.N )
     $      CALL SGEMV( 'Transpose', N-I, I-1, -ONE, A( I+1, 1 ), LDA,
     $                  A( I+1, I ), 1, ONE, A( I, 1 ), LDA )
         CALL SSCAL( I-1, ONE/A( I, I ), A( I, 1 ), LDA )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = I
      RETURN
*
*     End of SLTLI
*
      END
      SUBROUTINE SLTLBJ( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLTLBJ computes an L**T*L factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked j-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J, JB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SLTLJ, SSYRK, STRSM
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SLTLJ( N, A, LDA, INFO )
      ELSE
         DO 10 J = ( ( N-1 )/NB )*NB + 1, 1, -NB
            JB = MIN( NB, N-J+1 )
            IF( J+JB.LE.N ) THEN
               CALL STRSM( 'Left', 'Lower', 'Transpose', 'Non-unit',
     $                     N-J-JB+1, JB, ONE, A( J+JB, J+JB ), LDA,
     $                     A( J+JB, J ), LDA )
               CALL SSYRK( 'Lower', 'Transpose', JB, N-J-JB+1, -ONE,
     $                     A( J+JB, J ), LDA, ONE, A( J, J ), LDA )
            END IF
            CALL SLTLJ( JB, A( J, J ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + J - 1
      RETURN
*
*     End of SLTLBJ
*
      END
      SUBROUTINE SLTLJ( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLTLJ  computes an L**T*L factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked j-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ZERO
      PARAMETER          ( ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           STRSV
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 J = N, 1, -1
         IF( J.LT.N ) THEN
            CALL STRSV( 'Lower', 'Transpose', 'Non-unit', N-J,
     $                  A( J+1, J+1 ), LDA, A( J+1, J ), 1 )
            A( J, J ) = A( J, J ) - SDOT( N-J, A( J+1, J ), 1,
     $                  A( J+1, J ), 1 )
         END IF
         IF( A( J, J ).LE.ZERO )
     $      GO TO 20
         A( J, J ) = SQRT( A( J, J ) )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = J
      RETURN
*
*     End of SLTLJ
*
      END
      SUBROUTINE SLTLBK( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLTLBK computes an L**T*L factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked k-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K, KB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SLTLK, SSYRK, STRSM
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SLTLK( N, A, LDA, INFO )
      ELSE
         DO 10 K = ( ( N-1 )/NB )*NB + 1, 1, -NB
            KB = MIN( NB, N-K+1 )
            CALL SLTLK( KB, A( K, K ), LDA, INFO )
            CALL STRSM( 'Left', 'Lower', 'Transpose', 'Non-unit', KB,
     $                  K-1, ONE, A( K, K ), LDA, A( K, 1 ), LDA )
            CALL SSYRK( 'Lower', 'Transpose', K-1, KB, -ONE, A( K, 1 ),
     $                  LDA, ONE, A, LDA )
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + K - 1
      RETURN
*
*     End of SLTLBK
*
      END
      SUBROUTINE SLTLK( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SLTLK  computes an L**T*L factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked k-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the lower triangular
*          part of A is referenced.
*          On exit, the factor L from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSCAL, SSYR
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 K = N, 1, -1
         IF( A( K, K ).LE.ZERO )
     $      GO TO 20
         A( K, K ) = SQRT( A( K, K ) )
         CALL SSCAL( K-1, ONE/A( K, K ), A( K, 1 ), LDA )
         CALL SSYR( 'Lower', K-1, -ONE, A( K, 1 ), LDA, A, LDA )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = K
      RETURN
*
*     End of SLTLK
*
      END
      SUBROUTINE SUTUBI( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUTUBI computes a U**T*U factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked top-looking (i-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, IB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMM, SSYRK, STRSM, SUTUI
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SUTUI( N, A, LDA, INFO )
      ELSE
         DO 10 I = 1, N, NB
            IB = MIN( NB, N-I+1 )
            CALL SSYRK( 'Upper', 'Transpose', IB, I-1, -ONE, A( 1, I ),
     $                  LDA, ONE, A( I, I ), LDA )
            CALL SUTUI( IB, A( I, I ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
            IF( I+IB.LE.N ) THEN
               CALL SGEMM( 'Transpose', 'No transpose', IB, N-I-IB+1,
     $                     I-1, -ONE, A( 1, I ), LDA, A( 1, I+IB ), LDA,
     $                     ONE, A( I, I+IB ), LDA )
               CALL STRSM( 'Left', 'Upper', 'Transpose', 'Non-unit', IB,
     $                     N-I-IB+1, ONE, A( I, I ), LDA, A( I, I+IB ),
     $                     LDA )
            END IF
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + I - 1
      RETURN
*
*     End of SUTUBI
*
      END
      SUBROUTINE SUTUI( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUTUI  computes a U**T*U factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked top-looking (i-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMV, SSCAL
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 I = 1, N
         A( I, I ) = A( I, I ) - SDOT( I-1, A( 1, I ), 1, A( 1, I ), 1 )
         IF( A( I, I ).LE.ZERO )
     $      GO TO 20
         A( I, I ) = SQRT( A( I, I ) )
         IF( I.LT.N ) THEN
            CALL SGEMV( 'Transpose', I-1, N-I, -ONE, A( 1, I+1 ), LDA,
     $                  A( 1, I ), 1, ONE, A( I, I+1 ), LDA )
            CALL SSCAL( N-I, ONE/A( I, I ), A( I, I+1 ), LDA )
         END IF
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = I
      RETURN
*
*     End of SUTUI
*
      END
      SUBROUTINE SUTUBJ( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUTUBJ computes a U**T*U factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked left-looking (j-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J, JB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSYRK, STRSM, SUTUJ
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SUTUJ( N, A, LDA, INFO )
      ELSE
         DO 10 J = 1, N, NB
            JB = MIN( NB, N-J+1 )
            CALL STRSM( 'Left', 'Upper', 'Transpose', 'Non-unit', J-1,
     $                  JB, ONE, A, LDA, A( 1, J ), LDA )
            CALL SSYRK( 'Upper', 'Transpose', JB, J-1, -ONE, A( 1, J ),
     $                  LDA, ONE, A( J, J ), LDA )
            CALL SUTUJ( JB, A( J, J ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + J - 1
      RETURN
*
*     End of SUTUBJ
*
      END
      SUBROUTINE SUTUJ( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUTUJ  computes a U**T*U factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked left-looking (j-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ZERO
      PARAMETER          ( ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           STRSV
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 J = 1, N
         CALL STRSV( 'Upper', 'Transpose', 'Non-unit', J-1, A, LDA,
     $               A( 1, J ), 1 )
         A( J, J ) = A( J, J ) - SDOT( J-1, A( 1, J ), 1, A( 1, J ), 1 )
         IF( A( J, J ).LE.ZERO )
     $      GO TO 20
         A( J, J ) = SQRT( A( J, J ) )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = J
      RETURN
*
*     End of SUTUJ
*
      END
      SUBROUTINE SUTUBK( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUTUBK computes a U**T*U factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked right-looking (k-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K, KB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSYRK, STRSM, SUTUK
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SUTUK( N, A, LDA, INFO )
      ELSE
         DO 10 K = 1, N, NB
            KB = MIN( NB, N-K+1 )
            CALL SUTUK( KB, A( K, K ), LDA, INFO )
            IF( K+KB.LE.N ) THEN
               CALL STRSM( 'Left', 'Upper', 'Transpose', 'Non-unit', KB,
     $                     N-K-KB+1, ONE, A( K, K ), LDA, A( K, K+KB ),
     $                     LDA )
               CALL SSYRK( 'Upper', 'Transpose', N-K-KB+1, KB, -ONE,
     $                     A( K, K+KB ), LDA, ONE, A( K+KB, K+KB ),
     $                     LDA )
            END IF
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + K - 1
      RETURN
*
*     End of SUTUBK
*
      END
      SUBROUTINE SUTUK( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUTUK  computes a U**T*U factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked right-looking (k-variant) of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSCAL, SSYR
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 K = 1, N
         IF( A( K, K ).LE.ZERO )
     $      GO TO 20
         A( K, K ) = SQRT( A( K, K ) )
         IF( K.LT.N ) THEN
            CALL SSCAL( N-K, ONE/A( K, K ), A( K, K+1 ), LDA )
            CALL SSYR( 'Upper', N-K, -ONE, A( K, K+1 ), LDA,
     $                 A( K+1, K+1 ), LDA )
         END IF
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = K
      RETURN
*
*     End of SUTUK
*
      END
      SUBROUTINE SUUTBI( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUUTBI computes a U*U**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked i-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, IB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSYRK, STRSM, SUUTI
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SUUTI( N, A, LDA, INFO )
      ELSE
         DO 10 I = ( ( N-1 )/NB )*NB + 1, 1, -NB
            IB = MIN( NB, N-I+1 )
            IF( I+IB.LE.N ) THEN
               CALL STRSM( 'Right', 'Upper', 'Transpose', 'Non-unit',
     $                     IB, N-I-IB+1, ONE, A( I+IB, I+IB ), LDA,
     $                     A( I, I+IB ), LDA )
               CALL SSYRK( 'Upper', 'No transpose', IB, N-I-IB+1, -ONE,
     $                     A( I, I+IB ), LDA, ONE, A( I, I ), LDA )
            END IF
            CALL SUUTI( IB, A( I, I ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + I - 1
      RETURN
*
*     End of SUUTBI
*
      END
      SUBROUTINE SUUTI( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUUTI  computes a U*U**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked i-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ZERO
      PARAMETER          ( ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           STRSV
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 I = N, 1, -1
         IF( I.LT.N ) THEN
            CALL STRSV( 'Upper', 'No transpose', 'Non-unit', N-I,
     $                  A( I+1, I+1 ), LDA, A( I, I+1 ), LDA )
            A( I, I ) = A( I, I ) - SDOT( N-I, A( I, I+1 ), LDA,
     $                  A( I, I+1 ), LDA )
         END IF
         IF( A( I, I ).LE.ZERO )
     $      GO TO 20
         A( I, I ) = SQRT( A( I, I ) )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = I
      RETURN
*
*     End of SUUTI
*
      END
      SUBROUTINE SUUTBJ( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUUTBJ computes a U*U**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked j-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J, JB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMM, SSYRK, STRSM, SUUTJ
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SUUTJ( N, A, LDA, INFO )
      ELSE
         DO 10 J = ( ( N-1 )/NB )*NB + 1, 1, -NB
            JB = MIN( NB, N-J+1 )
            IF( J+JB.LE.N )
     $         CALL SSYRK( 'Upper', 'No transpose', JB, N-J-JB+1, -ONE,
     $                     A( J, J+JB ), LDA, ONE, A( J, J ), LDA )
            CALL SUUTJ( JB, A( J, J ), LDA, INFO )
            IF( INFO.NE.0 )
     $         GO TO 20
            IF( J+JB.LE.N )
     $         CALL SGEMM( 'No transpose', 'Transpose', J-1, JB,
     $                     N-J-JB+1, -ONE, A( 1, J+JB ), LDA,
     $                     A( J, J+JB ), LDA, ONE, A( 1, J ), LDA )
            CALL STRSM( 'Right', 'Upper', 'Transpose', 'Non-unit', J-1,
     $                  JB, ONE, A( J, J ), LDA, A( 1, J ), LDA )
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + J - 1
      RETURN
*
*     End of SUUTBJ
*
      END
      SUBROUTINE SUUTJ( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUUTJ  computes a U*U**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked j-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            J
*     ..
*     .. External Functions ..
      REAL               SDOT
      EXTERNAL           SDOT
*     ..
*     .. External Subroutines ..
      EXTERNAL           SGEMV, SSCAL
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 J = N, 1, -1
         IF( J.LT.N )
     $      A( J, J ) = A( J, J ) - SDOT( N-J, A( J, J+1 ), LDA,
     $                  A( J, J+1 ), LDA )
         IF( A( J, J ).LE.ZERO )
     $      GO TO 20
         A( J, J ) = SQRT( A( J, J ) )
         IF( J.LT.N )
     $      CALL SGEMV( 'No transpose', J-1, N-J, -ONE, A( 1, J+1 ),
     $                  LDA, A( J, J+1 ), LDA, ONE, A( 1, J ), 1 )
         CALL SSCAL( J-1, ONE/A( J, J ), A( 1, J ), 1 )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = J
      RETURN
*
*     End of SUUTJ
*
      END
      SUBROUTINE SUUTBK( N, A, LDA, INFO, NB )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N, NB
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUUTBK computes a U*U**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the blocked k-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  NB      (input) INTEGER
*          The blocksize for the blocked algorithm.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE
      PARAMETER          ( ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K, KB
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSYRK, STRSM, SUUTK
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          MIN
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      IF( NB.EQ.1 ) THEN
         CALL SUUTK( N, A, LDA, INFO )
      ELSE
         DO 10 K = ( ( N-1 )/NB )*NB + 1, 1, -NB
            KB = MIN( NB, N-K+1 )
            CALL SUUTK( KB, A( K, K ), LDA, INFO )
            CALL STRSM( 'Right', 'Upper', 'Transpose', 'Non-unit', K-1,
     $                  KB, ONE, A( K, K ), LDA, A( 1, K ), LDA )
            CALL SSYRK( 'Upper', 'No transpose', K-1, KB, -ONE,
     $                  A( 1, K ), LDA, ONE, A, LDA )
            IF( INFO.NE.0 )
     $         GO TO 20
   10    CONTINUE
      END IF
      RETURN
   20 CONTINUE
      INFO = INFO + K - 1
      RETURN
*
*     End of SUUTBK
*
      END
      SUBROUTINE SUUTK( N, A, LDA, INFO )
*
*  -- LAPACK variant --
*     Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab,
*     Courant Institute, NAG Ltd., and Rice University
*     January 17, 1990
*
*     .. Scalar Arguments ..
      INTEGER            INFO, LDA, N
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * )
*     ..
*
*  Purpose
*  =======
*
*  SUUTK  computes a U*U**T factorization of a real symmetric
*  positive definite matrix A.
*
*  This is the unblocked k-variant of the algorithm.
*
*  Arguments
*  =========
*
*  N       (input) INTEGER
*          The number of rows and columns of the matrix A.  N >= 0.
*
*  A       (input/output) REAL array, dimension (LDA,N)
*          On entry, the symmetric matrix A.  Only the upper triangular
*          part of A is referenced.
*          On exit, the factor U from the Cholesky factorization.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,N).
*
*  INFO    (output) INTEGER
*          = 0: successful exit
*          < 0: if INFO = -k, the k-th argument had an illegal value
*          > 0: if INFO = k, the leading minor of order K is not
*               positive definite, and the factorization could not be
*               completed.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            K
*     ..
*     .. External Subroutines ..
      EXTERNAL           SSCAL, SSYR
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          SQRT
*     ..
*     .. Executable Statements ..
*
      INFO = 0
      DO 10 K = N, 1, -1
         IF( A( K, K ).LE.ZERO )
     $      GO TO 20
         A( K, K ) = SQRT( A( K, K ) )
         CALL SSCAL( K-1, ONE/A( K, K ), A( 1, K ), 1 )
         CALL SSYR( 'Upper', K-1, -ONE, A( 1, K ), 1, A, LDA )
   10 CONTINUE
      RETURN
   20 CONTINUE
      INFO = K
      RETURN
*
*     End of SUUTK
*
      END