1        2        3        4        5        6        7        8        9       10       11       12       13       14       15       16       17       18       19       20       21       22       23       24       25       26       27       28       29       30       31       32       33       34       35       36       37       38       39       40       41       42       43       44       45       46       47       48       49       50       51       52       53       54       55       56       57       58       59       60       61       62       63       64       65       66       67       68       69       70       71       72       73       74       75       76       77       78       79       80       81       82       83       84       85       86       87       88       89       90       91       92       93       94       95       96       97       98       99      100      101      102      103      104      105      106      107      108      109      110      111      112      113      114      115      116      117      118      119      120      121      122      123      124      125      126      127      128      129      130      131      132      133      134      135      136      137      138      139      140      141      142      143      144      145      146      147      148      149      150      151      152 SUBROUTINE DPBSV( UPLO, N, KD, NRHS, AB, LDAB, B, LDB, INFO ) * *  -- LAPACK driver routine (version 3.3.1) -- *  -- LAPACK is a software package provided by Univ. of Tennessee,    -- *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- *  -- April 2011                                                      -- * *     .. Scalar Arguments ..       CHARACTER          UPLO       INTEGER            INFO, KD, LDAB, LDB, N, NRHS *     .. *     .. Array Arguments ..       DOUBLE PRECISION   AB( LDAB, * ), B( LDB, * ) *     .. * *  Purpose *  ======= * *  DPBSV computes the solution to a real system of linear equations *     A * X = B, *  where A is an N-by-N symmetric positive definite band matrix and X *  and B are N-by-NRHS matrices. * *  The Cholesky decomposition is used to factor A as *     A = U**T * U,  if UPLO = 'U', or *     A = L * L**T,  if UPLO = 'L', *  where U is an upper triangular band matrix, and L is a lower *  triangular band matrix, with the same number of superdiagonals or *  subdiagonals as A.  The factored form of A is then used to solve the *  system of equations A * X = B. * *  Arguments *  ========= * *  UPLO    (input) CHARACTER*1 *          = 'U':  Upper triangle of A is stored; *          = 'L':  Lower triangle of A is stored. * *  N       (input) INTEGER *          The number of linear equations, i.e., the order of the *          matrix A.  N >= 0. * *  KD      (input) INTEGER *          The number of superdiagonals of the matrix A if UPLO = 'U', *          or the number of subdiagonals if UPLO = 'L'.  KD >= 0. * *  NRHS    (input) INTEGER *          The number of right hand sides, i.e., the number of columns *          of the matrix B.  NRHS >= 0. * *  AB      (input/output) DOUBLE PRECISION array, dimension (LDAB,N) *          On entry, the upper or lower triangle of the symmetric band *          matrix A, stored in the first KD+1 rows of the array.  The *          j-th column of A is stored in the j-th column of the array AB *          as follows: *          if UPLO = 'U', AB(KD+1+i-j,j) = A(i,j) for max(1,j-KD)<=i<=j; *          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(N,j+KD). *          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 of the band *          matrix A, in the same storage format as A. * *  LDAB    (input) INTEGER *          The leading dimension of the array AB.  LDAB >= KD+1. * *  B       (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS) *          On entry, the N-by-NRHS right hand side matrix B. *          On exit, if INFO = 0, the N-by-NRHS solution matrix X. * *  LDB     (input) INTEGER *          The leading dimension of the array B.  LDB >= max(1,N). * *  INFO    (output) INTEGER *          = 0:  successful exit *          < 0:  if INFO = -i, the i-th argument had an illegal value *          > 0:  if INFO = i, the leading minor of order i of A is not *                positive definite, so the factorization could not be *                completed, and the solution has not been computed. * *  Further Details *  =============== * *  The band storage scheme is illustrated by the following example, when *  N = 6, KD = 2, and UPLO = 'U': * *  On entry:                       On exit: * *      *    *   a13  a24  a35  a46      *    *   u13  u24  u35  u46 *      *   a12  a23  a34  a45  a56      *   u12  u23  u34  u45  u56 *     a11  a22  a33  a44  a55  a66     u11  u22  u33  u44  u55  u66 * *  Similarly, if UPLO = 'L' the format of A is as follows: * *  On entry:                       On exit: * *     a11  a22  a33  a44  a55  a66     l11  l22  l33  l44  l55  l66 *     a21  a32  a43  a54  a65   *      l21  l32  l43  l54  l65   * *     a31  a42  a53  a64   *    *      l31  l42  l53  l64   *    * * *  Array elements marked * are not used by the routine. * *  ===================================================================== * *     .. External Functions ..       LOGICAL            LSAME       EXTERNAL           LSAME *     .. *     .. External Subroutines ..       EXTERNAL           DPBTRF, DPBTRS, XERBLA *     .. *     .. Intrinsic Functions ..       INTRINSIC          MAX *     .. *     .. Executable Statements .. * *     Test the input parameters. *       INFO = 0       IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN          INFO = -1       ELSE IF( N.LT.0 ) THEN          INFO = -2       ELSE IF( KD.LT.0 ) THEN          INFO = -3       ELSE IF( NRHS.LT.0 ) THEN          INFO = -4       ELSE IF( LDAB.LT.KD+1 ) THEN          INFO = -6       ELSE IF( LDB.LT.MAX( 1, N ) ) THEN          INFO = -8       END IF       IF( INFO.NE.0 ) THEN          CALL XERBLA( 'DPBSV ', -INFO )          RETURN       END IF * *     Compute the Cholesky factorization A = U**T*U or A = L*L**T. *       CALL DPBTRF( UPLO, N, KD, AB, LDAB, INFO )       IF( INFO.EQ.0 ) THEN * *        Solve the system A*X = B, overwriting B with X. *          CALL DPBTRS( UPLO, N, KD, NRHS, AB, LDAB, B, LDB, INFO ) *       END IF       RETURN * *     End of DPBSV *       END