128
129
130
131
132
133
134 CHARACTER UPLO
135 INTEGER INFO, ITYPE, LDA, LDB, N
136
137
138 COMPLEX A( LDA, * ), B( LDB, * )
139
140
141
142
143
144 REAL ONE
145 parameter( one = 1.0e+0 )
146 COMPLEX CONE, HALF
147 parameter( cone = ( 1.0e+0, 0.0e+0 ),
148 $ half = ( 0.5e+0, 0.0e+0 ) )
149
150
151 LOGICAL UPPER
152 INTEGER K, KB, NB
153
154
156
157
158 INTRINSIC max, min
159
160
161 LOGICAL LSAME
162 INTEGER ILAENV
164
165
166
167
168
169 info = 0
170 upper =
lsame( uplo,
'U' )
171 IF( itype.LT.1 .OR. itype.GT.3 ) THEN
172 info = -1
173 ELSE IF( .NOT.upper .AND. .NOT.
lsame( uplo,
'L' ) )
THEN
174 info = -2
175 ELSE IF( n.LT.0 ) THEN
176 info = -3
177 ELSE IF( lda.LT.max( 1, n ) ) THEN
178 info = -5
179 ELSE IF( ldb.LT.max( 1, n ) ) THEN
180 info = -7
181 END IF
182 IF( info.NE.0 ) THEN
183 CALL xerbla(
'CHEGST', -info )
184 RETURN
185 END IF
186
187
188
189 IF( n.EQ.0 )
190 $ RETURN
191
192
193
194 nb =
ilaenv( 1,
'CHEGST', uplo, n, -1, -1, -1 )
195
196 IF( nb.LE.1 .OR. nb.GE.n ) THEN
197
198
199
200 CALL chegs2( itype, uplo, n, a, lda, b, ldb, info )
201 ELSE
202
203
204
205 IF( itype.EQ.1 ) THEN
206 IF( upper ) THEN
207
208
209
210 DO 10 k = 1, n, nb
211 kb = min( n-k+1, nb )
212
213
214
215 CALL chegs2( itype, uplo, kb, a( k, k ), lda,
216 $ b( k, k ), ldb, info )
217 IF( k+kb.LE.n ) THEN
218 CALL ctrsm(
'Left', uplo,
'Conjugate transpose',
219 $ 'Non-unit', kb, n-k-kb+1, cone,
220 $ b( k, k ), ldb, a( k, k+kb ), lda )
221 CALL chemm(
'Left', uplo, kb, n-k-kb+1, -half,
222 $ a( k, k ), lda, b( k, k+kb ), ldb,
223 $ cone, a( k, k+kb ), lda )
224 CALL cher2k( uplo,
'Conjugate transpose', n-k-kb+1,
225 $ kb, -cone, a( k, k+kb ), lda,
226 $ b( k, k+kb ), ldb, one,
227 $ a( k+kb, k+kb ), lda )
228 CALL chemm(
'Left', uplo, kb, n-k-kb+1, -half,
229 $ a( k, k ), lda, b( k, k+kb ), ldb,
230 $ cone, a( k, k+kb ), lda )
231 CALL ctrsm(
'Right', uplo,
'No transpose',
232 $ 'Non-unit', kb, n-k-kb+1, cone,
233 $ b( k+kb, k+kb ), ldb, a( k, k+kb ),
234 $ lda )
235 END IF
236 10 CONTINUE
237 ELSE
238
239
240
241 DO 20 k = 1, n, nb
242 kb = min( n-k+1, nb )
243
244
245
246 CALL chegs2( itype, uplo, kb, a( k, k ), lda,
247 $ b( k, k ), ldb, info )
248 IF( k+kb.LE.n ) THEN
249 CALL ctrsm(
'Right', uplo,
'Conjugate transpose',
250 $ 'Non-unit', n-k-kb+1, kb, cone,
251 $ b( k, k ), ldb, a( k+kb, k ), lda )
252 CALL chemm(
'Right', uplo, n-k-kb+1, kb, -half,
253 $ a( k, k ), lda, b( k+kb, k ), ldb,
254 $ cone, a( k+kb, k ), lda )
255 CALL cher2k( uplo,
'No transpose', n-k-kb+1, kb,
256 $ -cone, a( k+kb, k ), lda,
257 $ b( k+kb, k ), ldb, one,
258 $ a( k+kb, k+kb ), lda )
259 CALL chemm(
'Right', uplo, n-k-kb+1, kb, -half,
260 $ a( k, k ), lda, b( k+kb, k ), ldb,
261 $ cone, a( k+kb, k ), lda )
262 CALL ctrsm(
'Left', uplo,
'No transpose',
263 $ 'Non-unit', n-k-kb+1, kb, cone,
264 $ b( k+kb, k+kb ), ldb, a( k+kb, k ),
265 $ lda )
266 END IF
267 20 CONTINUE
268 END IF
269 ELSE
270 IF( upper ) THEN
271
272
273
274 DO 30 k = 1, n, nb
275 kb = min( n-k+1, nb )
276
277
278
279 CALL ctrmm(
'Left', uplo,
'No transpose',
'Non-unit',
280 $ k-1, kb, cone, b, ldb, a( 1, k ), lda )
281 CALL chemm(
'Right', uplo, k-1, kb, half, a( k, k ),
282 $ lda, b( 1, k ), ldb, cone, a( 1, k ),
283 $ lda )
284 CALL cher2k( uplo,
'No transpose', k-1, kb, cone,
285 $ a( 1, k ), lda, b( 1, k ), ldb, one, a,
286 $ lda )
287 CALL chemm(
'Right', uplo, k-1, kb, half, a( k, k ),
288 $ lda, b( 1, k ), ldb, cone, a( 1, k ),
289 $ lda )
290 CALL ctrmm(
'Right', uplo,
'Conjugate transpose',
291 $ 'Non-unit', k-1, kb, cone, b( k, k ), ldb,
292 $ a( 1, k ), lda )
293 CALL chegs2( itype, uplo, kb, a( k, k ), lda,
294 $ b( k, k ), ldb, info )
295 30 CONTINUE
296 ELSE
297
298
299
300 DO 40 k = 1, n, nb
301 kb = min( n-k+1, nb )
302
303
304
305 CALL ctrmm(
'Right', uplo,
'No transpose',
'Non-unit',
306 $ kb, k-1, cone, b, ldb, a( k, 1 ), lda )
307 CALL chemm(
'Left', uplo, kb, k-1, half, a( k, k ),
308 $ lda, b( k, 1 ), ldb, cone, a( k, 1 ),
309 $ lda )
310 CALL cher2k( uplo,
'Conjugate transpose', k-1, kb,
311 $ cone, a( k, 1 ), lda, b( k, 1 ), ldb,
312 $ one, a, lda )
313 CALL chemm(
'Left', uplo, kb, k-1, half, a( k, k ),
314 $ lda, b( k, 1 ), ldb, cone, a( k, 1 ),
315 $ lda )
316 CALL ctrmm(
'Left', uplo,
'Conjugate transpose',
317 $ 'Non-unit', kb, k-1, cone, b( k, k ), ldb,
318 $ a( k, 1 ), lda )
319 CALL chegs2( itype, uplo, kb, a( k, k ), lda,
320 $ b( k, k ), ldb, info )
321 40 CONTINUE
322 END IF
323 END IF
324 END IF
325 RETURN
326
327
328
subroutine xerbla(srname, info)
subroutine chegs2(itype, uplo, n, a, lda, b, ldb, info)
CHEGS2 reduces a Hermitian definite generalized eigenproblem to standard form, using the factorizatio...
subroutine chemm(side, uplo, m, n, alpha, a, lda, b, ldb, beta, c, ldc)
CHEMM
subroutine cher2k(uplo, trans, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
CHER2K
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
logical function lsame(ca, cb)
LSAME
subroutine ctrmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
CTRMM
subroutine ctrsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
CTRSM