LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
slansf.f
Go to the documentation of this file.
1*> \brief \b SLANSF
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download SLANSF + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slansf.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slansf.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slansf.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* REAL FUNCTION SLANSF( NORM, TRANSR, UPLO, N, A, WORK )
22*
23* .. Scalar Arguments ..
24* CHARACTER NORM, TRANSR, UPLO
25* INTEGER N
26* ..
27* .. Array Arguments ..
28* REAL A( 0: * ), WORK( 0: * )
29* ..
30*
31*
32*> \par Purpose:
33* =============
34*>
35*> \verbatim
36*>
37*> SLANSF returns the value of the one norm, or the Frobenius norm, or
38*> the infinity norm, or the element of largest absolute value of a
39*> real symmetric matrix A in RFP format.
40*> \endverbatim
41*>
42*> \return SLANSF
43*> \verbatim
44*>
45*> SLANSF = ( max(abs(A(i,j))), NORM = 'M' or 'm'
46*> (
47*> ( norm1(A), NORM = '1', 'O' or 'o'
48*> (
49*> ( normI(A), NORM = 'I' or 'i'
50*> (
51*> ( normF(A), NORM = 'F', 'f', 'E' or 'e'
52*>
53*> where norm1 denotes the one norm of a matrix (maximum column sum),
54*> normI denotes the infinity norm of a matrix (maximum row sum) and
55*> normF denotes the Frobenius norm of a matrix (square root of sum of
56*> squares). Note that max(abs(A(i,j))) is not a matrix norm.
57*> \endverbatim
58*
59* Arguments:
60* ==========
61*
62*> \param[in] NORM
63*> \verbatim
64*> NORM is CHARACTER*1
65*> Specifies the value to be returned in SLANSF as described
66*> above.
67*> \endverbatim
68*>
69*> \param[in] TRANSR
70*> \verbatim
71*> TRANSR is CHARACTER*1
72*> Specifies whether the RFP format of A is normal or
73*> transposed format.
74*> = 'N': RFP format is Normal;
75*> = 'T': RFP format is Transpose.
76*> \endverbatim
77*>
78*> \param[in] UPLO
79*> \verbatim
80*> UPLO is CHARACTER*1
81*> On entry, UPLO specifies whether the RFP matrix A came from
82*> an upper or lower triangular matrix as follows:
83*> = 'U': RFP A came from an upper triangular matrix;
84*> = 'L': RFP A came from a lower triangular matrix.
85*> \endverbatim
86*>
87*> \param[in] N
88*> \verbatim
89*> N is INTEGER
90*> The order of the matrix A. N >= 0. When N = 0, SLANSF is
91*> set to zero.
92*> \endverbatim
93*>
94*> \param[in] A
95*> \verbatim
96*> A is REAL array, dimension ( N*(N+1)/2 );
97*> On entry, the upper (if UPLO = 'U') or lower (if UPLO = 'L')
98*> part of the symmetric matrix A stored in RFP format. See the
99*> "Notes" below for more details.
100*> Unchanged on exit.
101*> \endverbatim
102*>
103*> \param[out] WORK
104*> \verbatim
105*> WORK is REAL array, dimension (MAX(1,LWORK)),
106*> where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
107*> WORK is not referenced.
108*> \endverbatim
109*
110* Authors:
111* ========
112*
113*> \author Univ. of Tennessee
114*> \author Univ. of California Berkeley
115*> \author Univ. of Colorado Denver
116*> \author NAG Ltd.
117*
118*> \ingroup lanhf
119*
120*> \par Further Details:
121* =====================
122*>
123*> \verbatim
124*>
125*> We first consider Rectangular Full Packed (RFP) Format when N is
126*> even. We give an example where N = 6.
127*>
128*> AP is Upper AP is Lower
129*>
130*> 00 01 02 03 04 05 00
131*> 11 12 13 14 15 10 11
132*> 22 23 24 25 20 21 22
133*> 33 34 35 30 31 32 33
134*> 44 45 40 41 42 43 44
135*> 55 50 51 52 53 54 55
136*>
137*>
138*> Let TRANSR = 'N'. RFP holds AP as follows:
139*> For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
140*> three columns of AP upper. The lower triangle A(4:6,0:2) consists of
141*> the transpose of the first three columns of AP upper.
142*> For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
143*> three columns of AP lower. The upper triangle A(0:2,0:2) consists of
144*> the transpose of the last three columns of AP lower.
145*> This covers the case N even and TRANSR = 'N'.
146*>
147*> RFP A RFP A
148*>
149*> 03 04 05 33 43 53
150*> 13 14 15 00 44 54
151*> 23 24 25 10 11 55
152*> 33 34 35 20 21 22
153*> 00 44 45 30 31 32
154*> 01 11 55 40 41 42
155*> 02 12 22 50 51 52
156*>
157*> Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
158*> transpose of RFP A above. One therefore gets:
159*>
160*>
161*> RFP A RFP A
162*>
163*> 03 13 23 33 00 01 02 33 00 10 20 30 40 50
164*> 04 14 24 34 44 11 12 43 44 11 21 31 41 51
165*> 05 15 25 35 45 55 22 53 54 55 22 32 42 52
166*>
167*>
168*> We then consider Rectangular Full Packed (RFP) Format when N is
169*> odd. We give an example where N = 5.
170*>
171*> AP is Upper AP is Lower
172*>
173*> 00 01 02 03 04 00
174*> 11 12 13 14 10 11
175*> 22 23 24 20 21 22
176*> 33 34 30 31 32 33
177*> 44 40 41 42 43 44
178*>
179*>
180*> Let TRANSR = 'N'. RFP holds AP as follows:
181*> For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
182*> three columns of AP upper. The lower triangle A(3:4,0:1) consists of
183*> the transpose of the first two columns of AP upper.
184*> For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
185*> three columns of AP lower. The upper triangle A(0:1,1:2) consists of
186*> the transpose of the last two columns of AP lower.
187*> This covers the case N odd and TRANSR = 'N'.
188*>
189*> RFP A RFP A
190*>
191*> 02 03 04 00 33 43
192*> 12 13 14 10 11 44
193*> 22 23 24 20 21 22
194*> 00 33 34 30 31 32
195*> 01 11 44 40 41 42
196*>
197*> Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
198*> transpose of RFP A above. One therefore gets:
199*>
200*> RFP A RFP A
201*>
202*> 02 12 22 00 01 00 10 20 30 40 50
203*> 03 13 23 33 11 33 11 21 31 41 51
204*> 04 14 24 34 44 43 44 22 32 42 52
205*> \endverbatim
206*
207* =====================================================================
208 REAL function slansf( norm, transr, uplo, n, a, work )
209*
210* -- LAPACK computational routine --
211* -- LAPACK is a software package provided by Univ. of Tennessee, --
212* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
213*
214* .. Scalar Arguments ..
215 CHARACTER norm, transr, uplo
216 INTEGER n
217* ..
218* .. Array Arguments ..
219 REAL a( 0: * ), work( 0: * )
220* ..
221*
222* =====================================================================
223*
224* ..
225* .. Parameters ..
226 REAL one, zero
227 parameter( one = 1.0e+0, zero = 0.0e+0 )
228* ..
229* .. Local Scalars ..
230 INTEGER i, j, ifm, ilu, noe, n1, k, l, lda
231 REAL scale, s, VALUE, aa, temp
232* ..
233* .. External Functions ..
234 LOGICAL lsame, sisnan
235 EXTERNAL lsame, sisnan
236* ..
237* .. External Subroutines ..
238 EXTERNAL slassq
239* ..
240* .. Intrinsic Functions ..
241 INTRINSIC abs, sqrt
242* ..
243* .. Executable Statements ..
244*
245 IF( n.EQ.0 ) THEN
246 slansf = zero
247 RETURN
248 ELSE IF( n.EQ.1 ) THEN
249 slansf = abs( a(0) )
250 RETURN
251 END IF
252*
253* set noe = 1 if n is odd. if n is even set noe=0
254*
255 noe = 1
256 IF( mod( n, 2 ).EQ.0 )
257 $ noe = 0
258*
259* set ifm = 0 when form='T or 't' and 1 otherwise
260*
261 ifm = 1
262 IF( lsame( transr, 'T' ) )
263 $ ifm = 0
264*
265* set ilu = 0 when uplo='U or 'u' and 1 otherwise
266*
267 ilu = 1
268 IF( lsame( uplo, 'U' ) )
269 $ ilu = 0
270*
271* set lda = (n+1)/2 when ifm = 0
272* set lda = n when ifm = 1 and noe = 1
273* set lda = n+1 when ifm = 1 and noe = 0
274*
275 IF( ifm.EQ.1 ) THEN
276 IF( noe.EQ.1 ) THEN
277 lda = n
278 ELSE
279* noe=0
280 lda = n + 1
281 END IF
282 ELSE
283* ifm=0
284 lda = ( n+1 ) / 2
285 END IF
286*
287 IF( lsame( norm, 'M' ) ) THEN
288*
289* Find max(abs(A(i,j))).
290*
291 k = ( n+1 ) / 2
292 VALUE = zero
293 IF( noe.EQ.1 ) THEN
294* n is odd
295 IF( ifm.EQ.1 ) THEN
296* A is n by k
297 DO j = 0, k - 1
298 DO i = 0, n - 1
299 temp = abs( a( i+j*lda ) )
300 IF( VALUE .LT. temp .OR. sisnan( temp ) )
301 $ VALUE = temp
302 END DO
303 END DO
304 ELSE
305* xpose case; A is k by n
306 DO j = 0, n - 1
307 DO i = 0, k - 1
308 temp = abs( a( i+j*lda ) )
309 IF( VALUE .LT. temp .OR. sisnan( temp ) )
310 $ VALUE = temp
311 END DO
312 END DO
313 END IF
314 ELSE
315* n is even
316 IF( ifm.EQ.1 ) THEN
317* A is n+1 by k
318 DO j = 0, k - 1
319 DO i = 0, n
320 temp = abs( a( i+j*lda ) )
321 IF( VALUE .LT. temp .OR. sisnan( temp ) )
322 $ VALUE = temp
323 END DO
324 END DO
325 ELSE
326* xpose case; A is k by n+1
327 DO j = 0, n
328 DO i = 0, k - 1
329 temp = abs( a( i+j*lda ) )
330 IF( VALUE .LT. temp .OR. sisnan( temp ) )
331 $ VALUE = temp
332 END DO
333 END DO
334 END IF
335 END IF
336 ELSE IF( ( lsame( norm, 'I' ) ) .OR. ( lsame( norm, 'O' ) ) .OR.
337 $ ( norm.EQ.'1' ) ) THEN
338*
339* Find normI(A) ( = norm1(A), since A is symmetric).
340*
341 IF( ifm.EQ.1 ) THEN
342 k = n / 2
343 IF( noe.EQ.1 ) THEN
344* n is odd
345 IF( ilu.EQ.0 ) THEN
346 DO i = 0, k - 1
347 work( i ) = zero
348 END DO
349 DO j = 0, k
350 s = zero
351 DO i = 0, k + j - 1
352 aa = abs( a( i+j*lda ) )
353* -> A(i,j+k)
354 s = s + aa
355 work( i ) = work( i ) + aa
356 END DO
357 aa = abs( a( i+j*lda ) )
358* -> A(j+k,j+k)
359 work( j+k ) = s + aa
360 IF( i.EQ.k+k )
361 $ GO TO 10
362 i = i + 1
363 aa = abs( a( i+j*lda ) )
364* -> A(j,j)
365 work( j ) = work( j ) + aa
366 s = zero
367 DO l = j + 1, k - 1
368 i = i + 1
369 aa = abs( a( i+j*lda ) )
370* -> A(l,j)
371 s = s + aa
372 work( l ) = work( l ) + aa
373 END DO
374 work( j ) = work( j ) + s
375 END DO
376 10 CONTINUE
377 VALUE = work( 0 )
378 DO i = 1, n-1
379 temp = work( i )
380 IF( VALUE .LT. temp .OR. sisnan( temp ) )
381 $ VALUE = temp
382 END DO
383 ELSE
384* ilu = 1
385 k = k + 1
386* k=(n+1)/2 for n odd and ilu=1
387 DO i = k, n - 1
388 work( i ) = zero
389 END DO
390 DO j = k - 1, 0, -1
391 s = zero
392 DO i = 0, j - 2
393 aa = abs( a( i+j*lda ) )
394* -> A(j+k,i+k)
395 s = s + aa
396 work( i+k ) = work( i+k ) + aa
397 END DO
398 IF( j.GT.0 ) THEN
399 aa = abs( a( i+j*lda ) )
400* -> A(j+k,j+k)
401 s = s + aa
402 work( i+k ) = work( i+k ) + s
403* i=j
404 i = i + 1
405 END IF
406 aa = abs( a( i+j*lda ) )
407* -> A(j,j)
408 work( j ) = aa
409 s = zero
410 DO l = j + 1, n - 1
411 i = i + 1
412 aa = abs( a( i+j*lda ) )
413* -> A(l,j)
414 s = s + aa
415 work( l ) = work( l ) + aa
416 END DO
417 work( j ) = work( j ) + s
418 END DO
419 VALUE = work( 0 )
420 DO i = 1, n-1
421 temp = work( i )
422 IF( VALUE .LT. temp .OR. sisnan( temp ) )
423 $ VALUE = temp
424 END DO
425 END IF
426 ELSE
427* n is even
428 IF( ilu.EQ.0 ) THEN
429 DO i = 0, k - 1
430 work( i ) = zero
431 END DO
432 DO j = 0, k - 1
433 s = zero
434 DO i = 0, k + j - 1
435 aa = abs( a( i+j*lda ) )
436* -> A(i,j+k)
437 s = s + aa
438 work( i ) = work( i ) + aa
439 END DO
440 aa = abs( a( i+j*lda ) )
441* -> A(j+k,j+k)
442 work( j+k ) = s + aa
443 i = i + 1
444 aa = abs( a( i+j*lda ) )
445* -> A(j,j)
446 work( j ) = work( j ) + aa
447 s = zero
448 DO l = j + 1, k - 1
449 i = i + 1
450 aa = abs( a( i+j*lda ) )
451* -> A(l,j)
452 s = s + aa
453 work( l ) = work( l ) + aa
454 END DO
455 work( j ) = work( j ) + s
456 END DO
457 VALUE = work( 0 )
458 DO i = 1, n-1
459 temp = work( i )
460 IF( VALUE .LT. temp .OR. sisnan( temp ) )
461 $ VALUE = temp
462 END DO
463 ELSE
464* ilu = 1
465 DO i = k, n - 1
466 work( i ) = zero
467 END DO
468 DO j = k - 1, 0, -1
469 s = zero
470 DO i = 0, j - 1
471 aa = abs( a( i+j*lda ) )
472* -> A(j+k,i+k)
473 s = s + aa
474 work( i+k ) = work( i+k ) + aa
475 END DO
476 aa = abs( a( i+j*lda ) )
477* -> A(j+k,j+k)
478 s = s + aa
479 work( i+k ) = work( i+k ) + s
480* i=j
481 i = i + 1
482 aa = abs( a( i+j*lda ) )
483* -> A(j,j)
484 work( j ) = aa
485 s = zero
486 DO l = j + 1, n - 1
487 i = i + 1
488 aa = abs( a( i+j*lda ) )
489* -> A(l,j)
490 s = s + aa
491 work( l ) = work( l ) + aa
492 END DO
493 work( j ) = work( j ) + s
494 END DO
495 VALUE = work( 0 )
496 DO i = 1, n-1
497 temp = work( i )
498 IF( VALUE .LT. temp .OR. sisnan( temp ) )
499 $ VALUE = temp
500 END DO
501 END IF
502 END IF
503 ELSE
504* ifm=0
505 k = n / 2
506 IF( noe.EQ.1 ) THEN
507* n is odd
508 IF( ilu.EQ.0 ) THEN
509 n1 = k
510* n/2
511 k = k + 1
512* k is the row size and lda
513 DO i = n1, n - 1
514 work( i ) = zero
515 END DO
516 DO j = 0, n1 - 1
517 s = zero
518 DO i = 0, k - 1
519 aa = abs( a( i+j*lda ) )
520* A(j,n1+i)
521 work( i+n1 ) = work( i+n1 ) + aa
522 s = s + aa
523 END DO
524 work( j ) = s
525 END DO
526* j=n1=k-1 is special
527 s = abs( a( 0+j*lda ) )
528* A(k-1,k-1)
529 DO i = 1, k - 1
530 aa = abs( a( i+j*lda ) )
531* A(k-1,i+n1)
532 work( i+n1 ) = work( i+n1 ) + aa
533 s = s + aa
534 END DO
535 work( j ) = work( j ) + s
536 DO j = k, n - 1
537 s = zero
538 DO i = 0, j - k - 1
539 aa = abs( a( i+j*lda ) )
540* A(i,j-k)
541 work( i ) = work( i ) + aa
542 s = s + aa
543 END DO
544* i=j-k
545 aa = abs( a( i+j*lda ) )
546* A(j-k,j-k)
547 s = s + aa
548 work( j-k ) = work( j-k ) + s
549 i = i + 1
550 s = abs( a( i+j*lda ) )
551* A(j,j)
552 DO l = j + 1, n - 1
553 i = i + 1
554 aa = abs( a( i+j*lda ) )
555* A(j,l)
556 work( l ) = work( l ) + aa
557 s = s + aa
558 END DO
559 work( j ) = work( j ) + s
560 END DO
561 VALUE = work( 0 )
562 DO i = 1, n-1
563 temp = work( i )
564 IF( VALUE .LT. temp .OR. sisnan( temp ) )
565 $ VALUE = temp
566 END DO
567 ELSE
568* ilu=1
569 k = k + 1
570* k=(n+1)/2 for n odd and ilu=1
571 DO i = k, n - 1
572 work( i ) = zero
573 END DO
574 DO j = 0, k - 2
575* process
576 s = zero
577 DO i = 0, j - 1
578 aa = abs( a( i+j*lda ) )
579* A(j,i)
580 work( i ) = work( i ) + aa
581 s = s + aa
582 END DO
583 aa = abs( a( i+j*lda ) )
584* i=j so process of A(j,j)
585 s = s + aa
586 work( j ) = s
587* is initialised here
588 i = i + 1
589* i=j process A(j+k,j+k)
590 aa = abs( a( i+j*lda ) )
591 s = aa
592 DO l = k + j + 1, n - 1
593 i = i + 1
594 aa = abs( a( i+j*lda ) )
595* A(l,k+j)
596 s = s + aa
597 work( l ) = work( l ) + aa
598 END DO
599 work( k+j ) = work( k+j ) + s
600 END DO
601* j=k-1 is special :process col A(k-1,0:k-1)
602 s = zero
603 DO i = 0, k - 2
604 aa = abs( a( i+j*lda ) )
605* A(k,i)
606 work( i ) = work( i ) + aa
607 s = s + aa
608 END DO
609* i=k-1
610 aa = abs( a( i+j*lda ) )
611* A(k-1,k-1)
612 s = s + aa
613 work( i ) = s
614* done with col j=k+1
615 DO j = k, n - 1
616* process col j of A = A(j,0:k-1)
617 s = zero
618 DO i = 0, k - 1
619 aa = abs( a( i+j*lda ) )
620* A(j,i)
621 work( i ) = work( i ) + aa
622 s = s + aa
623 END DO
624 work( j ) = work( j ) + s
625 END DO
626 VALUE = work( 0 )
627 DO i = 1, n-1
628 temp = work( i )
629 IF( VALUE .LT. temp .OR. sisnan( temp ) )
630 $ VALUE = temp
631 END DO
632 END IF
633 ELSE
634* n is even
635 IF( ilu.EQ.0 ) THEN
636 DO i = k, n - 1
637 work( i ) = zero
638 END DO
639 DO j = 0, k - 1
640 s = zero
641 DO i = 0, k - 1
642 aa = abs( a( i+j*lda ) )
643* A(j,i+k)
644 work( i+k ) = work( i+k ) + aa
645 s = s + aa
646 END DO
647 work( j ) = s
648 END DO
649* j=k
650 aa = abs( a( 0+j*lda ) )
651* A(k,k)
652 s = aa
653 DO i = 1, k - 1
654 aa = abs( a( i+j*lda ) )
655* A(k,k+i)
656 work( i+k ) = work( i+k ) + aa
657 s = s + aa
658 END DO
659 work( j ) = work( j ) + s
660 DO j = k + 1, n - 1
661 s = zero
662 DO i = 0, j - 2 - k
663 aa = abs( a( i+j*lda ) )
664* A(i,j-k-1)
665 work( i ) = work( i ) + aa
666 s = s + aa
667 END DO
668* i=j-1-k
669 aa = abs( a( i+j*lda ) )
670* A(j-k-1,j-k-1)
671 s = s + aa
672 work( j-k-1 ) = work( j-k-1 ) + s
673 i = i + 1
674 aa = abs( a( i+j*lda ) )
675* A(j,j)
676 s = aa
677 DO l = j + 1, n - 1
678 i = i + 1
679 aa = abs( a( i+j*lda ) )
680* A(j,l)
681 work( l ) = work( l ) + aa
682 s = s + aa
683 END DO
684 work( j ) = work( j ) + s
685 END DO
686* j=n
687 s = zero
688 DO i = 0, k - 2
689 aa = abs( a( i+j*lda ) )
690* A(i,k-1)
691 work( i ) = work( i ) + aa
692 s = s + aa
693 END DO
694* i=k-1
695 aa = abs( a( i+j*lda ) )
696* A(k-1,k-1)
697 s = s + aa
698 work( i ) = work( i ) + s
699 VALUE = work( 0 )
700 DO i = 1, n-1
701 temp = work( i )
702 IF( VALUE .LT. temp .OR. sisnan( temp ) )
703 $ VALUE = temp
704 END DO
705 ELSE
706* ilu=1
707 DO i = k, n - 1
708 work( i ) = zero
709 END DO
710* j=0 is special :process col A(k:n-1,k)
711 s = abs( a( 0 ) )
712* A(k,k)
713 DO i = 1, k - 1
714 aa = abs( a( i ) )
715* A(k+i,k)
716 work( i+k ) = work( i+k ) + aa
717 s = s + aa
718 END DO
719 work( k ) = work( k ) + s
720 DO j = 1, k - 1
721* process
722 s = zero
723 DO i = 0, j - 2
724 aa = abs( a( i+j*lda ) )
725* A(j-1,i)
726 work( i ) = work( i ) + aa
727 s = s + aa
728 END DO
729 aa = abs( a( i+j*lda ) )
730* i=j-1 so process of A(j-1,j-1)
731 s = s + aa
732 work( j-1 ) = s
733* is initialised here
734 i = i + 1
735* i=j process A(j+k,j+k)
736 aa = abs( a( i+j*lda ) )
737 s = aa
738 DO l = k + j + 1, n - 1
739 i = i + 1
740 aa = abs( a( i+j*lda ) )
741* A(l,k+j)
742 s = s + aa
743 work( l ) = work( l ) + aa
744 END DO
745 work( k+j ) = work( k+j ) + s
746 END DO
747* j=k is special :process col A(k,0:k-1)
748 s = zero
749 DO i = 0, k - 2
750 aa = abs( a( i+j*lda ) )
751* A(k,i)
752 work( i ) = work( i ) + aa
753 s = s + aa
754 END DO
755* i=k-1
756 aa = abs( a( i+j*lda ) )
757* A(k-1,k-1)
758 s = s + aa
759 work( i ) = s
760* done with col j=k+1
761 DO j = k + 1, n
762* process col j-1 of A = A(j-1,0:k-1)
763 s = zero
764 DO i = 0, k - 1
765 aa = abs( a( i+j*lda ) )
766* A(j-1,i)
767 work( i ) = work( i ) + aa
768 s = s + aa
769 END DO
770 work( j-1 ) = work( j-1 ) + s
771 END DO
772 VALUE = work( 0 )
773 DO i = 1, n-1
774 temp = work( i )
775 IF( VALUE .LT. temp .OR. sisnan( temp ) )
776 $ VALUE = temp
777 END DO
778 END IF
779 END IF
780 END IF
781 ELSE IF( ( lsame( norm, 'F' ) ) .OR. ( lsame( norm, 'E' ) ) ) THEN
782*
783* Find normF(A).
784*
785 k = ( n+1 ) / 2
786 scale = zero
787 s = one
788 IF( noe.EQ.1 ) THEN
789* n is odd
790 IF( ifm.EQ.1 ) THEN
791* A is normal
792 IF( ilu.EQ.0 ) THEN
793* A is upper
794 DO j = 0, k - 3
795 CALL slassq( k-j-2, a( k+j+1+j*lda ), 1, scale, s )
796* L at A(k,0)
797 END DO
798 DO j = 0, k - 1
799 CALL slassq( k+j-1, a( 0+j*lda ), 1, scale, s )
800* trap U at A(0,0)
801 END DO
802 s = s + s
803* double s for the off diagonal elements
804 CALL slassq( k-1, a( k ), lda+1, scale, s )
805* tri L at A(k,0)
806 CALL slassq( k, a( k-1 ), lda+1, scale, s )
807* tri U at A(k-1,0)
808 ELSE
809* ilu=1 & A is lower
810 DO j = 0, k - 1
811 CALL slassq( n-j-1, a( j+1+j*lda ), 1, scale, s )
812* trap L at A(0,0)
813 END DO
814 DO j = 0, k - 2
815 CALL slassq( j, a( 0+( 1+j )*lda ), 1, scale, s )
816* U at A(0,1)
817 END DO
818 s = s + s
819* double s for the off diagonal elements
820 CALL slassq( k, a( 0 ), lda+1, scale, s )
821* tri L at A(0,0)
822 CALL slassq( k-1, a( 0+lda ), lda+1, scale, s )
823* tri U at A(0,1)
824 END IF
825 ELSE
826* A is xpose
827 IF( ilu.EQ.0 ) THEN
828* A**T is upper
829 DO j = 1, k - 2
830 CALL slassq( j, a( 0+( k+j )*lda ), 1, scale, s )
831* U at A(0,k)
832 END DO
833 DO j = 0, k - 2
834 CALL slassq( k, a( 0+j*lda ), 1, scale, s )
835* k by k-1 rect. at A(0,0)
836 END DO
837 DO j = 0, k - 2
838 CALL slassq( k-j-1, a( j+1+( j+k-1 )*lda ), 1,
839 $ scale, s )
840* L at A(0,k-1)
841 END DO
842 s = s + s
843* double s for the off diagonal elements
844 CALL slassq( k-1, a( 0+k*lda ), lda+1, scale, s )
845* tri U at A(0,k)
846 CALL slassq( k, a( 0+( k-1 )*lda ), lda+1, scale, s )
847* tri L at A(0,k-1)
848 ELSE
849* A**T is lower
850 DO j = 1, k - 1
851 CALL slassq( j, a( 0+j*lda ), 1, scale, s )
852* U at A(0,0)
853 END DO
854 DO j = k, n - 1
855 CALL slassq( k, a( 0+j*lda ), 1, scale, s )
856* k by k-1 rect. at A(0,k)
857 END DO
858 DO j = 0, k - 3
859 CALL slassq( k-j-2, a( j+2+j*lda ), 1, scale, s )
860* L at A(1,0)
861 END DO
862 s = s + s
863* double s for the off diagonal elements
864 CALL slassq( k, a( 0 ), lda+1, scale, s )
865* tri U at A(0,0)
866 CALL slassq( k-1, a( 1 ), lda+1, scale, s )
867* tri L at A(1,0)
868 END IF
869 END IF
870 ELSE
871* n is even
872 IF( ifm.EQ.1 ) THEN
873* A is normal
874 IF( ilu.EQ.0 ) THEN
875* A is upper
876 DO j = 0, k - 2
877 CALL slassq( k-j-1, a( k+j+2+j*lda ), 1, scale, s )
878* L at A(k+1,0)
879 END DO
880 DO j = 0, k - 1
881 CALL slassq( k+j, a( 0+j*lda ), 1, scale, s )
882* trap U at A(0,0)
883 END DO
884 s = s + s
885* double s for the off diagonal elements
886 CALL slassq( k, a( k+1 ), lda+1, scale, s )
887* tri L at A(k+1,0)
888 CALL slassq( k, a( k ), lda+1, scale, s )
889* tri U at A(k,0)
890 ELSE
891* ilu=1 & A is lower
892 DO j = 0, k - 1
893 CALL slassq( n-j-1, a( j+2+j*lda ), 1, scale, s )
894* trap L at A(1,0)
895 END DO
896 DO j = 1, k - 1
897 CALL slassq( j, a( 0+j*lda ), 1, scale, s )
898* U at A(0,0)
899 END DO
900 s = s + s
901* double s for the off diagonal elements
902 CALL slassq( k, a( 1 ), lda+1, scale, s )
903* tri L at A(1,0)
904 CALL slassq( k, a( 0 ), lda+1, scale, s )
905* tri U at A(0,0)
906 END IF
907 ELSE
908* A is xpose
909 IF( ilu.EQ.0 ) THEN
910* A**T is upper
911 DO j = 1, k - 1
912 CALL slassq( j, a( 0+( k+1+j )*lda ), 1, scale, s )
913* U at A(0,k+1)
914 END DO
915 DO j = 0, k - 1
916 CALL slassq( k, a( 0+j*lda ), 1, scale, s )
917* k by k rect. at A(0,0)
918 END DO
919 DO j = 0, k - 2
920 CALL slassq( k-j-1, a( j+1+( j+k )*lda ), 1, scale,
921 $ s )
922* L at A(0,k)
923 END DO
924 s = s + s
925* double s for the off diagonal elements
926 CALL slassq( k, a( 0+( k+1 )*lda ), lda+1, scale, s )
927* tri U at A(0,k+1)
928 CALL slassq( k, a( 0+k*lda ), lda+1, scale, s )
929* tri L at A(0,k)
930 ELSE
931* A**T is lower
932 DO j = 1, k - 1
933 CALL slassq( j, a( 0+( j+1 )*lda ), 1, scale, s )
934* U at A(0,1)
935 END DO
936 DO j = k + 1, n
937 CALL slassq( k, a( 0+j*lda ), 1, scale, s )
938* k by k rect. at A(0,k+1)
939 END DO
940 DO j = 0, k - 2
941 CALL slassq( k-j-1, a( j+1+j*lda ), 1, scale, s )
942* L at A(0,0)
943 END DO
944 s = s + s
945* double s for the off diagonal elements
946 CALL slassq( k, a( lda ), lda+1, scale, s )
947* tri L at A(0,1)
948 CALL slassq( k, a( 0 ), lda+1, scale, s )
949* tri U at A(0,0)
950 END IF
951 END IF
952 END IF
953 VALUE = scale*sqrt( s )
954 END IF
955*
956 slansf = VALUE
957 RETURN
958*
959* End of SLANSF
960*
961 END
logical function sisnan(sin)
SISNAN tests input for NaN.
Definition sisnan.f:59
real function slansf(norm, transr, uplo, n, a, work)
SLANSF
Definition slansf.f:209
subroutine slassq(n, x, incx, scale, sumsq)
SLASSQ updates a sum of squares represented in scaled form.
Definition slassq.f90:124
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48