ScaLAPACK 2.1  2.1 ScaLAPACK: Scalable Linear Algebra PACKage
ddbtf2.f
Go to the documentation of this file.
1  SUBROUTINE ddbtf2( M, N, KL, KU, AB, LDAB, INFO )
2 *
3 * -- ScaLAPACK auxiliary routine (version 2.0) --
4 * Univ. of Tennessee, Univ. of California Berkeley, Univ. of Colorado Denver
5 *
6 * Modified by Andrew J. Cleary in November, 96 from:
7 * -- LAPACK auxiliary routine (preliminary version) --
8 * Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
9 * Courant Institute, Argonne National Lab, and Rice University
10 * August 6, 1991
11 *
12 *
13 * .. Scalar Arguments ..
14  INTEGER INFO, KL, KU, LDAB, M, N
15 * ..
16 * .. Array Arguments ..
17  DOUBLE PRECISION AB( LDAB, * )
18 * ..
19 *
20 * Purpose
21 * =======
22 *
23 * Ddbtrf computes an LU factorization of a real m-by-n band matrix A
24 * without using partial pivoting with row interchanges.
25 *
26 * This is the unblocked version of the algorithm, calling Level 2 BLAS.
27 *
28 * Arguments
29 * =========
30 *
31 * M (input) INTEGER
32 * The number of rows of the matrix A. M >= 0.
33 *
34 * N (input) INTEGER
35 * The number of columns 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 * AB (input/output) DOUBLE PRECISION array, dimension (LDAB,N)
44 * On entry, the matrix A in band storage, in rows KL+1 to
45 * 2*KL+KU+1; rows 1 to KL of the array need not be set.
46 * The j-th column of A is stored in the j-th column of the
47 * array AB as follows:
48 * AB(kl+ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(m,j+kl)
49 *
50 * On exit, details of the factorization: U is stored as an
51 * upper triangular band matrix with KL+KU superdiagonals in
52 * rows 1 to KL+KU+1, and the multipliers used during the
53 * factorization are stored in rows KL+KU+2 to 2*KL+KU+1.
54 * See below for further details.
55 *
56 * LDAB (input) INTEGER
57 * The leading dimension of the array AB. LDAB >= 2*KL+KU+1.
58 *
59 * INFO (output) INTEGER
60 * = 0: successful exit
61 * < 0: if INFO = -i, the i-th argument had an illegal value
62 * > 0: if INFO = +i, U(i,i) is exactly zero. The factorization
63 * has been completed, but the factor U is exactly
64 * singular, and division by zero will occur if it is used
65 * to solve a system of equations.
66 *
67 * Further Details
68 * ===============
69 *
70 * The band storage scheme is illustrated by the following example, when
71 * M = N = 6, KL = 2, KU = 1:
72 *
73 * On entry: On exit:
74 *
75 * * a12 a23 a34 a45 a56 * u12 u23 u34 u45 u56
76 * a11 a22 a33 a44 a55 a66 u11 u22 u33 u44 u55 u66
77 * a21 a32 a43 a54 a65 * m21 m32 m43 m54 m65 *
78 * a31 a42 a53 a64 * * m31 m42 m53 m64 * *
79 *
80 * Array elements marked * are not used by the routine; elements marked
81 * + need not be set on entry, but are required by the routine to store
82 * elements of U, because of fill-in resulting from the row
83 * interchanges.
84 *
85 * =====================================================================
86 *
87 * .. Parameters ..
88  DOUBLE PRECISION ONE, ZERO
89  parameter( one = 1.0d+0 )
90  parameter( zero = 0.0d+0 )
91 * ..
92 * .. Local Scalars ..
93  INTEGER J, JP, JU, KM, KV
94 * ..
95 * .. External Functions ..
96  INTEGER ISAMAX
97  EXTERNAL isamax
98 * ..
99 * .. External Subroutines ..
100  EXTERNAL dger, dscal, dswap
101 * ..
102 * .. Intrinsic Functions ..
103  INTRINSIC max, min
104 * ..
105 * .. Executable Statements ..
106 *
107 * KV is the number of superdiagonals in the factor U, allowing for
108 * fill-in.
109 *
110  kv = ku
111 *
112 * Test the input parameters.
113 *
114  info = 0
115 *ECA IF( M.LT.0 ) THEN
116 *ECA INFO = -1
117 *ECA ELSE IF( N.LT.0 ) THEN
118 *ECA INFO = -2
119 *ECA ELSE IF( KL.LT.0 ) THEN
120 *ECA INFO = -3
121 *ECA ELSE IF( KU.LT.0 ) THEN
122 *ECA INFO = -4
123 *ECA ELSE IF( LDAB.LT.KL+KV+1 ) THEN
124 *ECA INFO = -6
125 *ECA END IF
126 *ECA IF( INFO.NE.0 ) THEN
127 *ECA CALL XERBLA( 'DDBTF2', -INFO )
128 *ECA RETURN
129 *ECA END IF
130 *
131 * Quick return if possible
132 *
133  IF( m.EQ.0 .OR. n.EQ.0 )
134  \$ RETURN
135 *
136 * Gaussian elimination without partial pivoting
137 *
138 * JU is the index of the last column affected by the current stage
139 * of the factorization.
140 *
141  ju = 1
142 *
143  DO 40 j = 1, min( m, n )
144 *
145 * Test for singularity. KM is the number of
146 * subdiagonal elements in the current column.
147 *
148  km = min( kl, m-j )
149  jp = 1
150  IF( ab( kv+1, j ).NE.zero ) THEN
151  ju = max( ju, min( j+ku, n ) )
152 *
153  IF( km.GT.0 ) THEN
154 *
155 * Compute multipliers.
156 *
157  CALL dscal( km, one / ab( ku+1, j ), ab( ku+2, j ), 1 )
158 *
159 * Update trailing submatrix within the band.
160 *
161  IF( ju.GT.j ) THEN
162  CALL dger( km, ju-j, -one, ab( ku+2, j ), 1,
163  \$ ab( ku, j+1 ), ldab-1, ab( ku+1, j+1 ),
164  \$ ldab-1 )
165  END IF
166  END IF
167  ELSE
168 *
169  IF( info.EQ.0 )
170  \$ info = j
171  END IF
172  40 CONTINUE
173  RETURN
174 *
175 * End of DDBTF2
176 *
177  END
max
#define max(A, B)
Definition: pcgemr.c:180
ddbtf2
subroutine ddbtf2(M, N, KL, KU, AB, LDAB, INFO)
Definition: ddbtf2.f:2
min
#define min(A, B)
Definition: pcgemr.c:181