1: SUBROUTINE DPBTRS( UPLO, N, KD, NRHS, AB, LDAB, B, LDB, INFO )
2: *
3: * -- LAPACK routine (version 3.2) --
4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6: * November 2006
7: *
8: * .. Scalar Arguments ..
9: CHARACTER UPLO
10: INTEGER INFO, KD, LDAB, LDB, N, NRHS
11: * ..
12: * .. Array Arguments ..
13: DOUBLE PRECISION AB( LDAB, * ), B( LDB, * )
14: * ..
15: *
16: * Purpose
17: * =======
18: *
19: * DPBTRS solves a system of linear equations A*X = B with a symmetric
20: * positive definite band matrix A using the Cholesky factorization
21: * A = U**T*U or A = L*L**T computed by DPBTRF.
22: *
23: * Arguments
24: * =========
25: *
26: * UPLO (input) CHARACTER*1
27: * = 'U': Upper triangular factor stored in AB;
28: * = 'L': Lower triangular factor stored in AB.
29: *
30: * N (input) INTEGER
31: * The order of the matrix A. N >= 0.
32: *
33: * KD (input) INTEGER
34: * The number of superdiagonals of the matrix A if UPLO = 'U',
35: * or the number of subdiagonals if UPLO = 'L'. KD >= 0.
36: *
37: * NRHS (input) INTEGER
38: * The number of right hand sides, i.e., the number of columns
39: * of the matrix B. NRHS >= 0.
40: *
41: * AB (input) DOUBLE PRECISION array, dimension (LDAB,N)
42: * The triangular factor U or L from the Cholesky factorization
43: * A = U**T*U or A = L*L**T of the band matrix A, stored in the
44: * first KD+1 rows of the array. The j-th column of U or L is
45: * stored in the j-th column of the array AB as follows:
46: * if UPLO ='U', AB(kd+1+i-j,j) = U(i,j) for max(1,j-kd)<=i<=j;
47: * if UPLO ='L', AB(1+i-j,j) = L(i,j) for j<=i<=min(n,j+kd).
48: *
49: * LDAB (input) INTEGER
50: * The leading dimension of the array AB. LDAB >= KD+1.
51: *
52: * B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
53: * On entry, the right hand side matrix B.
54: * On exit, the solution matrix X.
55: *
56: * LDB (input) INTEGER
57: * The leading dimension of the array B. LDB >= max(1,N).
58: *
59: * INFO (output) INTEGER
60: * = 0: successful exit
61: * < 0: if INFO = -i, the i-th argument had an illegal value
62: *
63: * =====================================================================
64: *
65: * .. Local Scalars ..
66: LOGICAL UPPER
67: INTEGER J
68: * ..
69: * .. External Functions ..
70: LOGICAL LSAME
71: EXTERNAL LSAME
72: * ..
73: * .. External Subroutines ..
74: EXTERNAL DTBSV, XERBLA
75: * ..
76: * .. Intrinsic Functions ..
77: INTRINSIC MAX
78: * ..
79: * .. Executable Statements ..
80: *
81: * Test the input parameters.
82: *
83: INFO = 0
84: UPPER = LSAME( UPLO, 'U' )
85: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
86: INFO = -1
87: ELSE IF( N.LT.0 ) THEN
88: INFO = -2
89: ELSE IF( KD.LT.0 ) THEN
90: INFO = -3
91: ELSE IF( NRHS.LT.0 ) THEN
92: INFO = -4
93: ELSE IF( LDAB.LT.KD+1 ) THEN
94: INFO = -6
95: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
96: INFO = -8
97: END IF
98: IF( INFO.NE.0 ) THEN
99: CALL XERBLA( 'DPBTRS', -INFO )
100: RETURN
101: END IF
102: *
103: * Quick return if possible
104: *
105: IF( N.EQ.0 .OR. NRHS.EQ.0 )
106: $ RETURN
107: *
108: IF( UPPER ) THEN
109: *
110: * Solve A*X = B where A = U'*U.
111: *
112: DO 10 J = 1, NRHS
113: *
114: * Solve U'*X = B, overwriting B with X.
115: *
116: CALL DTBSV( 'Upper', 'Transpose', 'Non-unit', N, KD, AB,
117: $ LDAB, B( 1, J ), 1 )
118: *
119: * Solve U*X = B, overwriting B with X.
120: *
121: CALL DTBSV( 'Upper', 'No transpose', 'Non-unit', N, KD, AB,
122: $ LDAB, B( 1, J ), 1 )
123: 10 CONTINUE
124: ELSE
125: *
126: * Solve A*X = B where A = L*L'.
127: *
128: DO 20 J = 1, NRHS
129: *
130: * Solve L*X = B, overwriting B with X.
131: *
132: CALL DTBSV( 'Lower', 'No transpose', 'Non-unit', N, KD, AB,
133: $ LDAB, B( 1, J ), 1 )
134: *
135: * Solve L'*X = B, overwriting B with X.
136: *
137: CALL DTBSV( 'Lower', 'Transpose', 'Non-unit', N, KD, AB,
138: $ LDAB, B( 1, J ), 1 )
139: 20 CONTINUE
140: END IF
141: *
142: RETURN
143: *
144: * End of DPBTRS
145: *
146: END
CVSweb interface <joel.bertrand@systella.fr>