LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlansf.f
Go to the documentation of this file.
1*> \brief \b DLANSF returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a symmetric 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 DLANSF + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlansf.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlansf.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlansf.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* DOUBLE PRECISION FUNCTION DLANSF( NORM, TRANSR, UPLO, N, A, WORK )
22*
23* .. Scalar Arguments ..
24* CHARACTER NORM, TRANSR, UPLO
25* INTEGER N
26* ..
27* .. Array Arguments ..
28* DOUBLE PRECISION A( 0: * ), WORK( 0: * )
29* ..
30*
31*
32*> \par Purpose:
33* =============
34*>
35*> \verbatim
36*>
37*> DLANSF 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 DLANSF
43*> \verbatim
44*>
45*> DLANSF = ( 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 DLANSF 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, DLANSF is
91*> set to zero.
92*> \endverbatim
93*>
94*> \param[in] A
95*> \verbatim
96*> A is DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION FUNCTION dlansf( 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 DOUBLE PRECISION a( 0: * ), work( 0: * )
220* ..
221*
222* =====================================================================
223*
224* .. Parameters ..
225 DOUBLE PRECISION one, zero
226 parameter( one = 1.0d+0, zero = 0.0d+0 )
227* ..
228* .. Local Scalars ..
229 INTEGER i, j, ifm, ilu, noe, n1, k, l, lda
230 DOUBLE PRECISION scale, s, VALUE, aa, temp
231* ..
232* .. External Functions ..
233 LOGICAL lsame, disnan
234 EXTERNAL lsame, disnan
235* ..
236* .. External Subroutines ..
237 EXTERNAL dlassq
238* ..
239* .. Intrinsic Functions ..
240 INTRINSIC abs, max, sqrt
241* ..
242* .. Executable Statements ..
243*
244 IF( n.EQ.0 ) THEN
245 dlansf = zero
246 RETURN
247 ELSE IF( n.EQ.1 ) THEN
248 dlansf = abs( a(0) )
249 RETURN
250 END IF
251*
252* set noe = 1 if n is odd. if n is even set noe=0
253*
254 noe = 1
255 IF( mod( n, 2 ).EQ.0 )
256 $ noe = 0
257*
258* set ifm = 0 when form='T or 't' and 1 otherwise
259*
260 ifm = 1
261 IF( lsame( transr, 'T' ) )
262 $ ifm = 0
263*
264* set ilu = 0 when uplo='U or 'u' and 1 otherwise
265*
266 ilu = 1
267 IF( lsame( uplo, 'U' ) )
268 $ ilu = 0
269*
270* set lda = (n+1)/2 when ifm = 0
271* set lda = n when ifm = 1 and noe = 1
272* set lda = n+1 when ifm = 1 and noe = 0
273*
274 IF( ifm.EQ.1 ) THEN
275 IF( noe.EQ.1 ) THEN
276 lda = n
277 ELSE
278* noe=0
279 lda = n + 1
280 END IF
281 ELSE
282* ifm=0
283 lda = ( n+1 ) / 2
284 END IF
285*
286 IF( lsame( norm, 'M' ) ) THEN
287*
288* Find max(abs(A(i,j))).
289*
290 k = ( n+1 ) / 2
291 VALUE = zero
292 IF( noe.EQ.1 ) THEN
293* n is odd
294 IF( ifm.EQ.1 ) THEN
295* A is n by k
296 DO j = 0, k - 1
297 DO i = 0, n - 1
298 temp = abs( a( i+j*lda ) )
299 IF( VALUE .LT. temp .OR. disnan( temp ) )
300 $ VALUE = temp
301 END DO
302 END DO
303 ELSE
304* xpose case; A is k by n
305 DO j = 0, n - 1
306 DO i = 0, k - 1
307 temp = abs( a( i+j*lda ) )
308 IF( VALUE .LT. temp .OR. disnan( temp ) )
309 $ VALUE = temp
310 END DO
311 END DO
312 END IF
313 ELSE
314* n is even
315 IF( ifm.EQ.1 ) THEN
316* A is n+1 by k
317 DO j = 0, k - 1
318 DO i = 0, n
319 temp = abs( a( i+j*lda ) )
320 IF( VALUE .LT. temp .OR. disnan( temp ) )
321 $ VALUE = temp
322 END DO
323 END DO
324 ELSE
325* xpose case; A is k by n+1
326 DO j = 0, n
327 DO i = 0, k - 1
328 temp = abs( a( i+j*lda ) )
329 IF( VALUE .LT. temp .OR. disnan( temp ) )
330 $ VALUE = temp
331 END DO
332 END DO
333 END IF
334 END IF
335 ELSE IF( ( lsame( norm, 'I' ) ) .OR. ( lsame( norm, 'O' ) ) .OR.
336 $ ( norm.EQ.'1' ) ) THEN
337*
338* Find normI(A) ( = norm1(A), since A is symmetric).
339*
340 IF( ifm.EQ.1 ) THEN
341 k = n / 2
342 IF( noe.EQ.1 ) THEN
343* n is odd
344 IF( ilu.EQ.0 ) THEN
345 DO i = 0, k - 1
346 work( i ) = zero
347 END DO
348 DO j = 0, k
349 s = zero
350 DO i = 0, k + j - 1
351 aa = abs( a( i+j*lda ) )
352* -> A(i,j+k)
353 s = s + aa
354 work( i ) = work( i ) + aa
355 END DO
356 aa = abs( a( i+j*lda ) )
357* -> A(j+k,j+k)
358 work( j+k ) = s + aa
359 IF( i.EQ.k+k )
360 $ GO TO 10
361 i = i + 1
362 aa = abs( a( i+j*lda ) )
363* -> A(j,j)
364 work( j ) = work( j ) + aa
365 s = zero
366 DO l = j + 1, k - 1
367 i = i + 1
368 aa = abs( a( i+j*lda ) )
369* -> A(l,j)
370 s = s + aa
371 work( l ) = work( l ) + aa
372 END DO
373 work( j ) = work( j ) + s
374 END DO
375 10 CONTINUE
376 VALUE = work( 0 )
377 DO i = 1, n-1
378 temp = work( i )
379 IF( VALUE .LT. temp .OR. disnan( temp ) )
380 $ VALUE = temp
381 END DO
382 ELSE
383* ilu = 1
384 k = k + 1
385* k=(n+1)/2 for n odd and ilu=1
386 DO i = k, n - 1
387 work( i ) = zero
388 END DO
389 DO j = k - 1, 0, -1
390 s = zero
391 DO i = 0, j - 2
392 aa = abs( a( i+j*lda ) )
393* -> A(j+k,i+k)
394 s = s + aa
395 work( i+k ) = work( i+k ) + aa
396 END DO
397 IF( j.GT.0 ) THEN
398 aa = abs( a( i+j*lda ) )
399* -> A(j+k,j+k)
400 s = s + aa
401 work( i+k ) = work( i+k ) + s
402* i=j
403 i = i + 1
404 END IF
405 aa = abs( a( i+j*lda ) )
406* -> A(j,j)
407 work( j ) = aa
408 s = zero
409 DO l = j + 1, n - 1
410 i = i + 1
411 aa = abs( a( i+j*lda ) )
412* -> A(l,j)
413 s = s + aa
414 work( l ) = work( l ) + aa
415 END DO
416 work( j ) = work( j ) + s
417 END DO
418 VALUE = work( 0 )
419 DO i = 1, n-1
420 temp = work( i )
421 IF( VALUE .LT. temp .OR. disnan( temp ) )
422 $ VALUE = temp
423 END DO
424 END IF
425 ELSE
426* n is even
427 IF( ilu.EQ.0 ) THEN
428 DO i = 0, k - 1
429 work( i ) = zero
430 END DO
431 DO j = 0, k - 1
432 s = zero
433 DO i = 0, k + j - 1
434 aa = abs( a( i+j*lda ) )
435* -> A(i,j+k)
436 s = s + aa
437 work( i ) = work( i ) + aa
438 END DO
439 aa = abs( a( i+j*lda ) )
440* -> A(j+k,j+k)
441 work( j+k ) = s + aa
442 i = i + 1
443 aa = abs( a( i+j*lda ) )
444* -> A(j,j)
445 work( j ) = work( j ) + aa
446 s = zero
447 DO l = j + 1, k - 1
448 i = i + 1
449 aa = abs( a( i+j*lda ) )
450* -> A(l,j)
451 s = s + aa
452 work( l ) = work( l ) + aa
453 END DO
454 work( j ) = work( j ) + s
455 END DO
456 VALUE = work( 0 )
457 DO i = 1, n-1
458 temp = work( i )
459 IF( VALUE .LT. temp .OR. disnan( temp ) )
460 $ VALUE = temp
461 END DO
462 ELSE
463* ilu = 1
464 DO i = k, n - 1
465 work( i ) = zero
466 END DO
467 DO j = k - 1, 0, -1
468 s = zero
469 DO i = 0, j - 1
470 aa = abs( a( i+j*lda ) )
471* -> A(j+k,i+k)
472 s = s + aa
473 work( i+k ) = work( i+k ) + aa
474 END DO
475 aa = abs( a( i+j*lda ) )
476* -> A(j+k,j+k)
477 s = s + aa
478 work( i+k ) = work( i+k ) + s
479* i=j
480 i = i + 1
481 aa = abs( a( i+j*lda ) )
482* -> A(j,j)
483 work( j ) = aa
484 s = zero
485 DO l = j + 1, n - 1
486 i = i + 1
487 aa = abs( a( i+j*lda ) )
488* -> A(l,j)
489 s = s + aa
490 work( l ) = work( l ) + aa
491 END DO
492 work( j ) = work( j ) + s
493 END DO
494 VALUE = work( 0 )
495 DO i = 1, n-1
496 temp = work( i )
497 IF( VALUE .LT. temp .OR. disnan( temp ) )
498 $ VALUE = temp
499 END DO
500 END IF
501 END IF
502 ELSE
503* ifm=0
504 k = n / 2
505 IF( noe.EQ.1 ) THEN
506* n is odd
507 IF( ilu.EQ.0 ) THEN
508 n1 = k
509* n/2
510 k = k + 1
511* k is the row size and lda
512 DO i = n1, n - 1
513 work( i ) = zero
514 END DO
515 DO j = 0, n1 - 1
516 s = zero
517 DO i = 0, k - 1
518 aa = abs( a( i+j*lda ) )
519* A(j,n1+i)
520 work( i+n1 ) = work( i+n1 ) + aa
521 s = s + aa
522 END DO
523 work( j ) = s
524 END DO
525* j=n1=k-1 is special
526 s = abs( a( 0+j*lda ) )
527* A(k-1,k-1)
528 DO i = 1, k - 1
529 aa = abs( a( i+j*lda ) )
530* A(k-1,i+n1)
531 work( i+n1 ) = work( i+n1 ) + aa
532 s = s + aa
533 END DO
534 work( j ) = work( j ) + s
535 DO j = k, n - 1
536 s = zero
537 DO i = 0, j - k - 1
538 aa = abs( a( i+j*lda ) )
539* A(i,j-k)
540 work( i ) = work( i ) + aa
541 s = s + aa
542 END DO
543* i=j-k
544 aa = abs( a( i+j*lda ) )
545* A(j-k,j-k)
546 s = s + aa
547 work( j-k ) = work( j-k ) + s
548 i = i + 1
549 s = abs( a( i+j*lda ) )
550* A(j,j)
551 DO l = j + 1, n - 1
552 i = i + 1
553 aa = abs( a( i+j*lda ) )
554* A(j,l)
555 work( l ) = work( l ) + aa
556 s = s + aa
557 END DO
558 work( j ) = work( j ) + s
559 END DO
560 VALUE = work( 0 )
561 DO i = 1, n-1
562 temp = work( i )
563 IF( VALUE .LT. temp .OR. disnan( temp ) )
564 $ VALUE = temp
565 END DO
566 ELSE
567* ilu=1
568 k = k + 1
569* k=(n+1)/2 for n odd and ilu=1
570 DO i = k, n - 1
571 work( i ) = zero
572 END DO
573 DO j = 0, k - 2
574* process
575 s = zero
576 DO i = 0, j - 1
577 aa = abs( a( i+j*lda ) )
578* A(j,i)
579 work( i ) = work( i ) + aa
580 s = s + aa
581 END DO
582 aa = abs( a( i+j*lda ) )
583* i=j so process of A(j,j)
584 s = s + aa
585 work( j ) = s
586* is initialised here
587 i = i + 1
588* i=j process A(j+k,j+k)
589 aa = abs( a( i+j*lda ) )
590 s = aa
591 DO l = k + j + 1, n - 1
592 i = i + 1
593 aa = abs( a( i+j*lda ) )
594* A(l,k+j)
595 s = s + aa
596 work( l ) = work( l ) + aa
597 END DO
598 work( k+j ) = work( k+j ) + s
599 END DO
600* j=k-1 is special :process col A(k-1,0:k-1)
601 s = zero
602 DO i = 0, k - 2
603 aa = abs( a( i+j*lda ) )
604* A(k,i)
605 work( i ) = work( i ) + aa
606 s = s + aa
607 END DO
608* i=k-1
609 aa = abs( a( i+j*lda ) )
610* A(k-1,k-1)
611 s = s + aa
612 work( i ) = s
613* done with col j=k+1
614 DO j = k, n - 1
615* process col j of A = A(j,0:k-1)
616 s = zero
617 DO i = 0, k - 1
618 aa = abs( a( i+j*lda ) )
619* A(j,i)
620 work( i ) = work( i ) + aa
621 s = s + aa
622 END DO
623 work( j ) = work( j ) + s
624 END DO
625 VALUE = work( 0 )
626 DO i = 1, n-1
627 temp = work( i )
628 IF( VALUE .LT. temp .OR. disnan( temp ) )
629 $ VALUE = temp
630 END DO
631 END IF
632 ELSE
633* n is even
634 IF( ilu.EQ.0 ) THEN
635 DO i = k, n - 1
636 work( i ) = zero
637 END DO
638 DO j = 0, k - 1
639 s = zero
640 DO i = 0, k - 1
641 aa = abs( a( i+j*lda ) )
642* A(j,i+k)
643 work( i+k ) = work( i+k ) + aa
644 s = s + aa
645 END DO
646 work( j ) = s
647 END DO
648* j=k
649 aa = abs( a( 0+j*lda ) )
650* A(k,k)
651 s = aa
652 DO i = 1, k - 1
653 aa = abs( a( i+j*lda ) )
654* A(k,k+i)
655 work( i+k ) = work( i+k ) + aa
656 s = s + aa
657 END DO
658 work( j ) = work( j ) + s
659 DO j = k + 1, n - 1
660 s = zero
661 DO i = 0, j - 2 - k
662 aa = abs( a( i+j*lda ) )
663* A(i,j-k-1)
664 work( i ) = work( i ) + aa
665 s = s + aa
666 END DO
667* i=j-1-k
668 aa = abs( a( i+j*lda ) )
669* A(j-k-1,j-k-1)
670 s = s + aa
671 work( j-k-1 ) = work( j-k-1 ) + s
672 i = i + 1
673 aa = abs( a( i+j*lda ) )
674* A(j,j)
675 s = aa
676 DO l = j + 1, n - 1
677 i = i + 1
678 aa = abs( a( i+j*lda ) )
679* A(j,l)
680 work( l ) = work( l ) + aa
681 s = s + aa
682 END DO
683 work( j ) = work( j ) + s
684 END DO
685* j=n
686 s = zero
687 DO i = 0, k - 2
688 aa = abs( a( i+j*lda ) )
689* A(i,k-1)
690 work( i ) = work( i ) + aa
691 s = s + aa
692 END DO
693* i=k-1
694 aa = abs( a( i+j*lda ) )
695* A(k-1,k-1)
696 s = s + aa
697 work( i ) = work( i ) + s
698 VALUE = work( 0 )
699 DO i = 1, n-1
700 temp = work( i )
701 IF( VALUE .LT. temp .OR. disnan( temp ) )
702 $ VALUE = temp
703 END DO
704 ELSE
705* ilu=1
706 DO i = k, n - 1
707 work( i ) = zero
708 END DO
709* j=0 is special :process col A(k:n-1,k)
710 s = abs( a( 0 ) )
711* A(k,k)
712 DO i = 1, k - 1
713 aa = abs( a( i ) )
714* A(k+i,k)
715 work( i+k ) = work( i+k ) + aa
716 s = s + aa
717 END DO
718 work( k ) = work( k ) + s
719 DO j = 1, k - 1
720* process
721 s = zero
722 DO i = 0, j - 2
723 aa = abs( a( i+j*lda ) )
724* A(j-1,i)
725 work( i ) = work( i ) + aa
726 s = s + aa
727 END DO
728 aa = abs( a( i+j*lda ) )
729* i=j-1 so process of A(j-1,j-1)
730 s = s + aa
731 work( j-1 ) = s
732* is initialised here
733 i = i + 1
734* i=j process A(j+k,j+k)
735 aa = abs( a( i+j*lda ) )
736 s = aa
737 DO l = k + j + 1, n - 1
738 i = i + 1
739 aa = abs( a( i+j*lda ) )
740* A(l,k+j)
741 s = s + aa
742 work( l ) = work( l ) + aa
743 END DO
744 work( k+j ) = work( k+j ) + s
745 END DO
746* j=k is special :process col A(k,0:k-1)
747 s = zero
748 DO i = 0, k - 2
749 aa = abs( a( i+j*lda ) )
750* A(k,i)
751 work( i ) = work( i ) + aa
752 s = s + aa
753 END DO
754* i=k-1
755 aa = abs( a( i+j*lda ) )
756* A(k-1,k-1)
757 s = s + aa
758 work( i ) = s
759* done with col j=k+1
760 DO j = k + 1, n
761* process col j-1 of A = A(j-1,0:k-1)
762 s = zero
763 DO i = 0, k - 1
764 aa = abs( a( i+j*lda ) )
765* A(j-1,i)
766 work( i ) = work( i ) + aa
767 s = s + aa
768 END DO
769 work( j-1 ) = work( j-1 ) + s
770 END DO
771 VALUE = work( 0 )
772 DO i = 1, n-1
773 temp = work( i )
774 IF( VALUE .LT. temp .OR. disnan( temp ) )
775 $ VALUE = temp
776 END DO
777 END IF
778 END IF
779 END IF
780 ELSE IF( ( lsame( norm, 'F' ) ) .OR. ( lsame( norm, 'E' ) ) ) THEN
781*
782* Find normF(A).
783*
784 k = ( n+1 ) / 2
785 scale = zero
786 s = one
787 IF( noe.EQ.1 ) THEN
788* n is odd
789 IF( ifm.EQ.1 ) THEN
790* A is normal
791 IF( ilu.EQ.0 ) THEN
792* A is upper
793 DO j = 0, k - 3
794 CALL dlassq( k-j-2, a( k+j+1+j*lda ), 1, scale, s )
795* L at A(k,0)
796 END DO
797 DO j = 0, k - 1
798 CALL dlassq( k+j-1, a( 0+j*lda ), 1, scale, s )
799* trap U at A(0,0)
800 END DO
801 s = s + s
802* double s for the off diagonal elements
803 CALL dlassq( k-1, a( k ), lda+1, scale, s )
804* tri L at A(k,0)
805 CALL dlassq( k, a( k-1 ), lda+1, scale, s )
806* tri U at A(k-1,0)
807 ELSE
808* ilu=1 & A is lower
809 DO j = 0, k - 1
810 CALL dlassq( n-j-1, a( j+1+j*lda ), 1, scale, s )
811* trap L at A(0,0)
812 END DO
813 DO j = 0, k - 2
814 CALL dlassq( j, a( 0+( 1+j )*lda ), 1, scale, s )
815* U at A(0,1)
816 END DO
817 s = s + s
818* double s for the off diagonal elements
819 CALL dlassq( k, a( 0 ), lda+1, scale, s )
820* tri L at A(0,0)
821 CALL dlassq( k-1, a( 0+lda ), lda+1, scale, s )
822* tri U at A(0,1)
823 END IF
824 ELSE
825* A is xpose
826 IF( ilu.EQ.0 ) THEN
827* A**T is upper
828 DO j = 1, k - 2
829 CALL dlassq( j, a( 0+( k+j )*lda ), 1, scale, s )
830* U at A(0,k)
831 END DO
832 DO j = 0, k - 2
833 CALL dlassq( k, a( 0+j*lda ), 1, scale, s )
834* k by k-1 rect. at A(0,0)
835 END DO
836 DO j = 0, k - 2
837 CALL dlassq( k-j-1, a( j+1+( j+k-1 )*lda ), 1,
838 $ scale, s )
839* L at A(0,k-1)
840 END DO
841 s = s + s
842* double s for the off diagonal elements
843 CALL dlassq( k-1, a( 0+k*lda ), lda+1, scale, s )
844* tri U at A(0,k)
845 CALL dlassq( k, a( 0+( k-1 )*lda ), lda+1, scale, s )
846* tri L at A(0,k-1)
847 ELSE
848* A**T is lower
849 DO j = 1, k - 1
850 CALL dlassq( j, a( 0+j*lda ), 1, scale, s )
851* U at A(0,0)
852 END DO
853 DO j = k, n - 1
854 CALL dlassq( k, a( 0+j*lda ), 1, scale, s )
855* k by k-1 rect. at A(0,k)
856 END DO
857 DO j = 0, k - 3
858 CALL dlassq( k-j-2, a( j+2+j*lda ), 1, scale, s )
859* L at A(1,0)
860 END DO
861 s = s + s
862* double s for the off diagonal elements
863 CALL dlassq( k, a( 0 ), lda+1, scale, s )
864* tri U at A(0,0)
865 CALL dlassq( k-1, a( 1 ), lda+1, scale, s )
866* tri L at A(1,0)
867 END IF
868 END IF
869 ELSE
870* n is even
871 IF( ifm.EQ.1 ) THEN
872* A is normal
873 IF( ilu.EQ.0 ) THEN
874* A is upper
875 DO j = 0, k - 2
876 CALL dlassq( k-j-1, a( k+j+2+j*lda ), 1, scale, s )
877* L at A(k+1,0)
878 END DO
879 DO j = 0, k - 1
880 CALL dlassq( k+j, a( 0+j*lda ), 1, scale, s )
881* trap U at A(0,0)
882 END DO
883 s = s + s
884* double s for the off diagonal elements
885 CALL dlassq( k, a( k+1 ), lda+1, scale, s )
886* tri L at A(k+1,0)
887 CALL dlassq( k, a( k ), lda+1, scale, s )
888* tri U at A(k,0)
889 ELSE
890* ilu=1 & A is lower
891 DO j = 0, k - 1
892 CALL dlassq( n-j-1, a( j+2+j*lda ), 1, scale, s )
893* trap L at A(1,0)
894 END DO
895 DO j = 1, k - 1
896 CALL dlassq( j, a( 0+j*lda ), 1, scale, s )
897* U at A(0,0)
898 END DO
899 s = s + s
900* double s for the off diagonal elements
901 CALL dlassq( k, a( 1 ), lda+1, scale, s )
902* tri L at A(1,0)
903 CALL dlassq( k, a( 0 ), lda+1, scale, s )
904* tri U at A(0,0)
905 END IF
906 ELSE
907* A is xpose
908 IF( ilu.EQ.0 ) THEN
909* A**T is upper
910 DO j = 1, k - 1
911 CALL dlassq( j, a( 0+( k+1+j )*lda ), 1, scale, s )
912* U at A(0,k+1)
913 END DO
914 DO j = 0, k - 1
915 CALL dlassq( k, a( 0+j*lda ), 1, scale, s )
916* k by k rect. at A(0,0)
917 END DO
918 DO j = 0, k - 2
919 CALL dlassq( k-j-1, a( j+1+( j+k )*lda ), 1, scale,
920 $ s )
921* L at A(0,k)
922 END DO
923 s = s + s
924* double s for the off diagonal elements
925 CALL dlassq( k, a( 0+( k+1 )*lda ), lda+1, scale, s )
926* tri U at A(0,k+1)
927 CALL dlassq( k, a( 0+k*lda ), lda+1, scale, s )
928* tri L at A(0,k)
929 ELSE
930* A**T is lower
931 DO j = 1, k - 1
932 CALL dlassq( j, a( 0+( j+1 )*lda ), 1, scale, s )
933* U at A(0,1)
934 END DO
935 DO j = k + 1, n
936 CALL dlassq( k, a( 0+j*lda ), 1, scale, s )
937* k by k rect. at A(0,k+1)
938 END DO
939 DO j = 0, k - 2
940 CALL dlassq( k-j-1, a( j+1+j*lda ), 1, scale, s )
941* L at A(0,0)
942 END DO
943 s = s + s
944* double s for the off diagonal elements
945 CALL dlassq( k, a( lda ), lda+1, scale, s )
946* tri L at A(0,1)
947 CALL dlassq( k, a( 0 ), lda+1, scale, s )
948* tri U at A(0,0)
949 END IF
950 END IF
951 END IF
952 VALUE = scale*sqrt( s )
953 END IF
954*
955 dlansf = VALUE
956 RETURN
957*
958* End of DLANSF
959*
960 END
disnan
logical function disnan(din)
DISNAN tests input for NaN.
Definition disnan.f:59
dlansf
double precision function dlansf(norm, transr, uplo, n, a, work)
DLANSF returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition dlansf.f:209
dlassq
subroutine dlassq(n, x, incx, scale, sumsq)
DLASSQ updates a sum of squares represented in scaled form.
Definition dlassq.f90:124
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48