156
157
158
159
160
161
162 CHARACTER DIAG, TRANS, UPLO
163 INTEGER INFO, LDA, LDB, N, NRHS
164
165
166 INTEGER IPIV( * )
167 COMPLEX A( LDA, * ), B( LDB, * )
168
169
170
171
172
173 COMPLEX CONE
174 parameter( cone = ( 1.0e+0, 0.0e+0 ) )
175
176
177 LOGICAL NOUNIT
178 INTEGER J, K, KP
179 COMPLEX D11, D12, D21, D22, T1, T2
180
181
182 LOGICAL LSAME
184
185
187
188
189 INTRINSIC abs, conjg, max
190
191
192
193
194
195 info = 0
196 IF( .NOT.
lsame( uplo,
'U' ) .AND. .NOT.
lsame( uplo,
'L' ) )
THEN
197 info = -1
198 ELSE IF( .NOT.
lsame( trans,
'N' ) .AND. .NOT.
lsame( trans,
'C' ) )
199 $ THEN
200 info = -2
201 ELSE IF( .NOT.
lsame( diag,
'U' ) .AND. .NOT.
lsame( diag,
'N' ) )
202 $ THEN
203 info = -3
204 ELSE IF( n.LT.0 ) THEN
205 info = -4
206 ELSE IF( lda.LT.max( 1, n ) ) THEN
207 info = -6
208 ELSE IF( ldb.LT.max( 1, n ) ) THEN
209 info = -9
210 END IF
211 IF( info.NE.0 ) THEN
212 CALL xerbla(
'CLAVHE_ROOK ', -info )
213 RETURN
214 END IF
215
216
217
218 IF( n.EQ.0 )
219 $ RETURN
220
221 nounit =
lsame( diag,
'N' )
222
223
224
225
226
227 IF(
lsame( trans,
'N' ) )
THEN
228
229
230
231
232 IF(
lsame( uplo,
'U' ) )
THEN
233
234
235
236 k = 1
237 10 CONTINUE
238 IF( k.GT.n )
239 $ GO TO 30
240 IF( ipiv( k ).GT.0 ) THEN
241
242
243
244
245
246 IF( nounit )
247 $
CALL cscal( nrhs, a( k, k ), b( k, 1 ), ldb )
248
249
250
251 IF( k.GT.1 ) THEN
252
253
254
255 CALL cgeru( k-1, nrhs, cone, a( 1, k ), 1, b( k, 1 ),
256 $ ldb, b( 1, 1 ), ldb )
257
258
259
260 kp = ipiv( k )
261 IF( kp.NE.k )
262 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
263 END IF
264 k = k + 1
265 ELSE
266
267
268
269
270
271 IF( nounit ) THEN
272 d11 = a( k, k )
273 d22 = a( k+1, k+1 )
274 d12 = a( k, k+1 )
275 d21 = conjg( d12 )
276 DO 20 j = 1, nrhs
277 t1 = b( k, j )
278 t2 = b( k+1, j )
279 b( k, j ) = d11*t1 + d12*t2
280 b( k+1, j ) = d21*t1 + d22*t2
281 20 CONTINUE
282 END IF
283
284
285
286 IF( k.GT.1 ) THEN
287
288
289
290 CALL cgeru( k-1, nrhs, cone, a( 1, k ), 1, b( k, 1 ),
291 $ ldb, b( 1, 1 ), ldb )
292 CALL cgeru( k-1, nrhs, cone, a( 1, k+1 ), 1,
293 $ b( k+1, 1 ), ldb, b( 1, 1 ), ldb )
294
295
296
297
298
299
300 kp = abs( ipiv( k ) )
301 IF( kp.NE.k )
302 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
303
304
305
306 kp = abs( ipiv( k+1 ) )
307 IF( kp.NE.k+1 )
308 $
CALL cswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ),
309 $ ldb )
310 END IF
311 k = k + 2
312 END IF
313 GO TO 10
314 30 CONTINUE
315
316
317
318
319 ELSE
320
321
322
323 k = n
324 40 CONTINUE
325 IF( k.LT.1 )
326 $ GO TO 60
327
328
329
330
331 IF( ipiv( k ).GT.0 ) THEN
332
333
334
335
336
337 IF( nounit )
338 $
CALL cscal( nrhs, a( k, k ), b( k, 1 ), ldb )
339
340
341
342 IF( k.NE.n ) THEN
343 kp = ipiv( k )
344
345
346
347 CALL cgeru( n-k, nrhs, cone, a( k+1, k ), 1,
348 $ b( k, 1 ), ldb, b( k+1, 1 ), ldb )
349
350
351
352
353 IF( kp.NE.k )
354 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
355 END IF
356 k = k - 1
357
358 ELSE
359
360
361
362
363
364 IF( nounit ) THEN
365 d11 = a( k-1, k-1 )
366 d22 = a( k, k )
367 d21 = a( k, k-1 )
368 d12 = conjg( d21 )
369 DO 50 j = 1, nrhs
370 t1 = b( k-1, j )
371 t2 = b( k, j )
372 b( k-1, j ) = d11*t1 + d12*t2
373 b( k, j ) = d21*t1 + d22*t2
374 50 CONTINUE
375 END IF
376
377
378
379 IF( k.NE.n ) THEN
380
381
382
383 CALL cgeru( n-k, nrhs, cone, a( k+1, k ), 1,
384 $ b( k, 1 ), ldb, b( k+1, 1 ), ldb )
385 CALL cgeru( n-k, nrhs, cone, a( k+1, k-1 ), 1,
386 $ b( k-1, 1 ), ldb, b( k+1, 1 ), ldb )
387
388
389
390
391
392
393
394 kp = abs( ipiv( k ) )
395 IF( kp.NE.k )
396 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
397
398
399
400 kp = abs( ipiv( k-1 ) )
401 IF( kp.NE.k-1 )
402 $
CALL cswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ),
403 $ ldb )
404
405 END IF
406 k = k - 2
407 END IF
408 GO TO 40
409 60 CONTINUE
410 END IF
411
412
413
414
415
416 ELSE
417
418
419
420
421
422 IF(
lsame( uplo,
'U' ) )
THEN
423
424
425
426 k = n
427 70 IF( k.LT.1 )
428 $ GO TO 90
429
430
431
432 IF( ipiv( k ).GT.0 ) THEN
433 IF( k.GT.1 ) THEN
434
435
436
437 kp = ipiv( k )
438 IF( kp.NE.k )
439 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
440
441
442
443
444
445 CALL clacgv( nrhs, b( k, 1 ), ldb )
446 CALL cgemv(
'Conjugate', k-1, nrhs, cone, b, ldb,
447 $ a( 1, k ), 1, cone, b( k, 1 ), ldb )
448 CALL clacgv( nrhs, b( k, 1 ), ldb )
449 END IF
450 IF( nounit )
451 $
CALL cscal( nrhs, a( k, k ), b( k, 1 ), ldb )
452 k = k - 1
453
454
455
456 ELSE
457 IF( k.GT.2 ) THEN
458
459
460
461 kp = abs( ipiv( k ) )
462 IF( kp.NE.k )
463 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
464
465
466
467 kp = abs( ipiv( k-1 ) )
468 IF( kp.NE.k-1 )
469 $
CALL cswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ),
470 $ ldb )
471
472
473
474
475
476
477 CALL clacgv( nrhs, b( k, 1 ), ldb )
478 CALL cgemv(
'Conjugate', k-2, nrhs, cone, b, ldb,
479 $ a( 1, k ), 1, cone, b( k, 1 ), ldb )
480 CALL clacgv( nrhs, b( k, 1 ), ldb )
481
482 CALL clacgv( nrhs, b( k-1, 1 ), ldb )
483 CALL cgemv(
'Conjugate', k-2, nrhs, cone, b, ldb,
484 $ a( 1, k-1 ), 1, cone, b( k-1, 1 ), ldb )
485 CALL clacgv( nrhs, b( k-1, 1 ), ldb )
486 END IF
487
488
489
490 IF( nounit ) THEN
491 d11 = a( k-1, k-1 )
492 d22 = a( k, k )
493 d12 = a( k-1, k )
494 d21 = conjg( d12 )
495 DO 80 j = 1, nrhs
496 t1 = b( k-1, j )
497 t2 = b( k, j )
498 b( k-1, j ) = d11*t1 + d12*t2
499 b( k, j ) = d21*t1 + d22*t2
500 80 CONTINUE
501 END IF
502 k = k - 2
503 END IF
504 GO TO 70
505 90 CONTINUE
506
507
508
509
510
511 ELSE
512
513
514
515 k = 1
516 100 CONTINUE
517 IF( k.GT.n )
518 $ GO TO 120
519
520
521
522 IF( ipiv( k ).GT.0 ) THEN
523 IF( k.LT.n ) THEN
524
525
526
527 kp = ipiv( k )
528 IF( kp.NE.k )
529 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
530
531
532
533 CALL clacgv( nrhs, b( k, 1 ), ldb )
534 CALL cgemv(
'Conjugate', n-k, nrhs, cone, b( k+1, 1 ),
535 $ ldb, a( k+1, k ), 1, cone, b( k, 1 ), ldb )
536 CALL clacgv( nrhs, b( k, 1 ), ldb )
537 END IF
538 IF( nounit )
539 $
CALL cscal( nrhs, a( k, k ), b( k, 1 ), ldb )
540 k = k + 1
541
542
543
544 ELSE
545 IF( k.LT.n-1 ) THEN
546
547
548
549 kp = abs( ipiv( k ) )
550 IF( kp.NE.k )
551 $
CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
552
553
554
555 kp = abs( ipiv( k+1 ) )
556 IF( kp.NE.k+1 )
557 $
CALL cswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ),
558 $ ldb )
559
560
561
562 CALL clacgv( nrhs, b( k+1, 1 ), ldb )
563 CALL cgemv(
'Conjugate', n-k-1, nrhs, cone,
564 $ b( k+2, 1 ), ldb, a( k+2, k+1 ), 1, cone,
565 $ b( k+1, 1 ), ldb )
566 CALL clacgv( nrhs, b( k+1, 1 ), ldb )
567
568 CALL clacgv( nrhs, b( k, 1 ), ldb )
569 CALL cgemv(
'Conjugate', n-k-1, nrhs, cone,
570 $ b( k+2, 1 ), ldb, a( k+2, k ), 1, cone,
571 $ b( k, 1 ), ldb )
572 CALL clacgv( nrhs, b( k, 1 ), ldb )
573 END IF
574
575
576
577 IF( nounit ) THEN
578 d11 = a( k, k )
579 d22 = a( k+1, k+1 )
580 d21 = a( k+1, k )
581 d12 = conjg( d21 )
582 DO 110 j = 1, nrhs
583 t1 = b( k, j )
584 t2 = b( k+1, j )
585 b( k, j ) = d11*t1 + d12*t2
586 b( k+1, j ) = d21*t1 + d22*t2
587 110 CONTINUE
588 END IF
589 k = k + 2
590 END IF
591 GO TO 100
592 120 CONTINUE
593 END IF
594
595 END IF
596 RETURN
597
598
599
subroutine xerbla(srname, info)
subroutine cgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
CGEMV
subroutine cgeru(m, n, alpha, x, incx, y, incy, a, lda)
CGERU
subroutine clacgv(n, x, incx)
CLACGV conjugates a complex vector.
logical function lsame(ca, cb)
LSAME
subroutine cscal(n, ca, cx, incx)
CSCAL
subroutine cswap(n, cx, incx, cy, incy)
CSWAP