Mercurial > octave-nkf

comparison libcruft/lapack/sgbtrs.f @ 7789:82be108cc558

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
First attempt at single precision tyeps * * * corrections to qrupdate single precision routines * * * prefer demotion to single over promotion to double * * * Add single precision support to log2 function * * * Trivial PROJECT file update * * * Cache optimized hermitian/transpose methods * * * Add tests for tranpose/hermitian and ChangeLog entry for new transpose code
author David Bateman <dbateman@free.fr>
date Sun, 27 Apr 2008 22:34:17 +0200
parents
children
comparison
equal deleted inserted replaced
7788:45f5faba05a2 7789:82be108cc558
1 SUBROUTINE SGBTRS( TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB,
2 $ INFO )
3 *
4 * -- LAPACK routine (version 3.1) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2006
7 *
8 * .. Scalar Arguments ..
9 CHARACTER TRANS
10 INTEGER INFO, KL, KU, LDAB, LDB, N, NRHS
11 * ..
12 * .. Array Arguments ..
13 INTEGER IPIV( * )
14 REAL AB( LDAB, * ), B( LDB, * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * SGBTRS solves a system of linear equations
21 * A * X = B or A' * X = B
22 * with a general band matrix A using the LU factorization computed
23 * by SGBTRF.
24 *
25 * Arguments
26 * =========
27 *
28 * TRANS (input) CHARACTER*1
29 * Specifies the form of the system of equations.
30 * = 'N': A * X = B (No transpose)
31 * = 'T': A'* X = B (Transpose)
32 * = 'C': A'* X = B (Conjugate transpose = Transpose)
33 *
34 * N (input) INTEGER
35 * The order of the matrix A. N >= 0.
36 *
37 * KL (input) INTEGER
38 * The number of subdiagonals within the band of A. KL >= 0.
39 *
40 * KU (input) INTEGER
41 * The number of superdiagonals within the band of A. KU >= 0.
42 *
43 * NRHS (input) INTEGER
44 * The number of right hand sides, i.e., the number of columns
45 * of the matrix B. NRHS >= 0.
46 *
47 * AB (input) REAL array, dimension (LDAB,N)
48 * Details of the LU factorization of the band matrix A, as
49 * computed by SGBTRF. U is stored as an upper triangular band
50 * matrix with KL+KU superdiagonals in rows 1 to KL+KU+1, and
51 * the multipliers used during the factorization are stored in
52 * rows KL+KU+2 to 2*KL+KU+1.
53 *
54 * LDAB (input) INTEGER
55 * The leading dimension of the array AB. LDAB >= 2*KL+KU+1.
56 *
57 * IPIV (input) INTEGER array, dimension (N)
58 * The pivot indices; for 1 <= i <= N, row i of the matrix was
59 * interchanged with row IPIV(i).
60 *
61 * B (input/output) REAL array, dimension (LDB,NRHS)
62 * On entry, the right hand side matrix B.
63 * On exit, the solution matrix X.
64 *
65 * LDB (input) INTEGER
66 * The leading dimension of the array B. LDB >= max(1,N).
67 *
68 * INFO (output) INTEGER
69 * = 0: successful exit
70 * < 0: if INFO = -i, the i-th argument had an illegal value
71 *
72 * =====================================================================
73 *
74 * .. Parameters ..
75 REAL ONE
76 PARAMETER ( ONE = 1.0E+0 )
77 * ..
78 * .. Local Scalars ..
79 LOGICAL LNOTI, NOTRAN
80 INTEGER I, J, KD, L, LM
81 * ..
82 * .. External Functions ..
83 LOGICAL LSAME
84 EXTERNAL LSAME
85 * ..
86 * .. External Subroutines ..
87 EXTERNAL SGEMV, SGER, SSWAP, STBSV, XERBLA
88 * ..
89 * .. Intrinsic Functions ..
90 INTRINSIC MAX, MIN
91 * ..
92 * .. Executable Statements ..
93 *
94 * Test the input parameters.
95 *
96 INFO = 0
97 NOTRAN = LSAME( TRANS, 'N' )
98 IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'T' ) .AND. .NOT.
99 $ LSAME( TRANS, 'C' ) ) THEN
100 INFO = -1
101 ELSE IF( N.LT.0 ) THEN
102 INFO = -2
103 ELSE IF( KL.LT.0 ) THEN
104 INFO = -3
105 ELSE IF( KU.LT.0 ) THEN
106 INFO = -4
107 ELSE IF( NRHS.LT.0 ) THEN
108 INFO = -5
109 ELSE IF( LDAB.LT.( 2*KL+KU+1 ) ) THEN
110 INFO = -7
111 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
112 INFO = -10
113 END IF
114 IF( INFO.NE.0 ) THEN
115 CALL XERBLA( 'SGBTRS', -INFO )
116 RETURN
117 END IF
118 *
119 * Quick return if possible
120 *
121 IF( N.EQ.0 .OR. NRHS.EQ.0 )
122 $ RETURN
123 *
124 KD = KU + KL + 1
125 LNOTI = KL.GT.0
126 *
127 IF( NOTRAN ) THEN
128 *
129 * Solve A*X = B.
130 *
131 * Solve L*X = B, overwriting B with X.
132 *
133 * L is represented as a product of permutations and unit lower
134 * triangular matrices L = P(1) * L(1) * ... * P(n-1) * L(n-1),
135 * where each transformation L(i) is a rank-one modification of
136 * the identity matrix.
137 *
138 IF( LNOTI ) THEN
139 DO 10 J = 1, N - 1
140 LM = MIN( KL, N-J )
141 L = IPIV( J )
142 IF( L.NE.J )
143 $ CALL SSWAP( NRHS, B( L, 1 ), LDB, B( J, 1 ), LDB )
144 CALL SGER( LM, NRHS, -ONE, AB( KD+1, J ), 1, B( J, 1 ),
145 $ LDB, B( J+1, 1 ), LDB )
146 10 CONTINUE
147 END IF
148 *
149 DO 20 I = 1, NRHS
150 *
151 * Solve U*X = B, overwriting B with X.
152 *
153 CALL STBSV( 'Upper', 'No transpose', 'Non-unit', N, KL+KU,
154 $ AB, LDAB, B( 1, I ), 1 )
155 20 CONTINUE
156 *
157 ELSE
158 *
159 * Solve A'*X = B.
160 *
161 DO 30 I = 1, NRHS
162 *
163 * Solve U'*X = B, overwriting B with X.
164 *
165 CALL STBSV( 'Upper', 'Transpose', 'Non-unit', N, KL+KU, AB,
166 $ LDAB, B( 1, I ), 1 )
167 30 CONTINUE
168 *
169 * Solve L'*X = B, overwriting B with X.
170 *
171 IF( LNOTI ) THEN
172 DO 40 J = N - 1, 1, -1
173 LM = MIN( KL, N-J )
174 CALL SGEMV( 'Transpose', LM, NRHS, -ONE, B( J+1, 1 ),
175 $ LDB, AB( KD+1, J ), 1, ONE, B( J, 1 ), LDB )
176 L = IPIV( J )
177 IF( L.NE.J )
178 $ CALL SSWAP( NRHS, B( L, 1 ), LDB, B( J, 1 ), LDB )
179 40 CONTINUE
180 END IF
181 END IF
182 RETURN
183 *
184 * End of SGBTRS
185 *
186 END