LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ clanhf()

real function clanhf ( character  norm,
character  transr,
character  uplo,
integer  n,
complex, dimension( 0: * )  a,
real, dimension( 0: * )  work 
)

CLANHF returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a Hermitian matrix in RFP format.

Download CLANHF + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 CLANHF  returns the value of the one norm,  or the Frobenius norm, or
 the  infinity norm,  or the  element of  largest absolute value  of a
 complex Hermitian matrix A in RFP format.
Returns
CLANHF 
    CLANHF = ( max(abs(A(i,j))), NORM = 'M' or 'm'
             (
             ( norm1(A),         NORM = '1', 'O' or 'o'
             (
             ( normI(A),         NORM = 'I' or 'i'
             (
             ( normF(A),         NORM = 'F', 'f', 'E' or 'e'

 where  norm1  denotes the  one norm of a matrix (maximum column sum),
 normI  denotes the  infinity norm  of a matrix  (maximum row sum) and
 normF  denotes the  Frobenius norm of a matrix (square root of sum of
 squares).  Note that  max(abs(A(i,j)))  is not a  matrix norm.
Parameters
[in]NORM
          NORM is CHARACTER
            Specifies the value to be returned in CLANHF as described
            above.
[in]TRANSR
          TRANSR is CHARACTER
            Specifies whether the RFP format of A is normal or
            conjugate-transposed format.
            = 'N':  RFP format is Normal
            = 'C':  RFP format is Conjugate-transposed
[in]UPLO
          UPLO is CHARACTER
            On entry, UPLO specifies whether the RFP matrix A came from
            an upper or lower triangular matrix as follows:

            UPLO = 'U' or 'u' RFP A came from an upper triangular
            matrix

            UPLO = 'L' or 'l' RFP A came from a  lower triangular
            matrix
[in]N
          N is INTEGER
            The order of the matrix A.  N >= 0.  When N = 0, CLANHF is
            set to zero.
[in]A
          A is COMPLEX array, dimension ( N*(N+1)/2 );
            On entry, the matrix A in RFP Format.
            RFP Format is described by TRANSR, UPLO and N as follows:
            If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
            K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
            TRANSR = 'C' then RFP is the Conjugate-transpose of RFP A
            as defined when TRANSR = 'N'. The contents of RFP A are
            defined by UPLO as follows: If UPLO = 'U' the RFP A
            contains the ( N*(N+1)/2 ) elements of upper packed A
            either in normal or conjugate-transpose Format. If
            UPLO = 'L' the RFP A contains the ( N*(N+1) /2 ) elements
            of lower packed A either in normal or conjugate-transpose
            Format. The LDA of RFP A is (N+1)/2 when TRANSR = 'C'. When
            TRANSR is 'N' the LDA is N+1 when N is even and is N when
            is odd. See the Note below for more details.
            Unchanged on exit.
[out]WORK
          WORK is REAL array, dimension (LWORK),
            where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
            WORK is not referenced.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  We first consider Standard Packed Format when N is even.
  We give an example where N = 6.

      AP is Upper             AP is Lower

   00 01 02 03 04 05       00
      11 12 13 14 15       10 11
         22 23 24 25       20 21 22
            33 34 35       30 31 32 33
               44 45       40 41 42 43 44
                  55       50 51 52 53 54 55


  Let TRANSR = 'N'. RFP holds AP as follows:
  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
  three columns of AP upper. The lower triangle A(4:6,0:2) consists of
  conjugate-transpose of the first three columns of AP upper.
  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
  three columns of AP lower. The upper triangle A(0:2,0:2) consists of
  conjugate-transpose of the last three columns of AP lower.
  To denote conjugate we place -- above the element. This covers the
  case N even and TRANSR = 'N'.

         RFP A                   RFP A

                                -- -- --
        03 04 05                33 43 53
                                   -- --
        13 14 15                00 44 54
                                      --
        23 24 25                10 11 55

        33 34 35                20 21 22
        --
        00 44 45                30 31 32
        -- --
        01 11 55                40 41 42
        -- -- --
        02 12 22                50 51 52

  Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
  transpose of RFP A above. One therefore gets:


           RFP A                   RFP A

     -- -- -- --                -- -- -- -- -- --
     03 13 23 33 00 01 02    33 00 10 20 30 40 50
     -- -- -- -- --                -- -- -- -- --
     04 14 24 34 44 11 12    43 44 11 21 31 41 51
     -- -- -- -- -- --                -- -- -- --
     05 15 25 35 45 55 22    53 54 55 22 32 42 52


  We next  consider Standard Packed Format when N is odd.
  We give an example where N = 5.

     AP is Upper                 AP is Lower

   00 01 02 03 04              00
      11 12 13 14              10 11
         22 23 24              20 21 22
            33 34              30 31 32 33
               44              40 41 42 43 44


  Let TRANSR = 'N'. RFP holds AP as follows:
  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
  conjugate-transpose of the first two   columns of AP upper.
  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
  conjugate-transpose of the last two   columns of AP lower.
  To denote conjugate we place -- above the element. This covers the
  case N odd  and TRANSR = 'N'.

         RFP A                   RFP A

                                   -- --
        02 03 04                00 33 43
                                      --
        12 13 14                10 11 44

        22 23 24                20 21 22
        --
        00 33 34                30 31 32
        -- --
        01 11 44                40 41 42

  Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
  transpose of RFP A above. One therefore gets:


           RFP A                   RFP A

     -- -- --                   -- -- -- -- -- --
     02 12 22 00 01             00 10 20 30 40 50
     -- -- -- --                   -- -- -- -- --
     03 13 23 33 11             33 11 21 31 41 51
     -- -- -- -- --                   -- -- -- --
     04 14 24 34 44             43 44 22 32 42 52

Definition at line 245 of file clanhf.f.

246*
247* -- LAPACK computational routine --
248* -- LAPACK is a software package provided by Univ. of Tennessee, --
249* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
250*
251* .. Scalar Arguments ..
252 CHARACTER NORM, TRANSR, UPLO
253 INTEGER N
254* ..
255* .. Array Arguments ..
256 REAL WORK( 0: * )
257 COMPLEX A( 0: * )
258* ..
259*
260* =====================================================================
261*
262* .. Parameters ..
263 REAL ONE, ZERO
264 parameter( one = 1.0e+0, zero = 0.0e+0 )
265* ..
266* .. Local Scalars ..
267 INTEGER I, J, IFM, ILU, NOE, N1, K, L, LDA
268 REAL SCALE, S, VALUE, AA, TEMP
269* ..
270* .. External Functions ..
271 LOGICAL LSAME, SISNAN
272 EXTERNAL lsame, sisnan
273* ..
274* .. External Subroutines ..
275 EXTERNAL classq
276* ..
277* .. Intrinsic Functions ..
278 INTRINSIC abs, real, sqrt
279* ..
280* .. Executable Statements ..
281*
282 IF( n.EQ.0 ) THEN
283 clanhf = zero
284 RETURN
285 ELSE IF( n.EQ.1 ) THEN
286 clanhf = abs(real(a(0)))
287 RETURN
288 END IF
289*
290* set noe = 1 if n is odd. if n is even set noe=0
291*
292 noe = 1
293 IF( mod( n, 2 ).EQ.0 )
294 $ noe = 0
295*
296* set ifm = 0 when form='C' or 'c' and 1 otherwise
297*
298 ifm = 1
299 IF( lsame( transr, 'C' ) )
300 $ ifm = 0
301*
302* set ilu = 0 when uplo='U or 'u' and 1 otherwise
303*
304 ilu = 1
305 IF( lsame( uplo, 'U' ) )
306 $ ilu = 0
307*
308* set lda = (n+1)/2 when ifm = 0
309* set lda = n when ifm = 1 and noe = 1
310* set lda = n+1 when ifm = 1 and noe = 0
311*
312 IF( ifm.EQ.1 ) THEN
313 IF( noe.EQ.1 ) THEN
314 lda = n
315 ELSE
316* noe=0
317 lda = n + 1
318 END IF
319 ELSE
320* ifm=0
321 lda = ( n+1 ) / 2
322 END IF
323*
324 IF( lsame( norm, 'M' ) ) THEN
325*
326* Find max(abs(A(i,j))).
327*
328 k = ( n+1 ) / 2
329 VALUE = zero
330 IF( noe.EQ.1 ) THEN
331* n is odd & n = k + k - 1
332 IF( ifm.EQ.1 ) THEN
333* A is n by k
334 IF( ilu.EQ.1 ) THEN
335* uplo ='L'
336 j = 0
337* -> L(0,0)
338 temp = abs( real( a( j+j*lda ) ) )
339 IF( VALUE .LT. temp .OR. sisnan( temp ) )
340 $ VALUE = temp
341 DO i = 1, n - 1
342 temp = abs( a( i+j*lda ) )
343 IF( VALUE .LT. temp .OR. sisnan( temp ) )
344 $ VALUE = temp
345 END DO
346 DO j = 1, k - 1
347 DO i = 0, j - 2
348 temp = abs( a( i+j*lda ) )
349 IF( VALUE .LT. temp .OR. sisnan( temp ) )
350 $ VALUE = temp
351 END DO
352 i = j - 1
353* L(k+j,k+j)
354 temp = abs( real( a( i+j*lda ) ) )
355 IF( VALUE .LT. temp .OR. sisnan( temp ) )
356 $ VALUE = temp
357 i = j
358* -> L(j,j)
359 temp = abs( real( a( i+j*lda ) ) )
360 IF( VALUE .LT. temp .OR. sisnan( temp ) )
361 $ VALUE = temp
362 DO i = j + 1, n - 1
363 temp = abs( a( i+j*lda ) )
364 IF( VALUE .LT. temp .OR. sisnan( temp ) )
365 $ VALUE = temp
366 END DO
367 END DO
368 ELSE
369* uplo = 'U'
370 DO j = 0, k - 2
371 DO i = 0, k + j - 2
372 temp = abs( a( i+j*lda ) )
373 IF( VALUE .LT. temp .OR. sisnan( temp ) )
374 $ VALUE = temp
375 END DO
376 i = k + j - 1
377* -> U(i,i)
378 temp = abs( real( a( i+j*lda ) ) )
379 IF( VALUE .LT. temp .OR. sisnan( temp ) )
380 $ VALUE = temp
381 i = i + 1
382* =k+j; i -> U(j,j)
383 temp = abs( real( a( i+j*lda ) ) )
384 IF( VALUE .LT. temp .OR. sisnan( temp ) )
385 $ VALUE = temp
386 DO i = k + j + 1, n - 1
387 temp = abs( a( i+j*lda ) )
388 IF( VALUE .LT. temp .OR. sisnan( temp ) )
389 $ VALUE = temp
390 END DO
391 END DO
392 DO i = 0, n - 2
393 temp = abs( a( i+j*lda ) )
394 IF( VALUE .LT. temp .OR. sisnan( temp ) )
395 $ VALUE = temp
396* j=k-1
397 END DO
398* i=n-1 -> U(n-1,n-1)
399 temp = abs( real( a( i+j*lda ) ) )
400 IF( VALUE .LT. temp .OR. sisnan( temp ) )
401 $ VALUE = temp
402 END IF
403 ELSE
404* xpose case; A is k by n
405 IF( ilu.EQ.1 ) THEN
406* uplo ='L'
407 DO j = 0, k - 2
408 DO i = 0, j - 1
409 temp = abs( a( i+j*lda ) )
410 IF( VALUE .LT. temp .OR. sisnan( temp ) )
411 $ VALUE = temp
412 END DO
413 i = j
414* L(i,i)
415 temp = abs( real( a( i+j*lda ) ) )
416 IF( VALUE .LT. temp .OR. sisnan( temp ) )
417 $ VALUE = temp
418 i = j + 1
419* L(j+k,j+k)
420 temp = abs( real( a( i+j*lda ) ) )
421 IF( VALUE .LT. temp .OR. sisnan( temp ) )
422 $ VALUE = temp
423 DO i = j + 2, k - 1
424 temp = abs( a( i+j*lda ) )
425 IF( VALUE .LT. temp .OR. sisnan( temp ) )
426 $ VALUE = temp
427 END DO
428 END DO
429 j = k - 1
430 DO i = 0, k - 2
431 temp = abs( a( i+j*lda ) )
432 IF( VALUE .LT. temp .OR. sisnan( temp ) )
433 $ VALUE = temp
434 END DO
435 i = k - 1
436* -> L(i,i) is at A(i,j)
437 temp = abs( real( a( i+j*lda ) ) )
438 IF( VALUE .LT. temp .OR. sisnan( temp ) )
439 $ VALUE = temp
440 DO j = k, n - 1
441 DO i = 0, k - 1
442 temp = abs( a( i+j*lda ) )
443 IF( VALUE .LT. temp .OR. sisnan( temp ) )
444 $ VALUE = temp
445 END DO
446 END DO
447 ELSE
448* uplo = 'U'
449 DO j = 0, k - 2
450 DO i = 0, k - 1
451 temp = abs( a( i+j*lda ) )
452 IF( VALUE .LT. temp .OR. sisnan( temp ) )
453 $ VALUE = temp
454 END DO
455 END DO
456 j = k - 1
457* -> U(j,j) is at A(0,j)
458 temp = abs( real( a( 0+j*lda ) ) )
459 IF( VALUE .LT. temp .OR. sisnan( temp ) )
460 $ VALUE = temp
461 DO i = 1, k - 1
462 temp = abs( a( i+j*lda ) )
463 IF( VALUE .LT. temp .OR. sisnan( temp ) )
464 $ VALUE = temp
465 END DO
466 DO j = k, n - 1
467 DO i = 0, j - k - 1
468 temp = abs( a( i+j*lda ) )
469 IF( VALUE .LT. temp .OR. sisnan( temp ) )
470 $ VALUE = temp
471 END DO
472 i = j - k
473* -> U(i,i) at A(i,j)
474 temp = abs( real( a( i+j*lda ) ) )
475 IF( VALUE .LT. temp .OR. sisnan( temp ) )
476 $ VALUE = temp
477 i = j - k + 1
478* U(j,j)
479 temp = abs( real( a( i+j*lda ) ) )
480 IF( VALUE .LT. temp .OR. sisnan( temp ) )
481 $ VALUE = temp
482 DO i = j - k + 2, k - 1
483 temp = abs( a( i+j*lda ) )
484 IF( VALUE .LT. temp .OR. sisnan( temp ) )
485 $ VALUE = temp
486 END DO
487 END DO
488 END IF
489 END IF
490 ELSE
491* n is even & k = n/2
492 IF( ifm.EQ.1 ) THEN
493* A is n+1 by k
494 IF( ilu.EQ.1 ) THEN
495* uplo ='L'
496 j = 0
497* -> L(k,k) & j=1 -> L(0,0)
498 temp = abs( real( a( j+j*lda ) ) )
499 IF( VALUE .LT. temp .OR. sisnan( temp ) )
500 $ VALUE = temp
501 temp = abs( real( a( j+1+j*lda ) ) )
502 IF( VALUE .LT. temp .OR. sisnan( temp ) )
503 $ VALUE = temp
504 DO i = 2, n
505 temp = abs( a( i+j*lda ) )
506 IF( VALUE .LT. temp .OR. sisnan( temp ) )
507 $ VALUE = temp
508 END DO
509 DO j = 1, k - 1
510 DO i = 0, j - 1
511 temp = abs( a( i+j*lda ) )
512 IF( VALUE .LT. temp .OR. sisnan( temp ) )
513 $ VALUE = temp
514 END DO
515 i = j
516* L(k+j,k+j)
517 temp = abs( real( a( i+j*lda ) ) )
518 IF( VALUE .LT. temp .OR. sisnan( temp ) )
519 $ VALUE = temp
520 i = j + 1
521* -> L(j,j)
522 temp = abs( real( a( i+j*lda ) ) )
523 IF( VALUE .LT. temp .OR. sisnan( temp ) )
524 $ VALUE = temp
525 DO i = j + 2, n
526 temp = abs( a( i+j*lda ) )
527 IF( VALUE .LT. temp .OR. sisnan( temp ) )
528 $ VALUE = temp
529 END DO
530 END DO
531 ELSE
532* uplo = 'U'
533 DO j = 0, k - 2
534 DO i = 0, k + j - 1
535 temp = abs( a( i+j*lda ) )
536 IF( VALUE .LT. temp .OR. sisnan( temp ) )
537 $ VALUE = temp
538 END DO
539 i = k + j
540* -> U(i,i)
541 temp = abs( real( a( i+j*lda ) ) )
542 IF( VALUE .LT. temp .OR. sisnan( temp ) )
543 $ VALUE = temp
544 i = i + 1
545* =k+j+1; i -> U(j,j)
546 temp = abs( real( a( i+j*lda ) ) )
547 IF( VALUE .LT. temp .OR. sisnan( temp ) )
548 $ VALUE = temp
549 DO i = k + j + 2, n
550 temp = abs( a( i+j*lda ) )
551 IF( VALUE .LT. temp .OR. sisnan( temp ) )
552 $ VALUE = temp
553 END DO
554 END DO
555 DO i = 0, n - 2
556 temp = abs( a( i+j*lda ) )
557 IF( VALUE .LT. temp .OR. sisnan( temp ) )
558 $ VALUE = temp
559* j=k-1
560 END DO
561* i=n-1 -> U(n-1,n-1)
562 temp = abs( real( a( i+j*lda ) ) )
563 IF( VALUE .LT. temp .OR. sisnan( temp ) )
564 $ VALUE = temp
565 i = n
566* -> U(k-1,k-1)
567 temp = abs( real( a( i+j*lda ) ) )
568 IF( VALUE .LT. temp .OR. sisnan( temp ) )
569 $ VALUE = temp
570 END IF
571 ELSE
572* xpose case; A is k by n+1
573 IF( ilu.EQ.1 ) THEN
574* uplo ='L'
575 j = 0
576* -> L(k,k) at A(0,0)
577 temp = abs( real( a( j+j*lda ) ) )
578 IF( VALUE .LT. temp .OR. sisnan( temp ) )
579 $ VALUE = temp
580 DO i = 1, k - 1
581 temp = abs( a( i+j*lda ) )
582 IF( VALUE .LT. temp .OR. sisnan( temp ) )
583 $ VALUE = temp
584 END DO
585 DO j = 1, k - 1
586 DO i = 0, j - 2
587 temp = abs( a( i+j*lda ) )
588 IF( VALUE .LT. temp .OR. sisnan( temp ) )
589 $ VALUE = temp
590 END DO
591 i = j - 1
592* L(i,i)
593 temp = abs( real( a( i+j*lda ) ) )
594 IF( VALUE .LT. temp .OR. sisnan( temp ) )
595 $ VALUE = temp
596 i = j
597* L(j+k,j+k)
598 temp = abs( real( a( i+j*lda ) ) )
599 IF( VALUE .LT. temp .OR. sisnan( temp ) )
600 $ VALUE = temp
601 DO i = j + 1, k - 1
602 temp = abs( a( i+j*lda ) )
603 IF( VALUE .LT. temp .OR. sisnan( temp ) )
604 $ VALUE = temp
605 END DO
606 END DO
607 j = k
608 DO i = 0, k - 2
609 temp = abs( a( i+j*lda ) )
610 IF( VALUE .LT. temp .OR. sisnan( temp ) )
611 $ VALUE = temp
612 END DO
613 i = k - 1
614* -> L(i,i) is at A(i,j)
615 temp = abs( real( a( i+j*lda ) ) )
616 IF( VALUE .LT. temp .OR. sisnan( temp ) )
617 $ VALUE = temp
618 DO j = k + 1, n
619 DO i = 0, k - 1
620 temp = abs( a( i+j*lda ) )
621 IF( VALUE .LT. temp .OR. sisnan( temp ) )
622 $ VALUE = temp
623 END DO
624 END DO
625 ELSE
626* uplo = 'U'
627 DO j = 0, k - 1
628 DO i = 0, k - 1
629 temp = abs( a( i+j*lda ) )
630 IF( VALUE .LT. temp .OR. sisnan( temp ) )
631 $ VALUE = temp
632 END DO
633 END DO
634 j = k
635* -> U(j,j) is at A(0,j)
636 temp = abs( real( a( 0+j*lda ) ) )
637 IF( VALUE .LT. temp .OR. sisnan( temp ) )
638 $ VALUE = temp
639 DO i = 1, k - 1
640 temp = abs( a( i+j*lda ) )
641 IF( VALUE .LT. temp .OR. sisnan( temp ) )
642 $ VALUE = temp
643 END DO
644 DO j = k + 1, n - 1
645 DO i = 0, j - k - 2
646 temp = abs( a( i+j*lda ) )
647 IF( VALUE .LT. temp .OR. sisnan( temp ) )
648 $ VALUE = temp
649 END DO
650 i = j - k - 1
651* -> U(i,i) at A(i,j)
652 temp = abs( real( a( i+j*lda ) ) )
653 IF( VALUE .LT. temp .OR. sisnan( temp ) )
654 $ VALUE = temp
655 i = j - k
656* U(j,j)
657 temp = abs( real( a( i+j*lda ) ) )
658 IF( VALUE .LT. temp .OR. sisnan( temp ) )
659 $ VALUE = temp
660 DO i = j - k + 1, k - 1
661 temp = abs( a( i+j*lda ) )
662 IF( VALUE .LT. temp .OR. sisnan( temp ) )
663 $ VALUE = temp
664 END DO
665 END DO
666 j = n
667 DO i = 0, k - 2
668 temp = abs( a( i+j*lda ) )
669 IF( VALUE .LT. temp .OR. sisnan( temp ) )
670 $ VALUE = temp
671 END DO
672 i = k - 1
673* U(k,k) at A(i,j)
674 temp = abs( real( a( i+j*lda ) ) )
675 IF( VALUE .LT. temp .OR. sisnan( temp ) )
676 $ VALUE = temp
677 END IF
678 END IF
679 END IF
680 ELSE IF( ( lsame( norm, 'I' ) ) .OR. ( lsame( norm, 'O' ) ) .OR.
681 $ ( norm.EQ.'1' ) ) THEN
682*
683* Find normI(A) ( = norm1(A), since A is Hermitian).
684*
685 IF( ifm.EQ.1 ) THEN
686* A is 'N'
687 k = n / 2
688 IF( noe.EQ.1 ) THEN
689* n is odd & A is n by (n+1)/2
690 IF( ilu.EQ.0 ) THEN
691* uplo = 'U'
692 DO i = 0, k - 1
693 work( i ) = zero
694 END DO
695 DO j = 0, k
696 s = zero
697 DO i = 0, k + j - 1
698 aa = abs( a( i+j*lda ) )
699* -> A(i,j+k)
700 s = s + aa
701 work( i ) = work( i ) + aa
702 END DO
703 aa = abs( real( a( i+j*lda ) ) )
704* -> A(j+k,j+k)
705 work( j+k ) = s + aa
706 IF( i.EQ.k+k )
707 $ GO TO 10
708 i = i + 1
709 aa = abs( real( a( i+j*lda ) ) )
710* -> A(j,j)
711 work( j ) = work( j ) + aa
712 s = zero
713 DO l = j + 1, k - 1
714 i = i + 1
715 aa = abs( a( i+j*lda ) )
716* -> A(l,j)
717 s = s + aa
718 work( l ) = work( l ) + aa
719 END DO
720 work( j ) = work( j ) + s
721 END DO
722 10 CONTINUE
723 VALUE = work( 0 )
724 DO i = 1, n-1
725 temp = work( i )
726 IF( VALUE .LT. temp .OR. sisnan( temp ) )
727 $ VALUE = temp
728 END DO
729 ELSE
730* ilu = 1 & uplo = 'L'
731 k = k + 1
732* k=(n+1)/2 for n odd and ilu=1
733 DO i = k, n - 1
734 work( i ) = zero
735 END DO
736 DO j = k - 1, 0, -1
737 s = zero
738 DO i = 0, j - 2
739 aa = abs( a( i+j*lda ) )
740* -> A(j+k,i+k)
741 s = s + aa
742 work( i+k ) = work( i+k ) + aa
743 END DO
744 IF( j.GT.0 ) THEN
745 aa = abs( real( a( i+j*lda ) ) )
746* -> A(j+k,j+k)
747 s = s + aa
748 work( i+k ) = work( i+k ) + s
749* i=j
750 i = i + 1
751 END IF
752 aa = abs( real( a( i+j*lda ) ) )
753* -> A(j,j)
754 work( j ) = aa
755 s = zero
756 DO l = j + 1, n - 1
757 i = i + 1
758 aa = abs( a( i+j*lda ) )
759* -> A(l,j)
760 s = s + aa
761 work( l ) = work( l ) + aa
762 END DO
763 work( j ) = work( j ) + s
764 END DO
765 VALUE = work( 0 )
766 DO i = 1, n-1
767 temp = work( i )
768 IF( VALUE .LT. temp .OR. sisnan( temp ) )
769 $ VALUE = temp
770 END DO
771 END IF
772 ELSE
773* n is even & A is n+1 by k = n/2
774 IF( ilu.EQ.0 ) THEN
775* uplo = 'U'
776 DO i = 0, k - 1
777 work( i ) = zero
778 END DO
779 DO j = 0, k - 1
780 s = zero
781 DO i = 0, k + j - 1
782 aa = abs( a( i+j*lda ) )
783* -> A(i,j+k)
784 s = s + aa
785 work( i ) = work( i ) + aa
786 END DO
787 aa = abs( real( a( i+j*lda ) ) )
788* -> A(j+k,j+k)
789 work( j+k ) = s + aa
790 i = i + 1
791 aa = abs( real( a( i+j*lda ) ) )
792* -> A(j,j)
793 work( j ) = work( j ) + aa
794 s = zero
795 DO l = j + 1, k - 1
796 i = i + 1
797 aa = abs( a( i+j*lda ) )
798* -> A(l,j)
799 s = s + aa
800 work( l ) = work( l ) + aa
801 END DO
802 work( j ) = work( j ) + s
803 END DO
804 VALUE = work( 0 )
805 DO i = 1, n-1
806 temp = work( i )
807 IF( VALUE .LT. temp .OR. sisnan( temp ) )
808 $ VALUE = temp
809 END DO
810 ELSE
811* ilu = 1 & uplo = 'L'
812 DO i = k, n - 1
813 work( i ) = zero
814 END DO
815 DO j = k - 1, 0, -1
816 s = zero
817 DO i = 0, j - 1
818 aa = abs( a( i+j*lda ) )
819* -> A(j+k,i+k)
820 s = s + aa
821 work( i+k ) = work( i+k ) + aa
822 END DO
823 aa = abs( real( a( i+j*lda ) ) )
824* -> A(j+k,j+k)
825 s = s + aa
826 work( i+k ) = work( i+k ) + s
827* i=j
828 i = i + 1
829 aa = abs( real( a( i+j*lda ) ) )
830* -> A(j,j)
831 work( j ) = aa
832 s = zero
833 DO l = j + 1, n - 1
834 i = i + 1
835 aa = abs( a( i+j*lda ) )
836* -> A(l,j)
837 s = s + aa
838 work( l ) = work( l ) + aa
839 END DO
840 work( j ) = work( j ) + s
841 END DO
842 VALUE = work( 0 )
843 DO i = 1, n-1
844 temp = work( i )
845 IF( VALUE .LT. temp .OR. sisnan( temp ) )
846 $ VALUE = temp
847 END DO
848 END IF
849 END IF
850 ELSE
851* ifm=0
852 k = n / 2
853 IF( noe.EQ.1 ) THEN
854* n is odd & A is (n+1)/2 by n
855 IF( ilu.EQ.0 ) THEN
856* uplo = 'U'
857 n1 = k
858* n/2
859 k = k + 1
860* k is the row size and lda
861 DO i = n1, n - 1
862 work( i ) = zero
863 END DO
864 DO j = 0, n1 - 1
865 s = zero
866 DO i = 0, k - 1
867 aa = abs( a( i+j*lda ) )
868* A(j,n1+i)
869 work( i+n1 ) = work( i+n1 ) + aa
870 s = s + aa
871 END DO
872 work( j ) = s
873 END DO
874* j=n1=k-1 is special
875 s = abs( real( a( 0+j*lda ) ) )
876* A(k-1,k-1)
877 DO i = 1, k - 1
878 aa = abs( a( i+j*lda ) )
879* A(k-1,i+n1)
880 work( i+n1 ) = work( i+n1 ) + aa
881 s = s + aa
882 END DO
883 work( j ) = work( j ) + s
884 DO j = k, n - 1
885 s = zero
886 DO i = 0, j - k - 1
887 aa = abs( a( i+j*lda ) )
888* A(i,j-k)
889 work( i ) = work( i ) + aa
890 s = s + aa
891 END DO
892* i=j-k
893 aa = abs( real( a( i+j*lda ) ) )
894* A(j-k,j-k)
895 s = s + aa
896 work( j-k ) = work( j-k ) + s
897 i = i + 1
898 s = abs( real( a( i+j*lda ) ) )
899* A(j,j)
900 DO l = j + 1, n - 1
901 i = i + 1
902 aa = abs( a( i+j*lda ) )
903* A(j,l)
904 work( l ) = work( l ) + aa
905 s = s + aa
906 END DO
907 work( j ) = work( j ) + s
908 END DO
909 VALUE = work( 0 )
910 DO i = 1, n-1
911 temp = work( i )
912 IF( VALUE .LT. temp .OR. sisnan( temp ) )
913 $ VALUE = temp
914 END DO
915 ELSE
916* ilu=1 & uplo = 'L'
917 k = k + 1
918* k=(n+1)/2 for n odd and ilu=1
919 DO i = k, n - 1
920 work( i ) = zero
921 END DO
922 DO j = 0, k - 2
923* process
924 s = zero
925 DO i = 0, j - 1
926 aa = abs( a( i+j*lda ) )
927* A(j,i)
928 work( i ) = work( i ) + aa
929 s = s + aa
930 END DO
931 aa = abs( real( a( i+j*lda ) ) )
932* i=j so process of A(j,j)
933 s = s + aa
934 work( j ) = s
935* is initialised here
936 i = i + 1
937* i=j process A(j+k,j+k)
938 aa = abs( real( a( i+j*lda ) ) )
939 s = aa
940 DO l = k + j + 1, n - 1
941 i = i + 1
942 aa = abs( a( i+j*lda ) )
943* A(l,k+j)
944 s = s + aa
945 work( l ) = work( l ) + aa
946 END DO
947 work( k+j ) = work( k+j ) + s
948 END DO
949* j=k-1 is special :process col A(k-1,0:k-1)
950 s = zero
951 DO i = 0, k - 2
952 aa = abs( a( i+j*lda ) )
953* A(k,i)
954 work( i ) = work( i ) + aa
955 s = s + aa
956 END DO
957* i=k-1
958 aa = abs( real( a( i+j*lda ) ) )
959* A(k-1,k-1)
960 s = s + aa
961 work( i ) = s
962* done with col j=k+1
963 DO j = k, n - 1
964* process col j of A = A(j,0:k-1)
965 s = zero
966 DO i = 0, k - 1
967 aa = abs( a( i+j*lda ) )
968* A(j,i)
969 work( i ) = work( i ) + aa
970 s = s + aa
971 END DO
972 work( j ) = work( j ) + s
973 END DO
974 VALUE = work( 0 )
975 DO i = 1, n-1
976 temp = work( i )
977 IF( VALUE .LT. temp .OR. sisnan( temp ) )
978 $ VALUE = temp
979 END DO
980 END IF
981 ELSE
982* n is even & A is k=n/2 by n+1
983 IF( ilu.EQ.0 ) THEN
984* uplo = 'U'
985 DO i = k, n - 1
986 work( i ) = zero
987 END DO
988 DO j = 0, k - 1
989 s = zero
990 DO i = 0, k - 1
991 aa = abs( a( i+j*lda ) )
992* A(j,i+k)
993 work( i+k ) = work( i+k ) + aa
994 s = s + aa
995 END DO
996 work( j ) = s
997 END DO
998* j=k
999 aa = abs( real( a( 0+j*lda ) ) )
1000* A(k,k)
1001 s = aa
1002 DO i = 1, k - 1
1003 aa = abs( a( i+j*lda ) )
1004* A(k,k+i)
1005 work( i+k ) = work( i+k ) + aa
1006 s = s + aa
1007 END DO
1008 work( j ) = work( j ) + s
1009 DO j = k + 1, n - 1
1010 s = zero
1011 DO i = 0, j - 2 - k
1012 aa = abs( a( i+j*lda ) )
1013* A(i,j-k-1)
1014 work( i ) = work( i ) + aa
1015 s = s + aa
1016 END DO
1017* i=j-1-k
1018 aa = abs( real( a( i+j*lda ) ) )
1019* A(j-k-1,j-k-1)
1020 s = s + aa
1021 work( j-k-1 ) = work( j-k-1 ) + s
1022 i = i + 1
1023 aa = abs( real( a( i+j*lda ) ) )
1024* A(j,j)
1025 s = aa
1026 DO l = j + 1, n - 1
1027 i = i + 1
1028 aa = abs( a( i+j*lda ) )
1029* A(j,l)
1030 work( l ) = work( l ) + aa
1031 s = s + aa
1032 END DO
1033 work( j ) = work( j ) + s
1034 END DO
1035* j=n
1036 s = zero
1037 DO i = 0, k - 2
1038 aa = abs( a( i+j*lda ) )
1039* A(i,k-1)
1040 work( i ) = work( i ) + aa
1041 s = s + aa
1042 END DO
1043* i=k-1
1044 aa = abs( real( a( i+j*lda ) ) )
1045* A(k-1,k-1)
1046 s = s + aa
1047 work( i ) = work( i ) + s
1048 VALUE = work( 0 )
1049 DO i = 1, n-1
1050 temp = work( i )
1051 IF( VALUE .LT. temp .OR. sisnan( temp ) )
1052 $ VALUE = temp
1053 END DO
1054 ELSE
1055* ilu=1 & uplo = 'L'
1056 DO i = k, n - 1
1057 work( i ) = zero
1058 END DO
1059* j=0 is special :process col A(k:n-1,k)
1060 s = abs( real( a( 0 ) ) )
1061* A(k,k)
1062 DO i = 1, k - 1
1063 aa = abs( a( i ) )
1064* A(k+i,k)
1065 work( i+k ) = work( i+k ) + aa
1066 s = s + aa
1067 END DO
1068 work( k ) = work( k ) + s
1069 DO j = 1, k - 1
1070* process
1071 s = zero
1072 DO i = 0, j - 2
1073 aa = abs( a( i+j*lda ) )
1074* A(j-1,i)
1075 work( i ) = work( i ) + aa
1076 s = s + aa
1077 END DO
1078 aa = abs( real( a( i+j*lda ) ) )
1079* i=j-1 so process of A(j-1,j-1)
1080 s = s + aa
1081 work( j-1 ) = s
1082* is initialised here
1083 i = i + 1
1084* i=j process A(j+k,j+k)
1085 aa = abs( real( a( i+j*lda ) ) )
1086 s = aa
1087 DO l = k + j + 1, n - 1
1088 i = i + 1
1089 aa = abs( a( i+j*lda ) )
1090* A(l,k+j)
1091 s = s + aa
1092 work( l ) = work( l ) + aa
1093 END DO
1094 work( k+j ) = work( k+j ) + s
1095 END DO
1096* j=k is special :process col A(k,0:k-1)
1097 s = zero
1098 DO i = 0, k - 2
1099 aa = abs( a( i+j*lda ) )
1100* A(k,i)
1101 work( i ) = work( i ) + aa
1102 s = s + aa
1103 END DO
1104*
1105* i=k-1
1106 aa = abs( real( a( i+j*lda ) ) )
1107* A(k-1,k-1)
1108 s = s + aa
1109 work( i ) = s
1110* done with col j=k+1
1111 DO j = k + 1, n
1112*
1113* process col j-1 of A = A(j-1,0:k-1)
1114 s = zero
1115 DO i = 0, k - 1
1116 aa = abs( a( i+j*lda ) )
1117* A(j-1,i)
1118 work( i ) = work( i ) + aa
1119 s = s + aa
1120 END DO
1121 work( j-1 ) = work( j-1 ) + s
1122 END DO
1123 VALUE = work( 0 )
1124 DO i = 1, n-1
1125 temp = work( i )
1126 IF( VALUE .LT. temp .OR. sisnan( temp ) )
1127 $ VALUE = temp
1128 END DO
1129 END IF
1130 END IF
1131 END IF
1132 ELSE IF( ( lsame( norm, 'F' ) ) .OR. ( lsame( norm, 'E' ) ) ) THEN
1133*
1134* Find normF(A).
1135*
1136 k = ( n+1 ) / 2
1137 scale = zero
1138 s = one
1139 IF( noe.EQ.1 ) THEN
1140* n is odd
1141 IF( ifm.EQ.1 ) THEN
1142* A is normal & A is n by k
1143 IF( ilu.EQ.0 ) THEN
1144* A is upper
1145 DO j = 0, k - 3
1146 CALL classq( k-j-2, a( k+j+1+j*lda ), 1, scale, s )
1147* L at A(k,0)
1148 END DO
1149 DO j = 0, k - 1
1150 CALL classq( k+j-1, a( 0+j*lda ), 1, scale, s )
1151* trap U at A(0,0)
1152 END DO
1153 s = s + s
1154* double s for the off diagonal elements
1155 l = k - 1
1156* -> U(k,k) at A(k-1,0)
1157 DO i = 0, k - 2
1158 aa = real( a( l ) )
1159* U(k+i,k+i)
1160 IF( aa.NE.zero ) THEN
1161 IF( scale.LT.aa ) THEN
1162 s = one + s*( scale / aa )**2
1163 scale = aa
1164 ELSE
1165 s = s + ( aa / scale )**2
1166 END IF
1167 END IF
1168 aa = real( a( l+1 ) )
1169* U(i,i)
1170 IF( aa.NE.zero ) THEN
1171 IF( scale.LT.aa ) THEN
1172 s = one + s*( scale / aa )**2
1173 scale = aa
1174 ELSE
1175 s = s + ( aa / scale )**2
1176 END IF
1177 END IF
1178 l = l + lda + 1
1179 END DO
1180 aa = real( a( l ) )
1181* U(n-1,n-1)
1182 IF( aa.NE.zero ) THEN
1183 IF( scale.LT.aa ) THEN
1184 s = one + s*( scale / aa )**2
1185 scale = aa
1186 ELSE
1187 s = s + ( aa / scale )**2
1188 END IF
1189 END IF
1190 ELSE
1191* ilu=1 & A is lower
1192 DO j = 0, k - 1
1193 CALL classq( n-j-1, a( j+1+j*lda ), 1, scale, s )
1194* trap L at A(0,0)
1195 END DO
1196 DO j = 1, k - 2
1197 CALL classq( j, a( 0+( 1+j )*lda ), 1, scale, s )
1198* U at A(0,1)
1199 END DO
1200 s = s + s
1201* double s for the off diagonal elements
1202 aa = real( a( 0 ) )
1203* L(0,0) at A(0,0)
1204 IF( aa.NE.zero ) THEN
1205 IF( scale.LT.aa ) THEN
1206 s = one + s*( scale / aa )**2
1207 scale = aa
1208 ELSE
1209 s = s + ( aa / scale )**2
1210 END IF
1211 END IF
1212 l = lda
1213* -> L(k,k) at A(0,1)
1214 DO i = 1, k - 1
1215 aa = real( a( l ) )
1216* L(k-1+i,k-1+i)
1217 IF( aa.NE.zero ) THEN
1218 IF( scale.LT.aa ) THEN
1219 s = one + s*( scale / aa )**2
1220 scale = aa
1221 ELSE
1222 s = s + ( aa / scale )**2
1223 END IF
1224 END IF
1225 aa = real( a( l+1 ) )
1226* L(i,i)
1227 IF( aa.NE.zero ) THEN
1228 IF( scale.LT.aa ) THEN
1229 s = one + s*( scale / aa )**2
1230 scale = aa
1231 ELSE
1232 s = s + ( aa / scale )**2
1233 END IF
1234 END IF
1235 l = l + lda + 1
1236 END DO
1237 END IF
1238 ELSE
1239* A is xpose & A is k by n
1240 IF( ilu.EQ.0 ) THEN
1241* A**H is upper
1242 DO j = 1, k - 2
1243 CALL classq( j, a( 0+( k+j )*lda ), 1, scale, s )
1244* U at A(0,k)
1245 END DO
1246 DO j = 0, k - 2
1247 CALL classq( k, a( 0+j*lda ), 1, scale, s )
1248* k by k-1 rect. at A(0,0)
1249 END DO
1250 DO j = 0, k - 2
1251 CALL classq( k-j-1, a( j+1+( j+k-1 )*lda ), 1,
1252 $ scale, s )
1253* L at A(0,k-1)
1254 END DO
1255 s = s + s
1256* double s for the off diagonal elements
1257 l = 0 + k*lda - lda
1258* -> U(k-1,k-1) at A(0,k-1)
1259 aa = real( a( l ) )
1260* U(k-1,k-1)
1261 IF( aa.NE.zero ) THEN
1262 IF( scale.LT.aa ) THEN
1263 s = one + s*( scale / aa )**2
1264 scale = aa
1265 ELSE
1266 s = s + ( aa / scale )**2
1267 END IF
1268 END IF
1269 l = l + lda
1270* -> U(0,0) at A(0,k)
1271 DO j = k, n - 1
1272 aa = real( a( l ) )
1273* -> U(j-k,j-k)
1274 IF( aa.NE.zero ) THEN
1275 IF( scale.LT.aa ) THEN
1276 s = one + s*( scale / aa )**2
1277 scale = aa
1278 ELSE
1279 s = s + ( aa / scale )**2
1280 END IF
1281 END IF
1282 aa = real( a( l+1 ) )
1283* -> U(j,j)
1284 IF( aa.NE.zero ) THEN
1285 IF( scale.LT.aa ) THEN
1286 s = one + s*( scale / aa )**2
1287 scale = aa
1288 ELSE
1289 s = s + ( aa / scale )**2
1290 END IF
1291 END IF
1292 l = l + lda + 1
1293 END DO
1294 ELSE
1295* A**H is lower
1296 DO j = 1, k - 1
1297 CALL classq( j, a( 0+j*lda ), 1, scale, s )
1298* U at A(0,0)
1299 END DO
1300 DO j = k, n - 1
1301 CALL classq( k, a( 0+j*lda ), 1, scale, s )
1302* k by k-1 rect. at A(0,k)
1303 END DO
1304 DO j = 0, k - 3
1305 CALL classq( k-j-2, a( j+2+j*lda ), 1, scale, s )
1306* L at A(1,0)
1307 END DO
1308 s = s + s
1309* double s for the off diagonal elements
1310 l = 0
1311* -> L(0,0) at A(0,0)
1312 DO i = 0, k - 2
1313 aa = real( a( l ) )
1314* L(i,i)
1315 IF( aa.NE.zero ) THEN
1316 IF( scale.LT.aa ) THEN
1317 s = one + s*( scale / aa )**2
1318 scale = aa
1319 ELSE
1320 s = s + ( aa / scale )**2
1321 END IF
1322 END IF
1323 aa = real( a( l+1 ) )
1324* L(k+i,k+i)
1325 IF( aa.NE.zero ) THEN
1326 IF( scale.LT.aa ) THEN
1327 s = one + s*( scale / aa )**2
1328 scale = aa
1329 ELSE
1330 s = s + ( aa / scale )**2
1331 END IF
1332 END IF
1333 l = l + lda + 1
1334 END DO
1335* L-> k-1 + (k-1)*lda or L(k-1,k-1) at A(k-1,k-1)
1336 aa = real( a( l ) )
1337* L(k-1,k-1) at A(k-1,k-1)
1338 IF( aa.NE.zero ) THEN
1339 IF( scale.LT.aa ) THEN
1340 s = one + s*( scale / aa )**2
1341 scale = aa
1342 ELSE
1343 s = s + ( aa / scale )**2
1344 END IF
1345 END IF
1346 END IF
1347 END IF
1348 ELSE
1349* n is even
1350 IF( ifm.EQ.1 ) THEN
1351* A is normal
1352 IF( ilu.EQ.0 ) THEN
1353* A is upper
1354 DO j = 0, k - 2
1355 CALL classq( k-j-1, a( k+j+2+j*lda ), 1, scale, s )
1356* L at A(k+1,0)
1357 END DO
1358 DO j = 0, k - 1
1359 CALL classq( k+j, a( 0+j*lda ), 1, scale, s )
1360* trap U at A(0,0)
1361 END DO
1362 s = s + s
1363* double s for the off diagonal elements
1364 l = k
1365* -> U(k,k) at A(k,0)
1366 DO i = 0, k - 1
1367 aa = real( a( l ) )
1368* U(k+i,k+i)
1369 IF( aa.NE.zero ) THEN
1370 IF( scale.LT.aa ) THEN
1371 s = one + s*( scale / aa )**2
1372 scale = aa
1373 ELSE
1374 s = s + ( aa / scale )**2
1375 END IF
1376 END IF
1377 aa = real( a( l+1 ) )
1378* U(i,i)
1379 IF( aa.NE.zero ) THEN
1380 IF( scale.LT.aa ) THEN
1381 s = one + s*( scale / aa )**2
1382 scale = aa
1383 ELSE
1384 s = s + ( aa / scale )**2
1385 END IF
1386 END IF
1387 l = l + lda + 1
1388 END DO
1389 ELSE
1390* ilu=1 & A is lower
1391 DO j = 0, k - 1
1392 CALL classq( n-j-1, a( j+2+j*lda ), 1, scale, s )
1393* trap L at A(1,0)
1394 END DO
1395 DO j = 1, k - 1
1396 CALL classq( j, a( 0+j*lda ), 1, scale, s )
1397* U at A(0,0)
1398 END DO
1399 s = s + s
1400* double s for the off diagonal elements
1401 l = 0
1402* -> L(k,k) at A(0,0)
1403 DO i = 0, k - 1
1404 aa = real( a( l ) )
1405* L(k-1+i,k-1+i)
1406 IF( aa.NE.zero ) THEN
1407 IF( scale.LT.aa ) THEN
1408 s = one + s*( scale / aa )**2
1409 scale = aa
1410 ELSE
1411 s = s + ( aa / scale )**2
1412 END IF
1413 END IF
1414 aa = real( a( l+1 ) )
1415* L(i,i)
1416 IF( aa.NE.zero ) THEN
1417 IF( scale.LT.aa ) THEN
1418 s = one + s*( scale / aa )**2
1419 scale = aa
1420 ELSE
1421 s = s + ( aa / scale )**2
1422 END IF
1423 END IF
1424 l = l + lda + 1
1425 END DO
1426 END IF
1427 ELSE
1428* A is xpose
1429 IF( ilu.EQ.0 ) THEN
1430* A**H is upper
1431 DO j = 1, k - 1
1432 CALL classq( j, a( 0+( k+1+j )*lda ), 1, scale, s )
1433* U at A(0,k+1)
1434 END DO
1435 DO j = 0, k - 1
1436 CALL classq( k, a( 0+j*lda ), 1, scale, s )
1437* k by k rect. at A(0,0)
1438 END DO
1439 DO j = 0, k - 2
1440 CALL classq( k-j-1, a( j+1+( j+k )*lda ), 1, scale,
1441 $ s )
1442* L at A(0,k)
1443 END DO
1444 s = s + s
1445* double s for the off diagonal elements
1446 l = 0 + k*lda
1447* -> U(k,k) at A(0,k)
1448 aa = real( a( l ) )
1449* U(k,k)
1450 IF( aa.NE.zero ) THEN
1451 IF( scale.LT.aa ) THEN
1452 s = one + s*( scale / aa )**2
1453 scale = aa
1454 ELSE
1455 s = s + ( aa / scale )**2
1456 END IF
1457 END IF
1458 l = l + lda
1459* -> U(0,0) at A(0,k+1)
1460 DO j = k + 1, n - 1
1461 aa = real( a( l ) )
1462* -> U(j-k-1,j-k-1)
1463 IF( aa.NE.zero ) THEN
1464 IF( scale.LT.aa ) THEN
1465 s = one + s*( scale / aa )**2
1466 scale = aa
1467 ELSE
1468 s = s + ( aa / scale )**2
1469 END IF
1470 END IF
1471 aa = real( a( l+1 ) )
1472* -> U(j,j)
1473 IF( aa.NE.zero ) THEN
1474 IF( scale.LT.aa ) THEN
1475 s = one + s*( scale / aa )**2
1476 scale = aa
1477 ELSE
1478 s = s + ( aa / scale )**2
1479 END IF
1480 END IF
1481 l = l + lda + 1
1482 END DO
1483* L=k-1+n*lda
1484* -> U(k-1,k-1) at A(k-1,n)
1485 aa = real( a( l ) )
1486* U(k,k)
1487 IF( aa.NE.zero ) THEN
1488 IF( scale.LT.aa ) THEN
1489 s = one + s*( scale / aa )**2
1490 scale = aa
1491 ELSE
1492 s = s + ( aa / scale )**2
1493 END IF
1494 END IF
1495 ELSE
1496* A**H is lower
1497 DO j = 1, k - 1
1498 CALL classq( j, a( 0+( j+1 )*lda ), 1, scale, s )
1499* U at A(0,1)
1500 END DO
1501 DO j = k + 1, n
1502 CALL classq( k, a( 0+j*lda ), 1, scale, s )
1503* k by k rect. at A(0,k+1)
1504 END DO
1505 DO j = 0, k - 2
1506 CALL classq( k-j-1, a( j+1+j*lda ), 1, scale, s )
1507* L at A(0,0)
1508 END DO
1509 s = s + s
1510* double s for the off diagonal elements
1511 l = 0
1512* -> L(k,k) at A(0,0)
1513 aa = real( a( l ) )
1514* L(k,k) at A(0,0)
1515 IF( aa.NE.zero ) THEN
1516 IF( scale.LT.aa ) THEN
1517 s = one + s*( scale / aa )**2
1518 scale = aa
1519 ELSE
1520 s = s + ( aa / scale )**2
1521 END IF
1522 END IF
1523 l = lda
1524* -> L(0,0) at A(0,1)
1525 DO i = 0, k - 2
1526 aa = real( a( l ) )
1527* L(i,i)
1528 IF( aa.NE.zero ) THEN
1529 IF( scale.LT.aa ) THEN
1530 s = one + s*( scale / aa )**2
1531 scale = aa
1532 ELSE
1533 s = s + ( aa / scale )**2
1534 END IF
1535 END IF
1536 aa = real( a( l+1 ) )
1537* L(k+i+1,k+i+1)
1538 IF( aa.NE.zero ) THEN
1539 IF( scale.LT.aa ) THEN
1540 s = one + s*( scale / aa )**2
1541 scale = aa
1542 ELSE
1543 s = s + ( aa / scale )**2
1544 END IF
1545 END IF
1546 l = l + lda + 1
1547 END DO
1548* L-> k - 1 + k*lda or L(k-1,k-1) at A(k-1,k)
1549 aa = real( a( l ) )
1550* L(k-1,k-1) at A(k-1,k)
1551 IF( aa.NE.zero ) THEN
1552 IF( scale.LT.aa ) THEN
1553 s = one + s*( scale / aa )**2
1554 scale = aa
1555 ELSE
1556 s = s + ( aa / scale )**2
1557 END IF
1558 END IF
1559 END IF
1560 END IF
1561 END IF
1562 VALUE = scale*sqrt( s )
1563 END IF
1564*
1565 clanhf = VALUE
1566 RETURN
1567*
1568* End of CLANHF
1569*
sisnan
logical function sisnan(sin)
SISNAN tests input for NaN.
Definition sisnan.f:59
clanhf
real function clanhf(norm, transr, uplo, n, a, work)
CLANHF returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition clanhf.f:246
classq
subroutine classq(n, x, incx, scale, sumsq)
CLASSQ updates a sum of squares represented in scaled form.
Definition classq.f90:124
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the call graph for this function:
Here is the caller graph for this function: