OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
pzlanhs.f
Go to the documentation of this file.
1 DOUBLE PRECISION FUNCTION pzlanhs( NORM, N, A, IA, JA, DESCA,
2 $ WORK )
3*
4* -- ScaLAPACK auxiliary routine (version 1.7) --
5* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
6* and University of California, Berkeley.
7* May 1, 1997
8*
9* .. Scalar Arguments ..
10 CHARACTER norm
11 INTEGER ia, ja, n
12* ..
13* .. Array Arguments ..
14 INTEGER desca( * )
15 DOUBLE PRECISION work( * )
16 COMPLEX*16 a( * )
17* ..
18*
19* Purpose
20* =======
21*
22* PZLANHS returns the value of the one norm, or the Frobenius norm,
23* or the infinity norm, or the element of largest absolute value of a
24* Hessenberg distributed matrix sub( A ) = A(IA:IA+N-1,JA:JA+N-1).
25*
26* PZLANHS returns the value
27*
28* ( max(abs(A(i,j))), NORM = 'M' or 'm' with IA <= i <= IA+N-1,
29* ( and JA <= j <= JA+N-1,
30* (
31* ( norm1( sub( A ) ), NORM = '1', 'O' or 'o'
32* (
33* ( normI( sub( A ) ), NORM = 'I' or 'i'
34* (
35* ( normF( sub( A ) ), NORM = 'F', 'f', 'E' or 'e'
36*
37* where norm1 denotes the one norm of a matrix (maximum column sum),
38* normI denotes the infinity norm of a matrix (maximum row sum) and
39* normF denotes the Frobenius norm of a matrix (square root of sum of
40* squares). Note that max(abs(A(i,j))) is not a matrix norm.
41*
42* Notes
43* =====
44*
45* Each global data object is described by an associated description
46* vector. This vector stores the information required to establish
47* the mapping between an object element and its corresponding process
48* and memory location.
49*
50* Let A be a generic term for any 2D block cyclicly distributed array.
51* Such a global array has an associated description vector DESCA.
52* In the following comments, the character _ should be read as
53* "of the global array".
54*
55* NOTATION STORED IN EXPLANATION
56* --------------- -------------- --------------------------------------
57* DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case,
58* DTYPE_A = 1.
59* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
60* the BLACS process grid A is distribu-
61* ted over. The context itself is glo-
62* bal, but the handle (the integer
63* value) may vary.
64* M_A (global) DESCA( M_ ) The number of rows in the global
65* array A.
66* N_A (global) DESCA( N_ ) The number of columns in the global
67* array A.
68* MB_A (global) DESCA( MB_ ) The blocking factor used to distribute
69* the rows of the array.
70* NB_A (global) DESCA( NB_ ) The blocking factor used to distribute
71* the columns of the array.
72* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
73* row of the array A is distributed.
74* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
75* first column of the array A is
76* distributed.
77* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
78* array. LLD_A >= MAX(1,LOCr(M_A)).
79*
80* Let K be the number of rows or columns of a distributed matrix,
81* and assume that its process grid has dimension p x q.
82* LOCr( K ) denotes the number of elements of K that a process
83* would receive if K were distributed over the p processes of its
84* process column.
85* Similarly, LOCc( K ) denotes the number of elements of K that a
86* process would receive if K were distributed over the q processes of
87* its process row.
88* The values of LOCr() and LOCc() may be determined via a call to the
89* ScaLAPACK tool function, NUMROC:
90* LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
91* LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
92* An upper bound for these quantities may be computed by:
93* LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
94* LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
95*
96* Arguments
97* =========
98*
99* NORM (global input) CHARACTER
100* Specifies the value to be returned in PZLANHS as described
101* above.
102*
103* N (global input) INTEGER
104* The number of rows and columns to be operated on i.e the
105* number of rows and columns of the distributed submatrix
106* sub( A ). When N = 0, PZLANHS is set to zero. N >= 0.
107*
108* A (local input) COMPLEX*16 pointer into the local memory
109* to an array of dimension (LLD_A, LOCc(JA+N-1) ) containing
110* the local pieces of sub( A ).
111*
112* IA (global input) INTEGER
113* The row index in the global array A indicating the first
114* row of sub( A ).
115*
116* JA (global input) INTEGER
117* The column index in the global array A indicating the
118* first column of sub( A ).
119*
120* DESCA (global and local input) INTEGER array of dimension DLEN_.
121* The array descriptor for the distributed matrix A.
122*
123* WORK (local workspace) DOUBLE PRECISION array dimension (LWORK)
124* LWORK >= 0 if NORM = 'M' or 'm' (not referenced),
125* Nq0 if NORM = '1', 'O' or 'o',
126* Mp0 if NORM = 'I' or 'i',
127* 0 if NORM = 'F', 'f', 'E' or 'e' (not referenced),
128* where
129*
130* IROFFA = MOD( IA-1, MB_A ), ICOFFA = MOD( JA-1, NB_A ),
131* IAROW = INDXG2P( IA, MB_A, MYROW, RSRC_A, NPROW ),
132* IACOL = INDXG2P( JA, NB_A, MYCOL, CSRC_A, NPCOL ),
133* Np0 = NUMROC( N+IROFFA, MB_A, MYROW, IAROW, NPROW ),
134* Nq0 = NUMROC( N+ICOFFA, NB_A, MYCOL, IACOL, NPCOL ),
135*
136* INDXG2P and NUMROC are ScaLAPACK tool functions; MYROW,
137* MYCOL, NPROW and NPCOL can be determined by calling the
138* subroutine BLACS_GRIDINFO.
139*
140* =====================================================================
141*
142* .. Parameters ..
143 INTEGER block_cyclic_2d, csrc_, ctxt_, dlen_, dtype_,
144 $ lld_, mb_, m_, nb_, n_, rsrc_
145 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
146 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
147 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
148 DOUBLE PRECISION one, zero
149 parameter( one = 1.0d+0, zero = 0.0d+0 )
150* ..
151* .. Local Scalars ..
152 INTEGER iacol, iarow, ictxt, ii, iia, icoff, inxtrow,
153 $ ioffa, iroff, j, jb, jj, jja, jn, kk, lda, ll,
154 $ mycol, myrow, np, npcol, nprow, nq
155 DOUBLE PRECISION scale, sum, value
156* ..
157* .. Local Arrays ..
158 DOUBLE PRECISION rwork( 2 )
159* ..
160* .. External Subroutines ..
161 EXTERNAL blacs_gridinfo, dcombssq, dgebr2d,
162 $ dgebs2d, dgamx2d, dgsum2d,
163 $ infog2l, pdtreecomb, zlassq
164* ..
165* .. External Functions ..
166 LOGICAL lsame
167 INTEGER iceil, IDAMAX, NUMROC
168 EXTERNAL lsame, iceil, idamax, numroc
169* ..
170* .. Intrinsic Functions ..
171 INTRINSIC abs, max, min, mod, sqrt
172* ..
173* .. Executable Statements ..
174*
175* Get grid parameters.
176*
177 ictxt = desca( ctxt_ )
178 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
179*
180 CALL infog2l( ia, ja, desca, nprow, npcol, myrow, mycol, iia, jja,
181 $ iarow, iacol )
182 iroff = mod( ia-1, desca( mb_ ) )
183 icoff = mod( ja-1, desca( nb_ ) )
184 np = numroc( n+iroff, desca( mb_ ), myrow, iarow, nprow )
185 nq = numroc( n+icoff, desca( nb_ ), mycol, iacol, npcol )
186 IF( myrow.EQ.iarow )
187 $ np = np - iroff
188 IF( mycol.EQ.iacol )
189 $ nq = nq - icoff
190 lda = desca( lld_ )
191 ioffa = ( jja - 1 ) * lda
192*
193 IF( n.EQ.0 ) THEN
194*
195 VALUE = zero
196*
197 ELSE IF( lsame( norm, 'M' ) ) THEN
198*
199 VALUE = zero
200*
201* Find max(abs(A(i,j))).
202*
203 ii = iia
204 jj = jja
205 jn = min( iceil( ja, desca( nb_ ) ) * desca( nb_ ), ja+n-1 )
206 jb = jn-ja+1
207*
208* Only one process row
209*
210 IF( nprow.EQ.1 ) THEN
211*
212* Handle first block of columns separately
213*
214 IF( mycol.EQ.iacol ) THEN
215 DO 20 ll = jj, jj+jb-1
216 DO 10 kk = iia, min( ii+ll-jj+1, iia+np-1 )
217 VALUE = max( VALUE, abs( a( ioffa+kk ) ) )
218 10 CONTINUE
219 ioffa = ioffa + lda
220 20 CONTINUE
221 jj = jj + jb
222 END IF
223*
224 iacol = mod( iacol+1, npcol )
225*
226* Loop over remaining block of columns
227*
228 DO 50 j = jn+1, ja+n-1, desca( nb_ )
229 jb = min( ja+n-j, desca( nb_ ) )
230*
231 IF( mycol.EQ.iacol ) THEN
232 DO 40 ll = jj, jj+jb-1
233 DO 30 kk = iia, min( ii+ll-jj+1, iia+np-1 )
234 VALUE = max( VALUE, abs( a( ioffa+kk ) ) )
235 30 CONTINUE
236 ioffa = ioffa + lda
237 40 CONTINUE
238 jj = jj + jb
239 END IF
240*
241 ii = ii + jb
242 iacol = mod( iacol+1, npcol )
243*
244 50 CONTINUE
245*
246 ELSE
247*
248* Handle first block of columns separately
249*
250 inxtrow = mod( iarow+1, nprow )
251 IF( mycol.EQ.iacol ) THEN
252 IF( myrow.EQ.iarow ) THEN
253 DO 70 ll = jj, jj + jb -1
254 DO 60 kk = iia, min( ii+ll-jj+1, iia+np-1 )
255 VALUE = max( VALUE, abs( a( ioffa+kk ) ) )
256 60 CONTINUE
257 ioffa = ioffa + lda
258 70 CONTINUE
259 ELSE
260 DO 90 ll = jj, jj+jb-1
261 DO 80 kk = iia, min( ii-1, iia+np-1 )
262 VALUE = max( VALUE, abs( a( ioffa+kk ) ) )
263 80 CONTINUE
264 ioffa = ioffa + lda
265 90 CONTINUE
266 IF( myrow.EQ.inxtrow .AND. ii.LE.iia+np-1 )
267 $ VALUE = max( VALUE, abs( a( ii+(jj+jb-2)*lda ) ) )
268 END IF
269 jj = jj + jb
270 END IF
271*
272 IF( myrow.EQ.iarow )
273 $ ii = ii + jb
274 iarow = inxtrow
275 iarow = mod( iarow+1, nprow )
276 iacol = mod( iacol+1, npcol )
277*
278* Loop over remaining block of columns
279*
280 DO 140 j = jn+1, ja+n-1, desca( nb_ )
281 jb = min( ja+n-j, desca( nb_ ) )
282*
283 IF( mycol.EQ.iacol ) THEN
284 IF( myrow.EQ.iarow ) THEN
285 DO 110 ll = jj, jj + jb -1
286 DO 100 kk = iia, min( ii+ll-jj+1, iia+np-1 )
287 VALUE = max( VALUE, abs( a( ioffa+kk ) ) )
288 100 CONTINUE
289 ioffa = ioffa + lda
290 110 CONTINUE
291 ELSE
292 DO 130 ll = jj, jj + jb -1
293 DO 120 kk = iia, min( ii-1, iia+np-1 )
294 VALUE = max( VALUE, abs( a( ioffa+kk ) ) )
295 120 CONTINUE
296 ioffa = ioffa + lda
297 130 CONTINUE
298 IF( myrow.EQ.inxtrow .AND. ii.LE.iia+np-1 )
299 $ VALUE = max( VALUE,
300 $ abs( a( ii+(jj+jb-2)*lda ) ) )
301 END IF
302 jj = jj + jb
303 END IF
304*
305 IF( myrow.EQ.iarow )
306 $ ii = ii + jb
307 iarow = inxtrow
308 iarow = mod( iarow+1, nprow )
309 iacol = mod( iacol+1, npcol )
310*
311 140 CONTINUE
312*
313 END IF
314*
315* Gather the intermediate results to process (0,0).
316*
317 CALL dgamx2d( ictxt, 'All', ' ', 1, 1, VALUE, 1, kk, ll, -1,
318 $ 0, 0 )
319*
320 ELSE IF( lsame( norm, 'O' ) .OR. norm.EQ.'1' ) THEN
321*
322 VALUE = zero
323 ii = iia
324 jj = jja
325 jn = min( iceil( ja, desca( nb_ ) ) * desca( nb_ ), ja+n-1 )
326 jb = jn-ja+1
327*
328* Only one process row
329*
330 IF( nprow.EQ.1 ) THEN
331*
332* Handle first block of columns separately
333*
334 IF( mycol.EQ.iacol ) THEN
335 DO 160 ll = jj, jj+jb-1
336 sum = zero
337 DO 150 kk = iia, min( ii+ll-jj+1, iia+np-1 )
338 sum = sum + abs( a( ioffa+kk ) )
339 150 CONTINUE
340 ioffa = ioffa + lda
341 work( ll-jja+1 ) = sum
342 160 CONTINUE
343 jj = jj + jb
344 END IF
345*
346 iacol = mod( iacol+1, npcol )
347*
348* Loop over remaining block of columns
349*
350 DO 190 j = jn+1, ja+n-1, desca( nb_ )
351 jb = min( ja+n-j, desca( nb_ ) )
352*
353 IF( mycol.EQ.iacol ) THEN
354 DO 180 ll = jj, jj+jb-1
355 sum = zero
356 DO 170 kk = iia, min( ii+ll-jj+1, iia+np-1 )
357 sum = sum + abs( a( ioffa+kk ) )
358 170 CONTINUE
359 ioffa = ioffa + lda
360 work( ll-jja+1 ) = sum
361 180 CONTINUE
362 jj = jj + jb
363 END IF
364*
365 ii = ii + jb
366 iacol = mod( iacol+1, npcol )
367*
368 190 CONTINUE
369*
370 ELSE
371*
372* Handle first block of columns separately
373*
374 inxtrow = mod( iarow+1, nprow )
375 IF( mycol.EQ.iacol ) THEN
376 IF( myrow.EQ.iarow ) THEN
377 DO 210 ll = jj, jj + jb -1
378 sum = zero
379 DO 200 kk = iia, min( ii+ll-jj+1, iia+np-1 )
380 sum = sum + abs( a( ioffa+kk ) )
381 200 CONTINUE
382 ioffa = ioffa + lda
383 work( ll-jja+1 ) = sum
384 210 CONTINUE
385 ELSE
386 DO 230 ll = jj, jj + jb -1
387 sum = zero
388 DO 220 kk = iia, min( ii-1, iia+np-1 )
389 sum = sum + abs( a( ioffa+kk ) )
390 220 CONTINUE
391 ioffa = ioffa + lda
392 work( ll-jja+1 ) = sum
393 230 CONTINUE
394 IF( myrow.EQ.inxtrow .AND. ii.LE.iia+np-1 )
395 $ work( jj+jb-jja ) = work( jj+jb-jja ) +
396 $ abs( a( ii+(jj+jb-2)*lda ) )
397 END IF
398 jj = jj + jb
399 END IF
400*
401 IF( myrow.EQ.iarow )
402 $ ii = ii + jb
403 iarow = inxtrow
404 iarow = mod( iarow+1, nprow )
405 iacol = mod( iacol+1, npcol )
406*
407* Loop over remaining block of columns
408*
409 DO 280 j = jn+1, ja+n-1, desca( nb_ )
410 jb = min( ja+n-j, desca( nb_ ) )
411*
412 IF( mycol.EQ.iacol ) THEN
413 IF( myrow.EQ.iarow ) THEN
414 DO 250 ll = jj, jj + jb -1
415 sum = zero
416 DO 240 kk = iia, min( ii+ll-jj+1, iia+np-1 )
417 sum = sum + abs( a( ioffa+kk ) )
418 240 CONTINUE
419 ioffa = ioffa + lda
420 work( ll-jja+1 ) = sum
421 250 CONTINUE
422 ELSE
423 DO 270 ll = jj, jj + jb -1
424 sum = zero
425 DO 260 kk = iia, min( ii-1, iia+np-1 )
426 sum = sum + abs( a( ioffa+kk ) )
427 260 CONTINUE
428 ioffa = ioffa + lda
429 work( ll-jja+1 ) = sum
430 270 CONTINUE
431 IF( myrow.EQ.inxtrow .AND. ii.LE.iia+np-1 )
432 $ work( jj+jb-jja ) = work( jj+jb-jja ) +
433 $ abs( a( ii+(jj+jb-2)*lda ) )
434 END IF
435 jj = jj + jb
436 END IF
437*
438 IF( myrow.EQ.iarow )
439 $ ii = ii + jb
440 iarow = inxtrow
441 iarow = mod( iarow+1, nprow )
442 iacol = mod( iacol+1, npcol )
443*
444 280 CONTINUE
445*
446 END IF
447*
448* Find sum of global matrix columns and store on row 0 of
449* process grid
450*
451 CALL dgsum2d( ictxt, 'Columnwise', ' ', 1, nq, work, 1,
452 $ 0, mycol )
453*
454* Find maximum sum of columns for 1-norm
455*
456 IF( myrow.EQ.0 ) THEN
457 IF( nq.GT.0 ) THEN
458 VALUE = work( idamax( nq, work, 1 ) )
459 ELSE
460 VALUE = zero
461 END IF
462 CALL dgamx2d( ictxt, 'Rowwise', ' ', 1, 1, VALUE, 1, KK, LL,
463 $ -1, 0, 0 )
464 END IF
465*
466 ELSE IF( LSAME( NORM, 'i' ) ) THEN
467*
468 DO 290 KK = IIA, IIA+NP-1
469 WORK( KK ) = ZERO
470 290 CONTINUE
471*
472 II = IIA
473 JJ = JJA
474 JN = MIN( ICEIL( JA, DESCA( NB_ ) ) * DESCA( NB_ ), JA+N-1 )
475 JB = JN-JA+1
476*
477* Only one process row
478*
479.EQ. IF( NPROW1 ) THEN
480*
481* Handle first block of columns separately
482*
483.EQ. IF( MYCOLIACOL ) THEN
484 DO 310 LL = JJ, JJ+JB-1
485 DO 300 KK = IIA, MIN( II+LL-JJ+1, IIA+NP-1 )
486 WORK( KK-IIA+1 ) = WORK( KK-IIA+1 ) +
487 $ ABS( A( IOFFA+KK ) )
488 300 CONTINUE
489 IOFFA = IOFFA + LDA
490 310 CONTINUE
491 JJ = JJ + JB
492 END IF
493*
494 IACOL = MOD( IACOL+1, NPCOL )
495*
496* Loop over remaining block of columns
497*
498 DO 340 J = JN+1, JA+N-1, DESCA( NB_ )
499 JB = MIN( JA+N-J, DESCA( NB_ ) )
500*
501.EQ. IF( MYCOLIACOL ) THEN
502 DO 330 LL = JJ, JJ+JB-1
503 DO 320 KK = IIA, MIN( II+LL-JJ+1, IIA+NP-1 )
504 WORK( KK-IIA+1 ) = WORK( KK-IIA+1 ) +
505 $ ABS( A( IOFFA+KK ) )
506 320 CONTINUE
507 IOFFA = IOFFA + LDA
508 330 CONTINUE
509 JJ = JJ + JB
510 END IF
511*
512 II = II + JB
513 IACOL = MOD( IACOL+1, NPCOL )
514*
515 340 CONTINUE
516*
517 ELSE
518*
519* Handle first block of columns separately
520*
521 INXTROW = MOD( IAROW+1, NPROW )
522.EQ. IF( MYCOLIACOL ) THEN
523.EQ. IF( MYROWIAROW ) THEN
524 DO 360 LL = JJ, JJ + JB -1
525 DO 350 KK = IIA, MIN( II+LL-JJ+1, IIA+NP-1 )
526 WORK( KK-IIA+1 ) = WORK( KK-IIA+1 ) +
527 $ ABS( A( IOFFA+KK ) )
528 350 CONTINUE
529 IOFFA = IOFFA + LDA
530 360 CONTINUE
531 ELSE
532 DO 380 LL = JJ, JJ + JB -1
533 DO 370 KK = IIA, MIN( II-1, IIA+NP-1 )
534 WORK( KK-IIA+1 ) = WORK( KK-IIA+1 ) +
535 $ ABS( A( IOFFA+KK ) )
536 370 CONTINUE
537 IOFFA = IOFFA + LDA
538 380 CONTINUE
539.EQ..AND..LE. IF( MYROWINXTROW IIIIA+NP-1 )
540 $ WORK( II-IIA+1 ) = WORK( II-IIA+1 ) +
541 $ ABS( A( II+(JJ+JB-2)*LDA ) )
542 END IF
543 JJ = JJ + JB
544 END IF
545*
546.EQ. IF( MYROWIAROW )
547 $ II = II + JB
548 IAROW = INXTROW
549 IAROW = MOD( IAROW+1, NPROW )
550 IACOL = MOD( IACOL+1, NPCOL )
551*
552* Loop over remaining block of columns
553*
554 DO 430 J = JN+1, JA+N-1, DESCA( NB_ )
555 JB = MIN( JA+N-J, DESCA( NB_ ) )
556*
557.EQ. IF( MYCOLIACOL ) THEN
558.EQ. IF( MYROWIAROW ) THEN
559 DO 400 LL = JJ, JJ + JB -1
560 DO 390 KK = IIA, MIN( II+LL-JJ+1, IIA+NP-1 )
561 WORK( KK-IIA+1 ) = WORK( KK-IIA+1 ) +
562 $ ABS( A( IOFFA+KK ) )
563 390 CONTINUE
564 IOFFA = IOFFA + LDA
565 400 CONTINUE
566 ELSE
567 DO 420 LL = JJ, JJ + JB -1
568 DO 410 KK = IIA, MIN( II-1, IIA+NP-1 )
569 WORK( KK-IIA+1 ) = WORK( KK-IIA+1 ) +
570 $ ABS(A(IOFFA+KK))
571 410 CONTINUE
572 IOFFA = IOFFA + LDA
573 420 CONTINUE
574.EQ..AND..LE. IF( MYROWINXTROW IIIIA+NP-1 )
575 $ WORK( II-IIA+1 ) = WORK( II-IIA+1 ) +
576 $ ABS( A( II+(JJ+JB-2)*LDA ) )
577 END IF
578 JJ = JJ + JB
579 END IF
580*
581.EQ. IF( MYROWIAROW )
582 $ II = II + JB
583 IAROW = INXTROW
584 IAROW = MOD( IAROW+1, NPROW )
585 IACOL = MOD( IACOL+1, NPCOL )
586*
587 430 CONTINUE
588*
589 END IF
590*
591* Find sum of global matrix rows and store on column 0 of
592* process grid
593*
594 CALL DGSUM2D( ICTXT, 'rowwise', ' ', NP, 1, WORK, MAX( 1, NP ),
595 $ MYROW, 0 )
596*
597* Find maximum sum of rows for Infinity-norm
598*
599.EQ. IF( MYCOL0 ) THEN
600.GT. IF( NP0 ) THEN
601 VALUE = WORK( IDAMAX( NP, WORK, 1 ) )
602 ELSE
603 VALUE = ZERO
604 END IF
605 CALL DGAMX2D( ICTXT, 'columnwise', ' ', 1, 1, VALUE, 1, KK,
606 $ LL, -1, 0, 0 )
607 END IF
608*
609 ELSE IF( LSAME( NORM, 'f.OR.' ) LSAME( NORM, 'e' ) ) THEN
610*
611 SCALE = ZERO
612 SUM = ONE
613 II = IIA
614 JJ = JJA
615 JN = MIN( ICEIL( JA, DESCA( NB_ ) ) * DESCA( NB_ ), JA+N-1 )
616 JB = JN-JA+1
617*
618* Only one process row
619*
620.EQ. IF( NPROW1 ) THEN
621*
622* Handle first block of columns separately
623*
624.EQ. IF( MYCOLIACOL ) THEN
625 DO 440 LL = JJ, JJ+JB-1
626 CALL ZLASSQ( MIN( II+LL-JJ+1, IIA+NP-1 )-IIA+1,
627 $ A( IIA+IOFFA ), 1, SCALE, SUM )
628 IOFFA = IOFFA + LDA
629 440 CONTINUE
630 JJ = JJ + JB
631 END IF
632*
633 IACOL = MOD( IACOL+1, NPCOL )
634*
635* Loop over remaining block of columns
636*
637 DO 460 J = JN+1, JA+N-1, DESCA( NB_ )
638 JB = MIN( JA+N-J, DESCA( NB_ ) )
639*
640.EQ. IF( MYCOLIACOL ) THEN
641 DO 450 LL = JJ, JJ+JB-1
642 CALL ZLASSQ( MIN( II+LL-JJ+1, IIA+NP-1 )-IIA+1,
643 $ A( IIA+IOFFA ), 1, SCALE, SUM )
644 IOFFA = IOFFA + LDA
645 450 CONTINUE
646 JJ = JJ + JB
647 END IF
648*
649 II = II + JB
650 IACOL = MOD( IACOL+1, NPCOL )
651*
652 460 CONTINUE
653*
654 ELSE
655*
656* Handle first block of columns separately
657*
658 INXTROW = MOD( IAROW+1, NPROW )
659.EQ. IF( MYCOLIACOL ) THEN
660.EQ. IF( MYROWIAROW ) THEN
661 DO 470 LL = JJ, JJ + JB -1
662 CALL ZLASSQ( MIN( II+LL-JJ+1, IIA+NP-1 )-IIA+1,
663 $ A( IIA+IOFFA ), 1, SCALE, SUM )
664 IOFFA = IOFFA + LDA
665 470 CONTINUE
666 ELSE
667 DO 480 LL = JJ, JJ + JB -1
668 CALL ZLASSQ( MIN( II-1, IIA+NP-1 )-IIA+1,
669 $ A( IIA+IOFFA ), 1, SCALE, SUM )
670 IOFFA = IOFFA + LDA
671 480 CONTINUE
672.EQ..AND..LE. IF( MYROWINXTROW IIIIA+NP-1 )
673 $ CALL ZLASSQ( 1, A( II+(JJ+JB-2)*LDA ), 1,
674 $ SCALE, SUM )
675 END IF
676 JJ = JJ + JB
677 END IF
678*
679.EQ. IF( MYROWIAROW )
680 $ II = II + JB
681 IAROW = INXTROW
682 IAROW = MOD( IAROW+1, NPROW )
683 IACOL = MOD( IACOL+1, NPCOL )
684*
685* Loop over remaining block of columns
686*
687 DO 510 J = JN+1, JA+N-1, DESCA( NB_ )
688 JB = MIN( JA+N-J, DESCA( NB_ ) )
689*
690.EQ. IF( MYCOLIACOL ) THEN
691.EQ. IF( MYROWIAROW ) THEN
692 DO 490 LL = JJ, JJ + JB -1
693 CALL ZLASSQ( MIN( II+LL-JJ+1, IIA+NP-1 )-IIA+1,
694 $ A( IIA+IOFFA ), 1, SCALE, SUM )
695 IOFFA = IOFFA + LDA
696 490 CONTINUE
697 ELSE
698 DO 500 LL = JJ, JJ + JB -1
699 CALL ZLASSQ( MIN( II-1, IIA+NP-1 )-IIA+1,
700 $ A( IIA+IOFFA ), 1, SCALE, SUM )
701 IOFFA = IOFFA + LDA
702 500 CONTINUE
703.EQ..AND..LE. IF( MYROWINXTROW IIIIA+NP-1 )
704 $ CALL ZLASSQ( 1, A( II+(JJ+JB-2)*LDA ), 1,
705 $ SCALE, SUM )
706 END IF
707 JJ = JJ + JB
708 END IF
709*
710.EQ. IF( MYROWIAROW )
711 $ II = II + JB
712 IAROW = INXTROW
713 IAROW = MOD( IAROW+1, NPROW )
714 IACOL = MOD( IACOL+1, NPCOL )
715*
716 510 CONTINUE
717*
718 END IF
719*
720* Perform the global scaled sum
721*
722 RWORK( 1 ) = SCALE
723 RWORK( 2 ) = SUM
724 CALL PDTREECOMB( ICTXT, 'all', 2, RWORK, 0, 0, DCOMBSSQ )
725 VALUE = RWORK( 1 ) * SQRT( RWORK( 2 ) )
726*
727 END IF
728*
729.EQ..AND..EQ. IF( MYROW0 MYCOL0 ) THEN
730 CALL DGEBS2D( ICTXT, 'all', ' ', 1, 1, VALUE, 1 )
731 ELSE
732 CALL DGEBR2D( ICTXT, 'all', ' ', 1, 1, VALUE, 1, 0, 0 )
733 END IF
734*
735 PZLANHS = VALUE
736*
737 RETURN
738*
739* End of PZLANHS
740*
741 END
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
zlassq
subroutine zlassq(n, x, incx, scl, sumsq)
ZLASSQ updates a sum of squares represented in scaled form.
Definition zlassq.f90:137
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:53
idamax
integer function idamax(n, dx, incx)
IDAMAX
Definition idamax.f:71
iceil
integer function iceil(inum, idenom)
Definition iceil.f:2
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
dgebs2d
subroutine dgebs2d(contxt, scope, top, m, n, a, lda)
Definition mpi.f:1082
dgebr2d
subroutine dgebr2d(contxt, scope, top, m, n, a, lda)
Definition mpi.f:1123
infog2l
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)
Definition mpi.f:937
blacs_gridinfo
subroutine blacs_gridinfo(cntxt, nprow, npcol, myrow, mycol)
Definition mpi.f:754
numroc
integer function numroc(n, nb, iproc, isrcproc, nprocs)
Definition mpi.f:786
dcombssq
subroutine dcombssq(v1, v2)
Definition pdtreecomb.f:259
pdtreecomb
subroutine pdtreecomb(ictxt, scope, n, mine, rdest0, cdest0, subptr)
Definition pdtreecomb.f:3
pzlanhs
double precision function pzlanhs(norm, n, a, ia, ja, desca, work)
Definition pzlanhs.f:3