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