Home
Browse Files
   1       DOUBLE PRECISION FUNCTION ZLANHF( NORM, TRANSR, UPLO, N, A, WORK )
   2 *
   3 *  -- LAPACK routine (version 3.3.1)                                    --
   4 *
   5 *  -- Contributed by Fred Gustavson of the IBM Watson Research Center --
   6 *  -- April 2009                                                      --
   7 *
   8 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
   9 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  10 *
  11 *     .. Scalar Arguments ..
  12       CHARACTER          NORM, TRANSR, UPLO
  13       INTEGER            N
  14 *     ..
  15 *     .. Array Arguments ..
  16       DOUBLE PRECISION   WORK( 0: * )
  17       COMPLEX*16         A( 0: * )
  18 *     ..
  19 *
  20 *  Purpose
  21 *  =======
  22 *
  23 *  ZLANHF  returns the value of the one norm,  or the Frobenius norm, or
  24 *  the  infinity norm,  or the  element of  largest absolute value  of a
  25 *  complex Hermitian matrix A in RFP format.
  26 *
  27 *  Description
  28 *  ===========
  29 *
  30 *  ZLANHF returns the value
  31 *
  32 *     ZLANHF = ( max(abs(A(i,j))), NORM = 'M' or 'm'
  33 *              (
  34 *              ( norm1(A),         NORM = '1', 'O' or 'o'
  35 *              (
  36 *              ( normI(A),         NORM = 'I' or 'i'
  37 *              (
  38 *              ( normF(A),         NORM = 'F', 'f', 'E' or 'e'
  39 *
  40 *  where  norm1  denotes the  one norm of a matrix (maximum column sum),
  41 *  normI  denotes the  infinity norm  of a matrix  (maximum row sum) and
  42 *  normF  denotes the  Frobenius norm of a matrix (square root of sum of
  43 *  squares).  Note that  max(abs(A(i,j)))  is not a  matrix norm.
  44 *
  45 *  Arguments
  46 *  =========
  47 *
  48 *  NORM      (input) CHARACTER
  49 *            Specifies the value to be returned in ZLANHF as described
  50 *            above.
  51 *
  52 *  TRANSR    (input) CHARACTER
  53 *            Specifies whether the RFP format of A is normal or
  54 *            conjugate-transposed format.
  55 *            = 'N':  RFP format is Normal
  56 *            = 'C':  RFP format is Conjugate-transposed
  57 *
  58 *  UPLO      (input) CHARACTER
  59 *            On entry, UPLO specifies whether the RFP matrix A came from
  60 *            an upper or lower triangular matrix as follows:
  61 *
  62 *            UPLO = 'U' or 'u' RFP A came from an upper triangular
  63 *            matrix
  64 *
  65 *            UPLO = 'L' or 'l' RFP A came from a  lower triangular
  66 *            matrix
  67 *
  68 *  N         (input) INTEGER
  69 *            The order of the matrix A.  N >= 0.  When N = 0, ZLANHF is
  70 *            set to zero.
  71 *
  72 *   A        (input) COMPLEX*16 array, dimension ( N*(N+1)/2 );
  73 *            On entry, the matrix A in RFP Format.
  74 *            RFP Format is described by TRANSR, UPLO and N as follows:
  75 *            If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
  76 *            K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
  77 *            TRANSR = 'C' then RFP is the Conjugate-transpose of RFP A
  78 *            as defined when TRANSR = 'N'. The contents of RFP A are
  79 *            defined by UPLO as follows: If UPLO = 'U' the RFP A
  80 *            contains the ( N*(N+1)/2 ) elements of upper packed A
  81 *            either in normal or conjugate-transpose Format. If
  82 *            UPLO = 'L' the RFP A contains the ( N*(N+1) /2 ) elements
  83 *            of lower packed A either in normal or conjugate-transpose
  84 *            Format. The LDA of RFP A is (N+1)/2 when TRANSR = 'C'. When
  85 *            TRANSR is 'N' the LDA is N+1 when N is even and is N when
  86 *            is odd. See the Note below for more details.
  87 *            Unchanged on exit.
  88 *
  89 *  WORK      (workspace) DOUBLE PRECISION array, dimension (LWORK),
  90 *            where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
  91 *            WORK is not referenced.
  92 *
  93 *  Further Details
  94 *  ===============
  95 *
  96 *  We first consider Standard Packed Format when N is even.
  97 *  We give an example where N = 6.
  98 *
  99 *      AP is Upper             AP is Lower
 100 *
 101 *   00 01 02 03 04 05       00
 102 *      11 12 13 14 15       10 11
 103 *         22 23 24 25       20 21 22
 104 *            33 34 35       30 31 32 33
 105 *               44 45       40 41 42 43 44
 106 *                  55       50 51 52 53 54 55
 107 *
 108 *
 109 *  Let TRANSR = 'N'. RFP holds AP as follows:
 110 *  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
 111 *  three columns of AP upper. The lower triangle A(4:6,0:2) consists of
 112 *  conjugate-transpose of the first three columns of AP upper.
 113 *  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
 114 *  three columns of AP lower. The upper triangle A(0:2,0:2) consists of
 115 *  conjugate-transpose of the last three columns of AP lower.
 116 *  To denote conjugate we place -- above the element. This covers the
 117 *  case N even and TRANSR = 'N'.
 118 *
 119 *         RFP A                   RFP A
 120 *
 121 *                                -- -- --
 122 *        03 04 05                33 43 53
 123 *                                   -- --
 124 *        13 14 15                00 44 54
 125 *                                      --
 126 *        23 24 25                10 11 55
 127 *
 128 *        33 34 35                20 21 22
 129 *        --
 130 *        00 44 45                30 31 32
 131 *        -- --
 132 *        01 11 55                40 41 42
 133 *        -- -- --
 134 *        02 12 22                50 51 52
 135 *
 136 *  Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
 137 *  transpose of RFP A above. One therefore gets:
 138 *
 139 *
 140 *           RFP A                   RFP A
 141 *
 142 *     -- -- -- --                -- -- -- -- -- --
 143 *     03 13 23 33 00 01 02    33 00 10 20 30 40 50
 144 *     -- -- -- -- --                -- -- -- -- --
 145 *     04 14 24 34 44 11 12    43 44 11 21 31 41 51
 146 *     -- -- -- -- -- --                -- -- -- --
 147 *     05 15 25 35 45 55 22    53 54 55 22 32 42 52
 148 *
 149 *
 150 *  We next  consider Standard Packed Format when N is odd.
 151 *  We give an example where N = 5.
 152 *
 153 *     AP is Upper                 AP is Lower
 154 *
 155 *   00 01 02 03 04              00
 156 *      11 12 13 14              10 11
 157 *         22 23 24              20 21 22
 158 *            33 34              30 31 32 33
 159 *               44              40 41 42 43 44
 160 *
 161 *
 162 *  Let TRANSR = 'N'. RFP holds AP as follows:
 163 *  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
 164 *  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
 165 *  conjugate-transpose of the first two   columns of AP upper.
 166 *  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
 167 *  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
 168 *  conjugate-transpose of the last two   columns of AP lower.
 169 *  To denote conjugate we place -- above the element. This covers the
 170 *  case N odd  and TRANSR = 'N'.
 171 *
 172 *         RFP A                   RFP A
 173 *
 174 *                                   -- --
 175 *        02 03 04                00 33 43
 176 *                                      --
 177 *        12 13 14                10 11 44
 178 *
 179 *        22 23 24                20 21 22
 180 *        --
 181 *        00 33 34                30 31 32
 182 *        -- --
 183 *        01 11 44                40 41 42
 184 *
 185 *  Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
 186 *  transpose of RFP A above. One therefore gets:
 187 *
 188 *
 189 *           RFP A                   RFP A
 190 *
 191 *     -- -- --                   -- -- -- -- -- --
 192 *     02 12 22 00 01             00 10 20 30 40 50
 193 *     -- -- -- --                   -- -- -- -- --
 194 *     03 13 23 33 11             33 11 21 31 41 51
 195 *     -- -- -- -- --                   -- -- -- --
 196 *     04 14 24 34 44             43 44 22 32 42 52
 197 *
 198 *  =====================================================================
 199 *
 200 *     .. Parameters ..
 201       DOUBLE PRECISION   ONE, ZERO
 202       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
 203 *     ..
 204 *     .. Local Scalars ..
 205       INTEGER            I, J, IFM, ILU, NOE, N1, K, L, LDA
 206       DOUBLE PRECISION   SCALE, S, VALUE, AA
 207 *     ..
 208 *     .. External Functions ..
 209       LOGICAL            LSAME
 210       INTEGER            IDAMAX
 211       EXTERNAL           LSAME, IDAMAX
 212 *     ..
 213 *     .. External Subroutines ..
 214       EXTERNAL           ZLASSQ
 215 *     ..
 216 *     .. Intrinsic Functions ..
 217       INTRINSIC          ABS, DBLE, MAX, SQRT
 218 *     ..
 219 *     .. Executable Statements ..
 220 *
 221       IF( N.EQ.0 ) THEN
 222          ZLANHF = ZERO
 223          RETURN
 224       END IF
 225 *
 226 *     set noe = 1 if n is odd. if n is even set noe=0
 227 *
 228       NOE = 1
 229       IF( MOD( N, 2 ).EQ.0 )
 230      $   NOE = 0
 231 *
 232 *     set ifm = 0 when form='C' or 'c' and 1 otherwise
 233 *
 234       IFM = 1
 235       IF( LSAME( TRANSR, 'C' ) )
 236      $   IFM = 0
 237 *
 238 *     set ilu = 0 when uplo='U or 'u' and 1 otherwise
 239 *
 240       ILU = 1
 241       IF( LSAME( UPLO, 'U' ) )
 242      $   ILU = 0
 243 *
 244 *     set lda = (n+1)/2 when ifm = 0
 245 *     set lda = n when ifm = 1 and noe = 1
 246 *     set lda = n+1 when ifm = 1 and noe = 0
 247 *
 248       IF( IFM.EQ.1 ) THEN
 249          IF( NOE.EQ.1 ) THEN
 250             LDA = N
 251          ELSE
 252 *           noe=0
 253             LDA = N + 1
 254          END IF
 255       ELSE
 256 *        ifm=0
 257          LDA = ( N+1 ) / 2
 258       END IF
 259 *
 260       IF( LSAME( NORM, 'M' ) ) THEN
 261 *
 262 *       Find max(abs(A(i,j))).
 263 *
 264          K = ( N+1 ) / 2
 265          VALUE = ZERO
 266          IF( NOE.EQ.1 ) THEN
 267 *           n is odd & n = k + k - 1
 268             IF( IFM.EQ.1 ) THEN
 269 *              A is n by k
 270                IF( ILU.EQ.1 ) THEN
 271 *                 uplo ='L'
 272                   J = 0
 273 *                 -> L(0,0)
 274                   VALUE = MAX( VALUE, ABS( DBLE( A( J+J*LDA ) ) ) )
 275                   DO I = 1, N - 1
 276                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 277                   END DO
 278                   DO J = 1, K - 1
 279                      DO I = 0, J - 2
 280                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 281                      END DO
 282                      I = J - 1
 283 *                    L(k+j,k+j)
 284                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 285                      I = J
 286 *                    -> L(j,j)
 287                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 288                      DO I = J + 1, N - 1
 289                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 290                      END DO
 291                   END DO
 292                ELSE
 293 *                 uplo = 'U'
 294                   DO J = 0, K - 2
 295                      DO I = 0, K + J - 2
 296                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 297                      END DO
 298                      I = K + J - 1
 299 *                    -> U(i,i)
 300                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 301                      I = I + 1
 302 *                    =k+j; i -> U(j,j)
 303                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 304                      DO I = K + J + 1, N - 1
 305                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 306                      END DO
 307                   END DO
 308                   DO I = 0, N - 2
 309                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 310 *                    j=k-1
 311                   END DO
 312 *                 i=n-1 -> U(n-1,n-1)
 313                   VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 314                END IF
 315             ELSE
 316 *              xpose case; A is k by n
 317                IF( ILU.EQ.1 ) THEN
 318 *                 uplo ='L'
 319                   DO J = 0, K - 2
 320                      DO I = 0, J - 1
 321                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 322                      END DO
 323                      I = J
 324 *                    L(i,i)
 325                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 326                      I = J + 1
 327 *                    L(j+k,j+k)
 328                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 329                      DO I = J + 2, K - 1
 330                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 331                      END DO
 332                   END DO
 333                   J = K - 1
 334                   DO I = 0, K - 2
 335                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 336                   END DO
 337                   I = K - 1
 338 *                 -> L(i,i) is at A(i,j)
 339                   VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 340                   DO J = K, N - 1
 341                      DO I = 0, K - 1
 342                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 343                      END DO
 344                   END DO
 345                ELSE
 346 *                 uplo = 'U'
 347                   DO J = 0, K - 2
 348                      DO I = 0, K - 1
 349                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 350                      END DO
 351                   END DO
 352                   J = K - 1
 353 *                 -> U(j,j) is at A(0,j)
 354                   VALUE = MAX( VALUE, ABS( DBLE( A( 0+J*LDA ) ) ) )
 355                   DO I = 1, K - 1
 356                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 357                   END DO
 358                   DO J = K, N - 1
 359                      DO I = 0, J - K - 1
 360                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 361                      END DO
 362                      I = J - K
 363 *                    -> U(i,i) at A(i,j)
 364                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 365                      I = J - K + 1
 366 *                    U(j,j)
 367                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 368                      DO I = J - K + 2, K - 1
 369                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 370                      END DO
 371                   END DO
 372                END IF
 373             END IF
 374          ELSE
 375 *           n is even & k = n/2
 376             IF( IFM.EQ.1 ) THEN
 377 *              A is n+1 by k
 378                IF( ILU.EQ.1 ) THEN
 379 *                 uplo ='L'
 380                   J = 0
 381 *                 -> L(k,k) & j=1 -> L(0,0)
 382                   VALUE = MAX( VALUE, ABS( DBLE( A( J+J*LDA ) ) ) )
 383                   VALUE = MAX( VALUE, ABS( DBLE( A( J+1+J*LDA ) ) ) )
 384                   DO I = 2, N
 385                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 386                   END DO
 387                   DO J = 1, K - 1
 388                      DO I = 0, J - 1
 389                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 390                      END DO
 391                      I = J
 392 *                    L(k+j,k+j)
 393                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 394                      I = J + 1
 395 *                    -> L(j,j)
 396                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 397                      DO I = J + 2, N
 398                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 399                      END DO
 400                   END DO
 401                ELSE
 402 *                 uplo = 'U'
 403                   DO J = 0, K - 2
 404                      DO I = 0, K + J - 1
 405                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 406                      END DO
 407                      I = K + J
 408 *                    -> U(i,i)
 409                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 410                      I = I + 1
 411 *                    =k+j+1; i -> U(j,j)
 412                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 413                      DO I = K + J + 2, N
 414                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 415                      END DO
 416                   END DO
 417                   DO I = 0, N - 2
 418                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 419 *                    j=k-1
 420                   END DO
 421 *                 i=n-1 -> U(n-1,n-1)
 422                   VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 423                   I = N
 424 *                 -> U(k-1,k-1)
 425                   VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 426                END IF
 427             ELSE
 428 *              xpose case; A is k by n+1
 429                IF( ILU.EQ.1 ) THEN
 430 *                 uplo ='L'
 431                   J = 0
 432 *                 -> L(k,k) at A(0,0)
 433                   VALUE = MAX( VALUE, ABS( DBLE( A( J+J*LDA ) ) ) )
 434                   DO I = 1, K - 1
 435                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 436                   END DO
 437                   DO J = 1, K - 1
 438                      DO I = 0, J - 2
 439                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 440                      END DO
 441                      I = J - 1
 442 *                    L(i,i)
 443                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 444                      I = J
 445 *                    L(j+k,j+k)
 446                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 447                      DO I = J + 1, K - 1
 448                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 449                      END DO
 450                   END DO
 451                   J = K
 452                   DO I = 0, K - 2
 453                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 454                   END DO
 455                   I = K - 1
 456 *                 -> L(i,i) is at A(i,j)
 457                   VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 458                   DO J = K + 1, N
 459                      DO I = 0, K - 1
 460                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 461                      END DO
 462                   END DO
 463                ELSE
 464 *                 uplo = 'U'
 465                   DO J = 0, K - 1
 466                      DO I = 0, K - 1
 467                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 468                      END DO
 469                   END DO
 470                   J = K
 471 *                 -> U(j,j) is at A(0,j)
 472                   VALUE = MAX( VALUE, ABS( DBLE( A( 0+J*LDA ) ) ) )
 473                   DO I = 1, K - 1
 474                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 475                   END DO
 476                   DO J = K + 1, N - 1
 477                      DO I = 0, J - K - 2
 478                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 479                      END DO
 480                      I = J - K - 1
 481 *                    -> U(i,i) at A(i,j)
 482                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 483                      I = J - K
 484 *                    U(j,j)
 485                      VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 486                      DO I = J - K + 1, K - 1
 487                         VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 488                      END DO
 489                   END DO
 490                   J = N
 491                   DO I = 0, K - 2
 492                      VALUE = MAX( VALUE, ABS( A( I+J*LDA ) ) )
 493                   END DO
 494                   I = K - 1
 495 *                 U(k,k) at A(i,j)
 496                   VALUE = MAX( VALUE, ABS( DBLE( A( I+J*LDA ) ) ) )
 497                END IF
 498             END IF
 499          END IF
 500       ELSE IF( ( LSAME( NORM, 'I' ) ) .OR. ( LSAME( NORM, 'O' ) ) .OR.
 501      $         ( NORM.EQ.'1' ) ) THEN
 502 *
 503 *       Find normI(A) ( = norm1(A), since A is Hermitian).
 504 *
 505          IF( IFM.EQ.1 ) THEN
 506 *           A is 'N'
 507             K = N / 2
 508             IF( NOE.EQ.1 ) THEN
 509 *              n is odd & A is n by (n+1)/2
 510                IF( ILU.EQ.0 ) THEN
 511 *                 uplo = 'U'
 512                   DO I = 0, K - 1
 513                      WORK( I ) = ZERO
 514                   END DO
 515                   DO J = 0, K
 516                      S = ZERO
 517                      DO I = 0, K + J - 1
 518                         AA = ABS( A( I+J*LDA ) )
 519 *                       -> A(i,j+k)
 520                         S = S + AA
 521                         WORK( I ) = WORK( I ) + AA
 522                      END DO
 523                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 524 *                    -> A(j+k,j+k)
 525                      WORK( J+K ) = S + AA
 526                      IF( I.EQ.K+K )
 527      $                  GO TO 10
 528                      I = I + 1
 529                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 530 *                    -> A(j,j)
 531                      WORK( J ) = WORK( J ) + AA
 532                      S = ZERO
 533                      DO L = J + 1, K - 1
 534                         I = I + 1
 535                         AA = ABS( A( I+J*LDA ) )
 536 *                       -> A(l,j)
 537                         S = S + AA
 538                         WORK( L ) = WORK( L ) + AA
 539                      END DO
 540                      WORK( J ) = WORK( J ) + S
 541                   END DO
 542    10             CONTINUE
 543                   I = IDAMAX( N, WORK, 1 )
 544                   VALUE = WORK( I-1 )
 545                ELSE
 546 *                 ilu = 1 & uplo = 'L'
 547                   K = K + 1
 548 *                 k=(n+1)/2 for n odd and ilu=1
 549                   DO I = K, N - 1
 550                      WORK( I ) = ZERO
 551                   END DO
 552                   DO J = K - 1, 0, -1
 553                      S = ZERO
 554                      DO I = 0, J - 2
 555                         AA = ABS( A( I+J*LDA ) )
 556 *                       -> A(j+k,i+k)
 557                         S = S + AA
 558                         WORK( I+K ) = WORK( I+K ) + AA
 559                      END DO
 560                      IF( J.GT.0 ) THEN
 561                         AA = ABS( DBLE( A( I+J*LDA ) ) )
 562 *                       -> A(j+k,j+k)
 563                         S = S + AA
 564                         WORK( I+K ) = WORK( I+K ) + S
 565 *                       i=j
 566                         I = I + 1
 567                      END IF
 568                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 569 *                    -> A(j,j)
 570                      WORK( J ) = AA
 571                      S = ZERO
 572                      DO L = J + 1, N - 1
 573                         I = I + 1
 574                         AA = ABS( A( I+J*LDA ) )
 575 *                       -> A(l,j)
 576                         S = S + AA
 577                         WORK( L ) = WORK( L ) + AA
 578                      END DO
 579                      WORK( J ) = WORK( J ) + S
 580                   END DO
 581                   I = IDAMAX( N, WORK, 1 )
 582                   VALUE = WORK( I-1 )
 583                END IF
 584             ELSE
 585 *              n is even & A is n+1 by k = n/2
 586                IF( ILU.EQ.0 ) THEN
 587 *                 uplo = 'U'
 588                   DO I = 0, K - 1
 589                      WORK( I ) = ZERO
 590                   END DO
 591                   DO J = 0, K - 1
 592                      S = ZERO
 593                      DO I = 0, K + J - 1
 594                         AA = ABS( A( I+J*LDA ) )
 595 *                       -> A(i,j+k)
 596                         S = S + AA
 597                         WORK( I ) = WORK( I ) + AA
 598                      END DO
 599                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 600 *                    -> A(j+k,j+k)
 601                      WORK( J+K ) = S + AA
 602                      I = I + 1
 603                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 604 *                    -> A(j,j)
 605                      WORK( J ) = WORK( J ) + AA
 606                      S = ZERO
 607                      DO L = J + 1, K - 1
 608                         I = I + 1
 609                         AA = ABS( A( I+J*LDA ) )
 610 *                       -> A(l,j)
 611                         S = S + AA
 612                         WORK( L ) = WORK( L ) + AA
 613                      END DO
 614                      WORK( J ) = WORK( J ) + S
 615                   END DO
 616                   I = IDAMAX( N, WORK, 1 )
 617                   VALUE = WORK( I-1 )
 618                ELSE
 619 *                 ilu = 1 & uplo = 'L'
 620                   DO I = K, N - 1
 621                      WORK( I ) = ZERO
 622                   END DO
 623                   DO J = K - 1, 0, -1
 624                      S = ZERO
 625                      DO I = 0, J - 1
 626                         AA = ABS( A( I+J*LDA ) )
 627 *                       -> A(j+k,i+k)
 628                         S = S + AA
 629                         WORK( I+K ) = WORK( I+K ) + AA
 630                      END DO
 631                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 632 *                    -> A(j+k,j+k)
 633                      S = S + AA
 634                      WORK( I+K ) = WORK( I+K ) + S
 635 *                    i=j
 636                      I = I + 1
 637                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 638 *                    -> A(j,j)
 639                      WORK( J ) = AA
 640                      S = ZERO
 641                      DO L = J + 1, N - 1
 642                         I = I + 1
 643                         AA = ABS( A( I+J*LDA ) )
 644 *                       -> A(l,j)
 645                         S = S + AA
 646                         WORK( L ) = WORK( L ) + AA
 647                      END DO
 648                      WORK( J ) = WORK( J ) + S
 649                   END DO
 650                   I = IDAMAX( N, WORK, 1 )
 651                   VALUE = WORK( I-1 )
 652                END IF
 653             END IF
 654          ELSE
 655 *           ifm=0
 656             K = N / 2
 657             IF( NOE.EQ.1 ) THEN
 658 *              n is odd & A is (n+1)/2 by n
 659                IF( ILU.EQ.0 ) THEN
 660 *                 uplo = 'U'
 661                   N1 = K
 662 *                 n/2
 663                   K = K + 1
 664 *                 k is the row size and lda
 665                   DO I = N1, N - 1
 666                      WORK( I ) = ZERO
 667                   END DO
 668                   DO J = 0, N1 - 1
 669                      S = ZERO
 670                      DO I = 0, K - 1
 671                         AA = ABS( A( I+J*LDA ) )
 672 *                       A(j,n1+i)
 673                         WORK( I+N1 ) = WORK( I+N1 ) + AA
 674                         S = S + AA
 675                      END DO
 676                      WORK( J ) = S
 677                   END DO
 678 *                 j=n1=k-1 is special
 679                   S = ABS( DBLE( A( 0+J*LDA ) ) )
 680 *                 A(k-1,k-1)
 681                   DO I = 1, K - 1
 682                      AA = ABS( A( I+J*LDA ) )
 683 *                    A(k-1,i+n1)
 684                      WORK( I+N1 ) = WORK( I+N1 ) + AA
 685                      S = S + AA
 686                   END DO
 687                   WORK( J ) = WORK( J ) + S
 688                   DO J = K, N - 1
 689                      S = ZERO
 690                      DO I = 0, J - K - 1
 691                         AA = ABS( A( I+J*LDA ) )
 692 *                       A(i,j-k)
 693                         WORK( I ) = WORK( I ) + AA
 694                         S = S + AA
 695                      END DO
 696 *                    i=j-k
 697                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 698 *                    A(j-k,j-k)
 699                      S = S + AA
 700                      WORK( J-K ) = WORK( J-K ) + S
 701                      I = I + 1
 702                      S = ABS( DBLE( A( I+J*LDA ) ) )
 703 *                    A(j,j)
 704                      DO L = J + 1, N - 1
 705                         I = I + 1
 706                         AA = ABS( A( I+J*LDA ) )
 707 *                       A(j,l)
 708                         WORK( L ) = WORK( L ) + AA
 709                         S = S + AA
 710                      END DO
 711                      WORK( J ) = WORK( J ) + S
 712                   END DO
 713                   I = IDAMAX( N, WORK, 1 )
 714                   VALUE = WORK( I-1 )
 715                ELSE
 716 *                 ilu=1 & uplo = 'L'
 717                   K = K + 1
 718 *                 k=(n+1)/2 for n odd and ilu=1
 719                   DO I = K, N - 1
 720                      WORK( I ) = ZERO
 721                   END DO
 722                   DO J = 0, K - 2
 723 *                    process
 724                      S = ZERO
 725                      DO I = 0, J - 1
 726                         AA = ABS( A( I+J*LDA ) )
 727 *                       A(j,i)
 728                         WORK( I ) = WORK( I ) + AA
 729                         S = S + AA
 730                      END DO
 731                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 732 *                    i=j so process of A(j,j)
 733                      S = S + AA
 734                      WORK( J ) = S
 735 *                    is initialised here
 736                      I = I + 1
 737 *                    i=j process A(j+k,j+k)
 738                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 739                      S = AA
 740                      DO L = K + J + 1, N - 1
 741                         I = I + 1
 742                         AA = ABS( A( I+J*LDA ) )
 743 *                       A(l,k+j)
 744                         S = S + AA
 745                         WORK( L ) = WORK( L ) + AA
 746                      END DO
 747                      WORK( K+J ) = WORK( K+J ) + S
 748                   END DO
 749 *                 j=k-1 is special :process col A(k-1,0:k-1)
 750                   S = ZERO
 751                   DO I = 0, K - 2
 752                      AA = ABS( A( I+J*LDA ) )
 753 *                    A(k,i)
 754                      WORK( I ) = WORK( I ) + AA
 755                      S = S + AA
 756                   END DO
 757 *                 i=k-1
 758                   AA = ABS( DBLE( A( I+J*LDA ) ) )
 759 *                 A(k-1,k-1)
 760                   S = S + AA
 761                   WORK( I ) = S
 762 *                 done with col j=k+1
 763                   DO J = K, N - 1
 764 *                    process col j of A = A(j,0:k-1)
 765                      S = ZERO
 766                      DO I = 0, K - 1
 767                         AA = ABS( A( I+J*LDA ) )
 768 *                       A(j,i)
 769                         WORK( I ) = WORK( I ) + AA
 770                         S = S + AA
 771                      END DO
 772                      WORK( J ) = WORK( J ) + S
 773                   END DO
 774                   I = IDAMAX( N, WORK, 1 )
 775                   VALUE = WORK( I-1 )
 776                END IF
 777             ELSE
 778 *              n is even & A is k=n/2 by n+1
 779                IF( ILU.EQ.0 ) THEN
 780 *                 uplo = 'U'
 781                   DO I = K, N - 1
 782                      WORK( I ) = ZERO
 783                   END DO
 784                   DO J = 0, K - 1
 785                      S = ZERO
 786                      DO I = 0, K - 1
 787                         AA = ABS( A( I+J*LDA ) )
 788 *                       A(j,i+k)
 789                         WORK( I+K ) = WORK( I+K ) + AA
 790                         S = S + AA
 791                      END DO
 792                      WORK( J ) = S
 793                   END DO
 794 *                 j=k
 795                   AA = ABS( DBLE( A( 0+J*LDA ) ) )
 796 *                 A(k,k)
 797                   S = AA
 798                   DO I = 1, K - 1
 799                      AA = ABS( A( I+J*LDA ) )
 800 *                    A(k,k+i)
 801                      WORK( I+K ) = WORK( I+K ) + AA
 802                      S = S + AA
 803                   END DO
 804                   WORK( J ) = WORK( J ) + S
 805                   DO J = K + 1, N - 1
 806                      S = ZERO
 807                      DO I = 0, J - 2 - K
 808                         AA = ABS( A( I+J*LDA ) )
 809 *                       A(i,j-k-1)
 810                         WORK( I ) = WORK( I ) + AA
 811                         S = S + AA
 812                      END DO
 813 *                    i=j-1-k
 814                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 815 *                    A(j-k-1,j-k-1)
 816                      S = S + AA
 817                      WORK( J-K-1 ) = WORK( J-K-1 ) + S
 818                      I = I + 1
 819                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 820 *                    A(j,j)
 821                      S = AA
 822                      DO L = J + 1, N - 1
 823                         I = I + 1
 824                         AA = ABS( A( I+J*LDA ) )
 825 *                       A(j,l)
 826                         WORK( L ) = WORK( L ) + AA
 827                         S = S + AA
 828                      END DO
 829                      WORK( J ) = WORK( J ) + S
 830                   END DO
 831 *                 j=n
 832                   S = ZERO
 833                   DO I = 0, K - 2
 834                      AA = ABS( A( I+J*LDA ) )
 835 *                    A(i,k-1)
 836                      WORK( I ) = WORK( I ) + AA
 837                      S = S + AA
 838                   END DO
 839 *                 i=k-1
 840                   AA = ABS( DBLE( A( I+J*LDA ) ) )
 841 *                 A(k-1,k-1)
 842                   S = S + AA
 843                   WORK( I ) = WORK( I ) + S
 844                   I = IDAMAX( N, WORK, 1 )
 845                   VALUE = WORK( I-1 )
 846                ELSE
 847 *                 ilu=1 & uplo = 'L'
 848                   DO I = K, N - 1
 849                      WORK( I ) = ZERO
 850                   END DO
 851 *                 j=0 is special :process col A(k:n-1,k)
 852                   S = ABS( DBLE( A( 0 ) ) )
 853 *                 A(k,k)
 854                   DO I = 1, K - 1
 855                      AA = ABS( A( I ) )
 856 *                    A(k+i,k)
 857                      WORK( I+K ) = WORK( I+K ) + AA
 858                      S = S + AA
 859                   END DO
 860                   WORK( K ) = WORK( K ) + S
 861                   DO J = 1, K - 1
 862 *                    process
 863                      S = ZERO
 864                      DO I = 0, J - 2
 865                         AA = ABS( A( I+J*LDA ) )
 866 *                       A(j-1,i)
 867                         WORK( I ) = WORK( I ) + AA
 868                         S = S + AA
 869                      END DO
 870                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 871 *                    i=j-1 so process of A(j-1,j-1)
 872                      S = S + AA
 873                      WORK( J-1 ) = S
 874 *                    is initialised here
 875                      I = I + 1
 876 *                    i=j process A(j+k,j+k)
 877                      AA = ABS( DBLE( A( I+J*LDA ) ) )
 878                      S = AA
 879                      DO L = K + J + 1, N - 1
 880                         I = I + 1
 881                         AA = ABS( A( I+J*LDA ) )
 882 *                       A(l,k+j)
 883                         S = S + AA
 884                         WORK( L ) = WORK( L ) + AA
 885                      END DO
 886                      WORK( K+J ) = WORK( K+J ) + S
 887                   END DO
 888 *                 j=k is special :process col A(k,0:k-1)
 889                   S = ZERO
 890                   DO I = 0, K - 2
 891                      AA = ABS( A( I+J*LDA ) )
 892 *                    A(k,i)
 893                      WORK( I ) = WORK( I ) + AA
 894                      S = S + AA
 895                   END DO
 896 *
 897 *                 i=k-1
 898                   AA = ABS( DBLE( A( I+J*LDA ) ) )
 899 *                 A(k-1,k-1)
 900                   S = S + AA
 901                   WORK( I ) = S
 902 *                 done with col j=k+1
 903                   DO J = K + 1, N
 904 *
 905 *                    process col j-1 of A = A(j-1,0:k-1)
 906                      S = ZERO
 907                      DO I = 0, K - 1
 908                         AA = ABS( A( I+J*LDA ) )
 909 *                       A(j-1,i)
 910                         WORK( I ) = WORK( I ) + AA
 911                         S = S + AA
 912                      END DO
 913                      WORK( J-1 ) = WORK( J-1 ) + S
 914                   END DO
 915                   I = IDAMAX( N, WORK, 1 )
 916                   VALUE = WORK( I-1 )
 917                END IF
 918             END IF
 919          END IF
 920       ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
 921 *
 922 *       Find normF(A).
 923 *
 924          K = ( N+1 ) / 2
 925          SCALE = ZERO
 926          S = ONE
 927          IF( NOE.EQ.1 ) THEN
 928 *           n is odd
 929             IF( IFM.EQ.1 ) THEN
 930 *              A is normal & A is n by k
 931                IF( ILU.EQ.0 ) THEN
 932 *                 A is upper
 933                   DO J = 0, K - 3
 934                      CALL ZLASSQ( K-J-2, A( K+J+1+J*LDA ), 1, SCALE, S )
 935 *                    L at A(k,0)
 936                   END DO
 937                   DO J = 0, K - 1
 938                      CALL ZLASSQ( K+J-1, A( 0+J*LDA ), 1, SCALE, S )
 939 *                    trap U at A(0,0)
 940                   END DO
 941                   S = S + S
 942 *                 double s for the off diagonal elements
 943                   L = K - 1
 944 *                 -> U(k,k) at A(k-1,0)
 945                   DO I = 0, K - 2
 946                      AA = DBLE( A( L ) )
 947 *                    U(k+i,k+i)
 948                      IF( AA.NE.ZERO ) THEN
 949                         IF( SCALE.LT.AA ) THEN
 950                            S = ONE + S*( SCALE / AA )**2
 951                            SCALE = AA
 952                         ELSE
 953                            S = S + ( AA / SCALE )**2
 954                         END IF
 955                      END IF
 956                      AA = DBLE( A( L+1 ) )
 957 *                    U(i,i)
 958                      IF( AA.NE.ZERO ) THEN
 959                         IF( SCALE.LT.AA ) THEN
 960                            S = ONE + S*( SCALE / AA )**2
 961                            SCALE = AA
 962                         ELSE
 963                            S = S + ( AA / SCALE )**2
 964                         END IF
 965                      END IF
 966                      L = L + LDA + 1
 967                   END DO
 968                   AA = DBLE( A( L ) )
 969 *                 U(n-1,n-1)
 970                   IF( AA.NE.ZERO ) THEN
 971                      IF( SCALE.LT.AA ) THEN
 972                         S = ONE + S*( SCALE / AA )**2
 973                         SCALE = AA
 974                      ELSE
 975                         S = S + ( AA / SCALE )**2
 976                      END IF
 977                   END IF
 978                ELSE
 979 *                 ilu=1 & A is lower
 980                   DO J = 0, K - 1
 981                      CALL ZLASSQ( N-J-1, A( J+1+J*LDA ), 1, SCALE, S )
 982 *                    trap L at A(0,0)
 983                   END DO
 984                   DO J = 1, K - 2
 985                      CALL ZLASSQ( J, A( 0+( 1+J )*LDA ), 1, SCALE, S )
 986 *                    U at A(0,1)
 987                   END DO
 988                   S = S + S
 989 *                 double s for the off diagonal elements
 990                   AA = DBLE( A( 0 ) )
 991 *                 L(0,0) at A(0,0)
 992                   IF( AA.NE.ZERO ) THEN
 993                      IF( SCALE.LT.AA ) THEN
 994                         S = ONE + S*( SCALE / AA )**2
 995                         SCALE = AA
 996                      ELSE
 997                         S = S + ( AA / SCALE )**2
 998                      END IF
 999                   END IF
1000                   L = LDA
1001 *                 -> L(k,k) at A(0,1)
1002                   DO I = 1, K - 1
1003                      AA = DBLE( A( L ) )
1004 *                    L(k-1+i,k-1+i)
1005                      IF( AA.NE.ZERO ) THEN
1006                         IF( SCALE.LT.AA ) THEN
1007                            S = ONE + S*( SCALE / AA )**2
1008                            SCALE = AA
1009                         ELSE
1010                            S = S + ( AA / SCALE )**2
1011                         END IF
1012                      END IF
1013                      AA = DBLE( A( L+1 ) )
1014 *                    L(i,i)
1015                      IF( AA.NE.ZERO ) THEN
1016                         IF( SCALE.LT.AA ) THEN
1017                            S = ONE + S*( SCALE / AA )**2
1018                            SCALE = AA
1019                         ELSE
1020                            S = S + ( AA / SCALE )**2
1021                         END IF
1022                      END IF
1023                      L = L + LDA + 1
1024                   END DO
1025                END IF
1026             ELSE
1027 *              A is xpose & A is k by n
1028                IF( ILU.EQ.0 ) THEN
1029 *                 A**H is upper
1030                   DO J = 1, K - 2
1031                      CALL ZLASSQ( J, A( 0+( K+J )*LDA ), 1, SCALE, S )
1032 *                    U at A(0,k)
1033                   END DO
1034                   DO J = 0, K - 2
1035                      CALL ZLASSQ( K, A( 0+J*LDA ), 1, SCALE, S )
1036 *                    k by k-1 rect. at A(0,0)
1037                   END DO
1038                   DO J = 0, K - 2
1039                      CALL ZLASSQ( K-J-1, A( J+1+( J+K-1 )*LDA ), 1,
1040      $                            SCALE, S )
1041 *                    L at A(0,k-1)
1042                   END DO
1043                   S = S + S
1044 *                 double s for the off diagonal elements
1045                   L = 0 + K*LDA - LDA
1046 *                 -> U(k-1,k-1) at A(0,k-1)
1047                   AA = DBLE( A( L ) )
1048 *                 U(k-1,k-1)
1049                   IF( AA.NE.ZERO ) THEN
1050                      IF( SCALE.LT.AA ) THEN
1051                         S = ONE + S*( SCALE / AA )**2
1052                         SCALE = AA
1053                      ELSE
1054                         S = S + ( AA / SCALE )**2
1055                      END IF
1056                   END IF
1057                   L = L + LDA
1058 *                 -> U(0,0) at A(0,k)
1059                   DO J = K, N - 1
1060                      AA = DBLE( A( L ) )
1061 *                    -> U(j-k,j-k)
1062                      IF( AA.NE.ZERO ) THEN
1063                         IF( SCALE.LT.AA ) THEN
1064                            S = ONE + S*( SCALE / AA )**2
1065                            SCALE = AA
1066                         ELSE
1067                            S = S + ( AA / SCALE )**2
1068                         END IF
1069                      END IF
1070                      AA = DBLE( A( L+1 ) )
1071 *                    -> U(j,j)
1072                      IF( AA.NE.ZERO ) THEN
1073                         IF( SCALE.LT.AA ) THEN
1074                            S = ONE + S*( SCALE / AA )**2
1075                            SCALE = AA
1076                         ELSE
1077                            S = S + ( AA / SCALE )**2
1078                         END IF
1079                      END IF
1080                      L = L + LDA + 1
1081                   END DO
1082                ELSE
1083 *                 A**H is lower
1084                   DO J = 1, K - 1
1085                      CALL ZLASSQ( J, A( 0+J*LDA ), 1, SCALE, S )
1086 *                    U at A(0,0)
1087                   END DO
1088                   DO J = K, N - 1
1089                      CALL ZLASSQ( K, A( 0+J*LDA ), 1, SCALE, S )
1090 *                    k by k-1 rect. at A(0,k)
1091                   END DO
1092                   DO J = 0, K - 3
1093                      CALL ZLASSQ( K-J-2, A( J+2+J*LDA ), 1, SCALE, S )
1094 *                    L at A(1,0)
1095                   END DO
1096                   S = S + S
1097 *                 double s for the off diagonal elements
1098                   L = 0
1099 *                 -> L(0,0) at A(0,0)
1100                   DO I = 0, K - 2
1101                      AA = DBLE( A( L ) )
1102 *                    L(i,i)
1103                      IF( AA.NE.ZERO ) THEN
1104                         IF( SCALE.LT.AA ) THEN
1105                            S = ONE + S*( SCALE / AA )**2
1106                            SCALE = AA
1107                         ELSE
1108                            S = S + ( AA / SCALE )**2
1109                         END IF
1110                      END IF
1111                      AA = DBLE( A( L+1 ) )
1112 *                    L(k+i,k+i)
1113                      IF( AA.NE.ZERO ) THEN
1114                         IF( SCALE.LT.AA ) THEN
1115                            S = ONE + S*( SCALE / AA )**2
1116                            SCALE = AA
1117                         ELSE
1118                            S = S + ( AA / SCALE )**2
1119                         END IF
1120                      END IF
1121                      L = L + LDA + 1
1122                   END DO
1123 *                 L-> k-1 + (k-1)*lda or L(k-1,k-1) at A(k-1,k-1)
1124                   AA = DBLE( A( L ) )
1125 *                 L(k-1,k-1) at A(k-1,k-1)
1126                   IF( AA.NE.ZERO ) THEN
1127                      IF( SCALE.LT.AA ) THEN
1128                         S = ONE + S*( SCALE / AA )**2
1129                         SCALE = AA
1130                      ELSE
1131                         S = S + ( AA / SCALE )**2
1132                      END IF
1133                   END IF
1134                END IF
1135             END IF
1136          ELSE
1137 *           n is even
1138             IF( IFM.EQ.1 ) THEN
1139 *              A is normal
1140                IF( ILU.EQ.0 ) THEN
1141 *                 A is upper
1142                   DO J = 0, K - 2
1143                      CALL ZLASSQ( K-J-1, A( K+J+2+J*LDA ), 1, SCALE, S )
1144 *                 L at A(k+1,0)
1145                   END DO
1146                   DO J = 0, K - 1
1147                      CALL ZLASSQ( K+J, A( 0+J*LDA ), 1, SCALE, S )
1148 *                 trap U at A(0,0)
1149                   END DO
1150                   S = S + S
1151 *                 double s for the off diagonal elements
1152                   L = K
1153 *                 -> U(k,k) at A(k,0)
1154                   DO I = 0, K - 1
1155                      AA = DBLE( A( L ) )
1156 *                    U(k+i,k+i)
1157                      IF( AA.NE.ZERO ) THEN
1158                         IF( SCALE.LT.AA ) THEN
1159                            S = ONE + S*( SCALE / AA )**2
1160                            SCALE = AA
1161                         ELSE
1162                            S = S + ( AA / SCALE )**2
1163                         END IF
1164                      END IF
1165                      AA = DBLE( A( L+1 ) )
1166 *                    U(i,i)
1167                      IF( AA.NE.ZERO ) THEN
1168                         IF( SCALE.LT.AA ) THEN
1169                            S = ONE + S*( SCALE / AA )**2
1170                            SCALE = AA
1171                         ELSE
1172                            S = S + ( AA / SCALE )**2
1173                         END IF
1174                      END IF
1175                      L = L + LDA + 1
1176                   END DO
1177                ELSE
1178 *                 ilu=1 & A is lower
1179                   DO J = 0, K - 1
1180                      CALL ZLASSQ( N-J-1, A( J+2+J*LDA ), 1, SCALE, S )
1181 *                    trap L at A(1,0)
1182                   END DO
1183                   DO J = 1, K - 1
1184                      CALL ZLASSQ( J, A( 0+J*LDA ), 1, SCALE, S )
1185 *                    U at A(0,0)
1186                   END DO
1187                   S = S + S
1188 *                 double s for the off diagonal elements
1189                   L = 0
1190 *                 -> L(k,k) at A(0,0)
1191                   DO I = 0, K - 1
1192                      AA = DBLE( A( L ) )
1193 *                    L(k-1+i,k-1+i)
1194                      IF( AA.NE.ZERO ) THEN
1195                         IF( SCALE.LT.AA ) THEN
1196                            S = ONE + S*( SCALE / AA )**2
1197                            SCALE = AA
1198                         ELSE
1199                            S = S + ( AA / SCALE )**2
1200                         END IF
1201                      END IF
1202                      AA = DBLE( A( L+1 ) )
1203 *                    L(i,i)
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 = L + LDA + 1
1213                   END DO
1214                END IF
1215             ELSE
1216 *              A is xpose
1217                IF( ILU.EQ.0 ) THEN
1218 *                 A**H is upper
1219                   DO J = 1, K - 1
1220                      CALL ZLASSQ( J, A( 0+( K+1+J )*LDA ), 1, SCALE, S )
1221 *                 U at A(0,k+1)
1222                   END DO
1223                   DO J = 0, K - 1
1224                      CALL ZLASSQ( K, A( 0+J*LDA ), 1, SCALE, S )
1225 *                 k by k rect. at A(0,0)
1226                   END DO
1227                   DO J = 0, K - 2
1228                      CALL ZLASSQ( K-J-1, A( J+1+( J+K )*LDA ), 1, SCALE,
1229      $                            S )
1230 *                 L at A(0,k)
1231                   END DO
1232                   S = S + S
1233 *                 double s for the off diagonal elements
1234                   L = 0 + K*LDA
1235 *                 -> U(k,k) at A(0,k)
1236                   AA = DBLE( A( L ) )
1237 *                 U(k,k)
1238                   IF( AA.NE.ZERO ) THEN
1239                      IF( SCALE.LT.AA ) THEN
1240                         S = ONE + S*( SCALE / AA )**2
1241                         SCALE = AA
1242                      ELSE
1243                         S = S + ( AA / SCALE )**2
1244                      END IF
1245                   END IF
1246                   L = L + LDA
1247 *                 -> U(0,0) at A(0,k+1)
1248                   DO J = K + 1, N - 1
1249                      AA = DBLE( A( L ) )
1250 *                    -> U(j-k-1,j-k-1)
1251                      IF( AA.NE.ZERO ) THEN
1252                         IF( SCALE.LT.AA ) THEN
1253                            S = ONE + S*( SCALE / AA )**2
1254                            SCALE = AA
1255                         ELSE
1256                            S = S + ( AA / SCALE )**2
1257                         END IF
1258                      END IF
1259                      AA = DBLE( A( L+1 ) )
1260 *                    -> U(j,j)
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 + 1
1270                   END DO
1271 *                 L=k-1+n*lda
1272 *                 -> U(k-1,k-1) at A(k-1,n)
1273                   AA = DBLE( A( L ) )
1274 *                 U(k,k)
1275                   IF( AA.NE.ZERO ) THEN
1276                      IF( SCALE.LT.AA ) THEN
1277                         S = ONE + S*( SCALE / AA )**2
1278                         SCALE = AA
1279                      ELSE
1280                         S = S + ( AA / SCALE )**2
1281                      END IF
1282                   END IF
1283                ELSE
1284 *                 A**H is lower
1285                   DO J = 1, K - 1
1286                      CALL ZLASSQ( J, A( 0+( J+1 )*LDA ), 1, SCALE, S )
1287 *                 U at A(0,1)
1288                   END DO
1289                   DO J = K + 1, N
1290                      CALL ZLASSQ( K, A( 0+J*LDA ), 1, SCALE, S )
1291 *                 k by k rect. at A(0,k+1)
1292                   END DO
1293                   DO J = 0, K - 2
1294                      CALL ZLASSQ( K-J-1, A( J+1+J*LDA ), 1, SCALE, S )
1295 *                 L at A(0,0)
1296                   END DO
1297                   S = S + S
1298 *                 double s for the off diagonal elements
1299                   L = 0
1300 *                 -> L(k,k) at A(0,0)
1301                   AA = DBLE( A( L ) )
1302 *                 L(k,k) at A(0,0)
1303                   IF( AA.NE.ZERO ) THEN
1304                      IF( SCALE.LT.AA ) THEN
1305                         S = ONE + S*( SCALE / AA )**2
1306                         SCALE = AA
1307                      ELSE
1308                         S = S + ( AA / SCALE )**2
1309                      END IF
1310                   END IF
1311                   L = LDA
1312 *                 -> L(0,0) at A(0,1)
1313                   DO I = 0, K - 2
1314                      AA = DBLE( A( L ) )
1315 *                    L(i,i)
1316                      IF( AA.NE.ZERO ) THEN
1317                         IF( SCALE.LT.AA ) THEN
1318                            S = ONE + S*( SCALE / AA )**2
1319                            SCALE = AA
1320                         ELSE
1321                            S = S + ( AA / SCALE )**2
1322                         END IF
1323                      END IF
1324                      AA = DBLE( A( L+1 ) )
1325 *                    L(k+i+1,k+i+1)
1326                      IF( AA.NE.ZERO ) THEN
1327                         IF( SCALE.LT.AA ) THEN
1328                            S = ONE + S*( SCALE / AA )**2
1329                            SCALE = AA
1330                         ELSE
1331                            S = S + ( AA / SCALE )**2
1332                         END IF
1333                      END IF
1334                      L = L + LDA + 1
1335                   END DO
1336 *                 L-> k - 1 + k*lda or L(k-1,k-1) at A(k-1,k)
1337                   AA = DBLE( A( L ) )
1338 *                 L(k-1,k-1) at A(k-1,k)
1339                   IF( AA.NE.ZERO ) THEN
1340                      IF( SCALE.LT.AA ) THEN
1341                         S = ONE + S*( SCALE / AA )**2
1342                         SCALE = AA
1343                      ELSE
1344                         S = S + ( AA / SCALE )**2
1345                      END IF
1346                   END IF
1347                END IF
1348             END IF
1349          END IF
1350          VALUE = SCALE*SQRT( S )
1351       END IF
1352 *
1353       ZLANHF = VALUE
1354       RETURN
1355 *
1356 *     End of ZLANHF
1357 *
1358       END
Home