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