LAPACK
3.5.0
LAPACK: Linear Algebra PACKage
Main Page
Data Types List
Files
File List
File Members
All
Classes
Files
Functions
Variables
Typedefs
Macros
zpbtf2.f
Go to the documentation of this file.
1
*> \brief \b ZPBTF2 computes the Cholesky factorization of a symmetric/Hermitian positive definite band matrix (unblocked algorithm).
2
*
3
* =========== DOCUMENTATION ===========
4
*
5
* Online html documentation available at
6
* http://www.netlib.org/lapack/explore-html/
7
*
8
*> \htmlonly
9
*> Download ZPBTF2 + dependencies
10
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zpbtf2.f">
11
*> [TGZ]</a>
12
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zpbtf2.f">
13
*> [ZIP]</a>
14
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zpbtf2.f">
15
*> [TXT]</a>
16
*> \endhtmlonly
17
*
18
* Definition:
19
* ===========
20
*
21
* SUBROUTINE ZPBTF2( UPLO, N, KD, AB, LDAB, INFO )
22
*
23
* .. Scalar Arguments ..
24
* CHARACTER UPLO
25
* INTEGER INFO, KD, LDAB, N
26
* ..
27
* .. Array Arguments ..
28
* COMPLEX*16 AB( LDAB, * )
29
* ..
30
*
31
*
32
*> \par Purpose:
33
* =============
34
*>
35
*> \verbatim
36
*>
37
*> ZPBTF2 computes the Cholesky factorization of a complex Hermitian
38
*> positive definite band matrix A.
39
*>
40
*> The factorization has the form
41
*> A = U**H * U , if UPLO = 'U', or
42
*> A = L * L**H, if UPLO = 'L',
43
*> where U is an upper triangular matrix, U**H is the conjugate transpose
44
*> of U, and L is lower triangular.
45
*>
46
*> This is the unblocked version of the algorithm, calling Level 2 BLAS.
47
*> \endverbatim
48
*
49
* Arguments:
50
* ==========
51
*
52
*> \param[in] UPLO
53
*> \verbatim
54
*> UPLO is CHARACTER*1
55
*> Specifies whether the upper or lower triangular part of the
56
*> Hermitian matrix A is stored:
57
*> = 'U': Upper triangular
58
*> = 'L': Lower triangular
59
*> \endverbatim
60
*>
61
*> \param[in] N
62
*> \verbatim
63
*> N is INTEGER
64
*> The order of the matrix A. N >= 0.
65
*> \endverbatim
66
*>
67
*> \param[in] KD
68
*> \verbatim
69
*> KD is INTEGER
70
*> The number of super-diagonals of the matrix A if UPLO = 'U',
71
*> or the number of sub-diagonals if UPLO = 'L'. KD >= 0.
72
*> \endverbatim
73
*>
74
*> \param[in,out] AB
75
*> \verbatim
76
*> AB is COMPLEX*16 array, dimension (LDAB,N)
77
*> On entry, the upper or lower triangle of the Hermitian band
78
*> matrix A, stored in the first KD+1 rows of the array. The
79
*> j-th column of A is stored in the j-th column of the array AB
80
*> as follows:
81
*> if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
82
*> if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd).
83
*>
84
*> On exit, if INFO = 0, the triangular factor U or L from the
85
*> Cholesky factorization A = U**H *U or A = L*L**H of the band
86
*> matrix A, in the same storage format as A.
87
*> \endverbatim
88
*>
89
*> \param[in] LDAB
90
*> \verbatim
91
*> LDAB is INTEGER
92
*> The leading dimension of the array AB. LDAB >= KD+1.
93
*> \endverbatim
94
*>
95
*> \param[out] INFO
96
*> \verbatim
97
*> INFO is INTEGER
98
*> = 0: successful exit
99
*> < 0: if INFO = -k, the k-th argument had an illegal value
100
*> > 0: if INFO = k, the leading minor of order k is not
101
*> positive definite, and the factorization could not be
102
*> completed.
103
*> \endverbatim
104
*
105
* Authors:
106
* ========
107
*
108
*> \author Univ. of Tennessee
109
*> \author Univ. of California Berkeley
110
*> \author Univ. of Colorado Denver
111
*> \author NAG Ltd.
112
*
113
*> \date September 2012
114
*
115
*> \ingroup complex16OTHERcomputational
116
*
117
*> \par Further Details:
118
* =====================
119
*>
120
*> \verbatim
121
*>
122
*> The band storage scheme is illustrated by the following example, when
123
*> N = 6, KD = 2, and UPLO = 'U':
124
*>
125
*> On entry: On exit:
126
*>
127
*> * * a13 a24 a35 a46 * * u13 u24 u35 u46
128
*> * a12 a23 a34 a45 a56 * u12 u23 u34 u45 u56
129
*> a11 a22 a33 a44 a55 a66 u11 u22 u33 u44 u55 u66
130
*>
131
*> Similarly, if UPLO = 'L' the format of A is as follows:
132
*>
133
*> On entry: On exit:
134
*>
135
*> a11 a22 a33 a44 a55 a66 l11 l22 l33 l44 l55 l66
136
*> a21 a32 a43 a54 a65 * l21 l32 l43 l54 l65 *
137
*> a31 a42 a53 a64 * * l31 l42 l53 l64 * *
138
*>
139
*> Array elements marked * are not used by the routine.
140
*> \endverbatim
141
*>
142
* =====================================================================
143
SUBROUTINE
zpbtf2
( UPLO, N, KD, AB, LDAB, INFO )
144
*
145
* -- LAPACK computational routine (version 3.4.2) --
146
* -- LAPACK is a software package provided by Univ. of Tennessee, --
147
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
148
* September 2012
149
*
150
* .. Scalar Arguments ..
151
CHARACTER
uplo
152
INTEGER
info, kd, ldab, n
153
* ..
154
* .. Array Arguments ..
155
COMPLEX*16
ab( ldab, * )
156
* ..
157
*
158
* =====================================================================
159
*
160
* .. Parameters ..
161
DOUBLE PRECISION
one, zero
162
parameter( one = 1.0d+0, zero = 0.0d+0 )
163
* ..
164
* .. Local Scalars ..
165
LOGICAL
upper
166
INTEGER
j
, kld, kn
167
DOUBLE PRECISION
ajj
168
* ..
169
* .. External Functions ..
170
LOGICAL
lsame
171
EXTERNAL
lsame
172
* ..
173
* .. External Subroutines ..
174
EXTERNAL
xerbla
,
zdscal
,
zher
,
zlacgv
175
* ..
176
* .. Intrinsic Functions ..
177
INTRINSIC
dble, max, min, sqrt
178
* ..
179
* .. Executable Statements ..
180
*
181
* Test the input parameters.
182
*
183
info = 0
184
upper =
lsame
( uplo,
'U'
)
185
IF
( .NOT.upper .AND. .NOT.
lsame
( uplo,
'L'
) )
THEN
186
info = -1
187
ELSE
IF
( n.LT.0 )
THEN
188
info = -2
189
ELSE
IF
( kd.LT.0 )
THEN
190
info = -3
191
ELSE
IF
( ldab.LT.kd+1 )
THEN
192
info = -5
193
END IF
194
IF
( info.NE.0 )
THEN
195
CALL
xerbla
(
'ZPBTF2'
, -info )
196
RETURN
197
END IF
198
*
199
* Quick return if possible
200
*
201
IF
( n.EQ.0 )
202
$
RETURN
203
*
204
kld = max( 1, ldab-1 )
205
*
206
IF
( upper )
THEN
207
*
208
* Compute the Cholesky factorization A = U**H * U.
209
*
210
DO
10
j
= 1, n
211
*
212
* Compute U(J,J) and test for non-positive-definiteness.
213
*
214
ajj = dble( ab( kd+1,
j
) )
215
IF
( ajj.LE.zero )
THEN
216
ab( kd+1,
j
) = ajj
217
go to 30
218
END IF
219
ajj = sqrt( ajj )
220
ab( kd+1,
j
) = ajj
221
*
222
* Compute elements J+1:J+KN of row J and update the
223
* trailing submatrix within the band.
224
*
225
kn = min( kd, n-
j
)
226
IF
( kn.GT.0 )
THEN
227
CALL
zdscal
( kn, one / ajj, ab( kd,
j
+1 ), kld )
228
CALL
zlacgv
( kn, ab( kd,
j
+1 ), kld )
229
CALL
zher
(
'Upper'
, kn, -one, ab( kd,
j
+1 ), kld,
230
$ ab( kd+1,
j
+1 ), kld )
231
CALL
zlacgv
( kn, ab( kd,
j
+1 ), kld )
232
END IF
233
10
CONTINUE
234
ELSE
235
*
236
* Compute the Cholesky factorization A = L*L**H.
237
*
238
DO
20
j
= 1, n
239
*
240
* Compute L(J,J) and test for non-positive-definiteness.
241
*
242
ajj = dble( ab( 1,
j
) )
243
IF
( ajj.LE.zero )
THEN
244
ab( 1,
j
) = ajj
245
go to 30
246
END IF
247
ajj = sqrt( ajj )
248
ab( 1,
j
) = ajj
249
*
250
* Compute elements J+1:J+KN of column J and update the
251
* trailing submatrix within the band.
252
*
253
kn = min( kd, n-
j
)
254
IF
( kn.GT.0 )
THEN
255
CALL
zdscal
( kn, one / ajj, ab( 2,
j
), 1 )
256
CALL
zher
(
'Lower'
, kn, -one, ab( 2,
j
), 1,
257
$ ab( 1,
j
+1 ), kld )
258
END IF
259
20
CONTINUE
260
END IF
261
RETURN
262
*
263
30
CONTINUE
264
info =
j
265
RETURN
266
*
267
* End of ZPBTF2
268
*
269
END
src
zpbtf2.f
Generated on Mon Dec 30 2013 16:10:53 for LAPACK by
1.8.1.2