OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
pdblastim.f
Go to the documentation of this file.
1 SUBROUTINE pdlascal( TYPE, M, N, ALPHA, A, IA, JA, DESCA )
2*
3* -- PBLAS test routine (version 2.0) --
4* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
5* and University of California, Berkeley.
6* April 1, 1998
7*
8* .. Scalar Arguments ..
9 CHARACTER*1 TYPE
10 INTEGER IA, JA, M, N
11 DOUBLE PRECISION ALPHA
12* ..
13* .. Array Arguments ..
14 INTEGER DESCA( * )
15 DOUBLE PRECISION A( * )
16* ..
17*
18* Purpose
19* =======
20*
21* PDLASCAL scales the m by n submatrix A(IA:IA+M-1,JA:JA+N-1) denoted
22* by sub( A ) by the scalar alpha. TYPE specifies if sub( A ) is full,
23* upper triangular, lower triangular or upper Hessenberg.
24*
25* Notes
26* =====
27*
28* A description vector is associated with each 2D block-cyclicly dis-
29* tributed matrix. This vector stores the information required to
30* establish the mapping between a matrix entry and its corresponding
31* process and memory location.
32*
33* In the following comments, the character _ should be read as
34* "of the distributed matrix". Let A be a generic term for any 2D
35* block cyclicly distributed matrix. Its description vector is DESCA:
36*
37* NOTATION STORED IN EXPLANATION
38* ---------------- --------------- ------------------------------------
39* DTYPE_A (global) DESCA( DTYPE_ ) The descriptor type.
40* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
41* the NPROW x NPCOL BLACS process grid
42* A is distributed over. The context
43* itself is global, but the handle
44* (the integer value) may vary.
45* M_A (global) DESCA( M_ ) The number of rows in the distribu-
46* ted matrix A, M_A >= 0.
47* N_A (global) DESCA( N_ ) The number of columns in the distri-
48* buted matrix A, N_A >= 0.
49* IMB_A (global) DESCA( IMB_ ) The number of rows of the upper left
50* block of the matrix A, IMB_A > 0.
51* INB_A (global) DESCA( INB_ ) The number of columns of the upper
52* left block of the matrix A,
53* INB_A > 0.
54* MB_A (global) DESCA( MB_ ) The blocking factor used to distri-
55* bute the last M_A-IMB_A rows of A,
56* MB_A > 0.
57* NB_A (global) DESCA( NB_ ) The blocking factor used to distri-
58* bute the last N_A-INB_A columns of
59* A, NB_A > 0.
60* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
61* row of the matrix A is distributed,
62* NPROW > RSRC_A >= 0.
63* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
64* first column of A is distributed.
65* NPCOL > CSRC_A >= 0.
66* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
67* array storing the local blocks of
68* the distributed matrix A,
69* IF( Lc( 1, N_A ) > 0 )
70* LLD_A >= MAX( 1, Lr( 1, M_A ) )
71* ELSE
72* LLD_A >= 1.
73*
74* Let K be the number of rows of a matrix A starting at the global in-
75* dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
76* that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
77* receive if these K rows were distributed over NPROW processes. If K
78* is the number of columns of a matrix A starting at the global index
79* JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
80* lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
81* these K columns were distributed over NPCOL processes.
82*
83* The values of Lr() and Lc() may be determined via a call to the func-
84* tion PB_NUMROC:
85* Lr( IA, K ) = PB_NUMROC( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
86* Lc( JA, K ) = PB_NUMROC( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
87*
88* Arguments
89* =========
90*
91* TYPE (global input) CHARACTER*1
92* On entry, TYPE specifies the type of the input submatrix as
93* follows:
94* = 'L' or 'l': sub( A ) is a lower triangular matrix,
95* = 'U' or 'u': sub( A ) is an upper triangular matrix,
96* = 'H' or 'h': sub( A ) is an upper Hessenberg matrix,
97* otherwise sub( A ) is a full matrix.
98*
99* M (global input) INTEGER
100* On entry, M specifies the number of rows of the submatrix
101* sub( A ). M must be at least zero.
102*
103* N (global input) INTEGER
104* On entry, N specifies the number of columns of the submatrix
105* sub( A ). N must be at least zero.
106*
107* ALPHA (global input) DOUBLE PRECISION
108* On entry, ALPHA specifies the scalar alpha.
109*
110* A (local input/local output) DOUBLE PRECISION array
111* On entry, A is an array of dimension (LLD_A, Ka), where Ka is
112* at least Lc( 1, JA+N-1 ). Before entry, this array contains
113* the local entries of the matrix A.
114* On exit, the local entries of this array corresponding to the
115* to the entries of the submatrix sub( A ) are overwritten by
116* the local entries of the m by n scaled submatrix.
117*
118* IA (global input) INTEGER
119* On entry, IA specifies A's global row index, which points to
120* the beginning of the submatrix sub( A ).
121*
122* JA (global input) INTEGER
123* On entry, JA specifies A's global column index, which points
124* to the beginning of the submatrix sub( A ).
125*
126* DESCA (global and local input) INTEGER array
127* On entry, DESCA is an integer array of dimension DLEN_. This
128* is the array descriptor for the matrix A.
129*
130* -- Written on April 1, 1998 by
131* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
132*
133* =====================================================================
134*
135* .. Parameters ..
136 INTEGER BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
137 $ DTYPE_, IMB_, INB_, LLD_, MB_, M_, NB_, N_,
138 $ RSRC_
139 parameter( block_cyclic_2d_inb = 2, dlen_ = 11,
140 $ dtype_ = 1, ctxt_ = 2, m_ = 3, n_ = 4,
141 $ imb_ = 5, inb_ = 6, mb_ = 7, nb_ = 8,
142 $ rsrc_ = 9, csrc_ = 10, lld_ = 11 )
143* ..
144* .. Local Scalars ..
145 CHARACTER*1 UPLO
146 LOGICAL GODOWN, GOLEFT, LOWER, UPPER
147 INTEGER IACOL, IAROW, ICTXT, IIA, IIMAX, ILOW, IMB1,
148 $ IMBLOC, INB1, INBLOC, IOFFA, IOFFD, ITYPE,
149 $ IUPP, JJA, JJMAX, JOFFA, JOFFD, LCMT, LCMT00,
150 $ LDA, LMBLOC, LNBLOC, LOW, M1, MB, MBLKD, MBLKS,
151 $ MBLOC, MP, MRCOL, MRROW, MYCOL, MYROW, N1, NB,
152 $ NBLKD, NBLKS, NBLOC, NPCOL, NPROW, NQ, PMB,
153 $ QNB, TMP1, UPP
154* ..
155* .. Local Arrays ..
156 INTEGER DESCA2( DLEN_ )
157* ..
158* .. External Subroutines ..
161* ..
162* .. External Functions ..
163 LOGICAL LSAME
164 INTEGER PB_NUMROC
165 EXTERNAL lsame, pb_numroc
166* ..
167* .. Intrinsic Functions ..
168 INTRINSIC min
169* ..
170* .. Executable Statements ..
171*
172* Convert descriptor
173*
174 CALL pb_desctrans( desca, desca2 )
175*
176* Get grid parameters
177*
178 ictxt = desca2( ctxt_ )
179 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
180*
181* Quick return if possible
182*
183 IF( m.EQ.0 .OR. n.EQ.0 )
184 $ RETURN
185*
186 IF( lsame( TYPE, 'L' ) ) then
187 itype = 1
188 uplo = TYPE
189 upper = .false.
190 lower = .true.
191 ioffd = 0
192 ELSE IF( lsame( TYPE, 'U' ) ) then
193 itype = 2
194 uplo = TYPE
195 upper = .true.
196 lower = .false.
197 ioffd = 0
198 ELSE IF( lsame( TYPE, 'H' ) ) then
199 itype = 3
200 uplo = 'U'
201 upper = .true.
202 lower = .false.
203 ioffd = 1
204 ELSE
205 itype = 0
206 uplo = 'A'
207 upper = .true.
208 lower = .true.
209 ioffd = 0
210 END IF
211*
212* Compute local indexes
213*
214 IF( itype.EQ.0 ) THEN
215*
216* Full matrix
217*
218 CALL pb_infog2l( ia, ja, desca2, nprow, npcol, myrow, mycol,
219 $ iia, jja, iarow, iacol )
220 mp = pb_numroc( m, ia, desca2( imb_ ), desca2( mb_ ), myrow,
221 $ desca2( rsrc_ ), nprow )
222 nq = pb_numroc( n, ja, desca2( inb_ ), desca2( nb_ ), mycol,
223 $ desca2( csrc_ ), npcol )
224*
225 IF( mp.LE.0 .OR. nq.LE.0 )
226 $ RETURN
227*
228 lda = desca2( lld_ )
229 ioffa = iia + ( jja - 1 ) * lda
230*
231 CALL pb_dlascal( 'All', mp, nq, 0, alpha, a( ioffa ), lda )
232*
233 ELSE
234*
235* Trapezoidal matrix
236*
237 CALL pb_ainfog2l( m, n, ia, ja, desca2, nprow, npcol, myrow,
238 $ mycol, imb1, inb1, mp, nq, iia, jja, iarow,
239 $ iacol, mrrow, mrcol )
240*
241 IF( mp.LE.0 .OR. nq.LE.0 )
242 $ RETURN
243*
244* Initialize LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC, LMBLOC,
245* LNBLOC, ILOW, LOW, IUPP, and UPP.
246*
247 mb = desca2( mb_ )
248 nb = desca2( nb_ )
249 lda = desca2( lld_ )
250*
251 CALL pb_binfo( ioffd, mp, nq, imb1, inb1, mb, nb, mrrow,
252 $ mrcol, lcmt00, mblks, nblks, imbloc, inbloc,
253 $ lmbloc, lnbloc, ilow, low, iupp, upp )
254*
255 m1 = mp
256 n1 = nq
257 ioffa = iia - 1
258 joffa = jja - 1
259 iimax = ioffa + mp
260 jjmax = joffa + nq
261*
262 IF( desca2( rsrc_ ).LT.0 ) THEN
263 pmb = mb
264 ELSE
265 pmb = nprow * mb
266 END IF
267 IF( desca2( csrc_ ).LT.0 ) THEN
268 qnb = nb
269 ELSE
270 qnb = npcol * nb
271 END IF
272*
273* Handle the first block of rows or columns separately, and
274* update LCMT00, MBLKS and NBLKS.
275*
276 godown = ( lcmt00.GT.iupp )
277 goleft = ( lcmt00.LT.ilow )
278*
279 IF( .NOT.godown .AND. .NOT.goleft ) THEN
280*
281* LCMT00 >= ILOW && LCMT00 <= IUPP
282*
283 goleft = ( ( lcmt00 - ( iupp - upp + pmb ) ).LT.ilow )
284 godown = .NOT.goleft
285*
286 CALL pb_dlascal( uplo, imbloc, inbloc, lcmt00, alpha,
287 $ a( iia+joffa*lda ), lda )
288 IF( godown ) THEN
289 IF( upper .AND. nq.GT.inbloc )
290 $ CALL pb_dlascal( 'All', imbloc, nq-inbloc, 0, alpha,
291 $ a( iia+(joffa+inbloc)*lda ), lda )
292 iia = iia + imbloc
293 m1 = m1 - imbloc
294 ELSE
295 IF( lower .AND. mp.GT.imbloc )
296 $ CALL pb_dlascal( 'All', mp-imbloc, inbloc, 0, alpha,
297 $ a( iia+imbloc+joffa*lda ), lda )
298 jja = jja + inbloc
299 n1 = n1 - inbloc
300 END IF
301*
302 END IF
303*
304 IF( godown ) THEN
305*
306 lcmt00 = lcmt00 - ( iupp - upp + pmb )
307 mblks = mblks - 1
308 ioffa = ioffa + imbloc
309*
310 10 CONTINUE
311 IF( mblks.GT.0 .AND. lcmt00.GT.upp ) THEN
312 lcmt00 = lcmt00 - pmb
313 mblks = mblks - 1
314 ioffa = ioffa + mb
315 GO TO 10
316 END IF
317*
318 tmp1 = min( ioffa, iimax ) - iia + 1
319 IF( upper .AND. tmp1.GT.0 ) THEN
320 CALL pb_dlascal( 'All', tmp1, n1, 0, alpha,
321 $ a( iia+joffa*lda ), lda )
322 iia = iia + tmp1
323 m1 = m1 - tmp1
324 END IF
325*
326 IF( mblks.LE.0 )
327 $ RETURN
328*
329 lcmt = lcmt00
330 mblkd = mblks
331 ioffd = ioffa
332*
333 mbloc = mb
334 20 CONTINUE
335 IF( mblkd.GT.0 .AND. lcmt.GE.ilow ) THEN
336 IF( mblkd.EQ.1 )
337 $ mbloc = lmbloc
338 CALL pb_dlascal( uplo, mbloc, inbloc, lcmt, alpha,
339 $ a( ioffd+1+joffa*lda ), lda )
340 lcmt00 = lcmt
341 lcmt = lcmt - pmb
342 mblks = mblkd
343 mblkd = mblkd - 1
344 ioffa = ioffd
345 ioffd = ioffd + mbloc
346 GO TO 20
347 END IF
348*
349 tmp1 = m1 - ioffd + iia - 1
350 IF( lower .AND. tmp1.GT.0 )
351 $ CALL pb_dlascal( 'All', tmp1, inbloc, 0, alpha,
352 $ a( ioffd+1+joffa*lda ), lda )
353*
354 tmp1 = ioffa - iia + 1
355 m1 = m1 - tmp1
356 n1 = n1 - inbloc
357 lcmt00 = lcmt00 + low - ilow + qnb
358 nblks = nblks - 1
359 joffa = joffa + inbloc
360*
361 IF( upper .AND. tmp1.GT.0 .AND. n1.GT.0 )
362 $ CALL pb_dlascal( 'All', tmp1, n1, 0, alpha,
363 $ a( iia+joffa*lda ), lda )
364*
365 iia = ioffa + 1
366 jja = joffa + 1
367*
368 ELSE IF( goleft ) THEN
369*
370 lcmt00 = lcmt00 + low - ilow + qnb
371 nblks = nblks - 1
372 joffa = joffa + inbloc
373*
374 30 CONTINUE
375 IF( nblks.GT.0 .AND. lcmt00.LT.low ) THEN
376 lcmt00 = lcmt00 + qnb
377 nblks = nblks - 1
378 joffa = joffa + nb
379 GO TO 30
380 END IF
381*
382 tmp1 = min( joffa, jjmax ) - jja + 1
383 IF( lower .AND. tmp1.GT.0 ) THEN
384 CALL pb_dlascal( 'All', m1, tmp1, 0, alpha,
385 $ a( iia+(jja-1)*lda ), lda )
386 jja = jja + tmp1
387 n1 = n1 - tmp1
388 END IF
389*
390 IF( nblks.LE.0 )
391 $ RETURN
392*
393 lcmt = lcmt00
394 nblkd = nblks
395 joffd = joffa
396*
397 nbloc = nb
398 40 CONTINUE
399 IF( nblkd.GT.0 .AND. lcmt.LE.iupp ) THEN
400 IF( nblkd.EQ.1 )
401 $ nbloc = lnbloc
402 CALL pb_dlascal( uplo, imbloc, nbloc, lcmt, alpha,
403 $ a( iia+joffd*lda ), lda )
404 lcmt00 = lcmt
405 lcmt = lcmt + qnb
406 nblks = nblkd
407 nblkd = nblkd - 1
408 joffa = joffd
409 joffd = joffd + nbloc
410 GO TO 40
411 END IF
412*
413 tmp1 = n1 - joffd + jja - 1
414 IF( upper .AND. tmp1.GT.0 )
415 $ CALL pb_dlascal( 'All', imbloc, tmp1, 0, alpha,
416 $ a( iia+joffd*lda ), lda )
417*
418 tmp1 = joffa - jja + 1
419 m1 = m1 - imbloc
420 n1 = n1 - tmp1
421 lcmt00 = lcmt00 - ( iupp - upp + pmb )
422 mblks = mblks - 1
423 ioffa = ioffa + imbloc
424*
425 IF( lower .AND. m1.GT.0 .AND. tmp1.GT.0 )
426 $ CALL pb_dlascal( 'All', m1, tmp1, 0, alpha,
427 $ a( ioffa+1+(jja-1)*lda ), lda )
428*
429 iia = ioffa + 1
430 jja = joffa + 1
431*
432 END IF
433*
434 nbloc = nb
435 50 CONTINUE
436 IF( nblks.GT.0 ) THEN
437 IF( nblks.EQ.1 )
438 $ nbloc = lnbloc
439 60 CONTINUE
440 IF( mblks.GT.0 .AND. lcmt00.GT.upp ) THEN
441 lcmt00 = lcmt00 - pmb
442 mblks = mblks - 1
443 ioffa = ioffa + mb
444 GO TO 60
445 END IF
446*
447 tmp1 = min( ioffa, iimax ) - iia + 1
448 IF( upper .AND. tmp1.GT.0 ) THEN
449 CALL pb_dlascal( 'all', TMP1, N1, 0, ALPHA,
450 $ A( IIA+JOFFA*LDA ), LDA )
451 IIA = IIA + TMP1
452 M1 = M1 - TMP1
453 END IF
454*
455.LE. IF( MBLKS0 )
456 $ RETURN
457*
458 LCMT = LCMT00
459 MBLKD = MBLKS
460 IOFFD = IOFFA
461*
462 MBLOC = MB
463 70 CONTINUE
464.GT..AND..GE. IF( MBLKD0 LCMTLOW ) THEN
465.EQ. IF( MBLKD1 )
466 $ MBLOC = LMBLOC
467 CALL PB_DLASCAL( UPLO, MBLOC, NBLOC, LCMT, ALPHA,
468 $ A( IOFFD+1+JOFFA*LDA ), LDA )
469 LCMT00 = LCMT
470 LCMT = LCMT - PMB
471 MBLKS = MBLKD
472 MBLKD = MBLKD - 1
473 IOFFA = IOFFD
474 IOFFD = IOFFD + MBLOC
475 GO TO 70
476 END IF
477*
478 TMP1 = M1 - IOFFD + IIA - 1
479.AND..GT. IF( LOWER TMP10 )
480 $ CALL PB_DLASCAL( 'all', TMP1, NBLOC, 0, ALPHA,
481 $ A( IOFFD+1+JOFFA*LDA ), LDA )
482*
483 TMP1 = MIN( IOFFA, IIMAX ) - IIA + 1
484 M1 = M1 - TMP1
485 N1 = N1 - NBLOC
486 LCMT00 = LCMT00 + QNB
487 NBLKS = NBLKS - 1
488 JOFFA = JOFFA + NBLOC
489*
490.AND..GT..AND..GT. IF( UPPER TMP10 N10 )
491 $ CALL PB_DLASCAL( 'all', TMP1, N1, 0, ALPHA,
492 $ A( IIA+JOFFA*LDA ), LDA )
493*
494 IIA = IOFFA + 1
495 JJA = JOFFA + 1
496*
497 GO TO 50
498*
499 END IF
500*
501 END IF
502*
503 RETURN
504*
505* End of PDLASCAL
506*
507 END
508 SUBROUTINE PDLAGEN( INPLACE, AFORM, DIAG, OFFA, M, N, IA, JA,
509 $ DESCA, IASEED, A, LDA )
510*
511* -- PBLAS test routine (version 2.0) --
512* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
513* and University of California, Berkeley.
514* April 1, 1998
515*
516* .. Scalar Arguments ..
517 LOGICAL INPLACE
518 CHARACTER*1 AFORM, DIAG
519 INTEGER IA, IASEED, JA, LDA, M, N, OFFA
520* ..
521* .. Array Arguments ..
522 INTEGER DESCA( * )
523 DOUBLE PRECISION A( LDA, * )
524* ..
525*
526* Purpose
527* =======
528*
529* PDLAGEN generates (or regenerates) a submatrix sub( A ) denoting
530* A(IA:IA+M-1,JA:JA+N-1).
531*
532* Notes
533* =====
534*
535* A description vector is associated with each 2D block-cyclicly dis-
536* tributed matrix. This vector stores the information required to
537* establish the mapping between a matrix entry and its corresponding
538* process and memory location.
539*
540* In the following comments, the character _ should be read as
541* "of the distributed matrix". Let A be a generic term for any 2D
542* block cyclicly distributed matrix. Its description vector is DESCA:
543*
544* NOTATION STORED IN EXPLANATION
545* ---------------- --------------- ------------------------------------
546* DTYPE_A (global) DESCA( DTYPE_ ) The descriptor type.
547* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
548* the NPROW x NPCOL BLACS process grid
549* A is distributed over. The context
550* itself is global, but the handle
551* (the integer value) may vary.
552* M_A (global) DESCA( M_ ) The number of rows in the distribu-
553* ted matrix A, M_A >= 0.
554* N_A (global) DESCA( N_ ) The number of columns in the distri-
555* buted matrix A, N_A >= 0.
556* IMB_A (global) DESCA( IMB_ ) The number of rows of the upper left
557* block of the matrix A, IMB_A > 0.
558* INB_A (global) DESCA( INB_ ) The number of columns of the upper
559* left block of the matrix A,
560* INB_A > 0.
561* MB_A (global) DESCA( MB_ ) The blocking factor used to distri-
562* bute the last M_A-IMB_A rows of A,
563* MB_A > 0.
564* NB_A (global) DESCA( NB_ ) The blocking factor used to distri-
565* bute the last N_A-INB_A columns of
566* A, NB_A > 0.
567* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
568* row of the matrix A is distributed,
569* NPROW > RSRC_A >= 0.
570* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
571* first column of A is distributed.
572* NPCOL > CSRC_A >= 0.
573* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
574* array storing the local blocks of
575* the distributed matrix A,
576* IF( Lc( 1, N_A ) > 0 )
577* LLD_A >= MAX( 1, Lr( 1, M_A ) )
578* ELSE
579* LLD_A >= 1.
580*
581* Let K be the number of rows of a matrix A starting at the global in-
582* dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
583* that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
584* receive if these K rows were distributed over NPROW processes. If K
585* is the number of columns of a matrix A starting at the global index
586* JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
587* lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
588* these K columns were distributed over NPCOL processes.
589*
590* The values of Lr() and Lc() may be determined via a call to the func-
591* tion PB_NUMROC:
592* Lr( IA, K ) = PB_NUMROC( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
593* Lc( JA, K ) = PB_NUMROC( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
594*
595* Arguments
596* =========
597*
598* INPLACE (global input) LOGICAL
599* On entry, INPLACE specifies if the matrix should be generated
600* in place or not. If INPLACE is .TRUE., the local random array
601* to be generated will start in memory at the local memory lo-
602* cation A( 1, 1 ), otherwise it will start at the local posi-
603* tion induced by IA and JA.
604*
605* AFORM (global input) CHARACTER*1
606* On entry, AFORM specifies the type of submatrix to be genera-
607* ted as follows:
608* AFORM = 'S', sub( A ) is a symmetric matrix,
609* AFORM = 'H', sub( A ) is a Hermitian matrix,
610* AFORM = 'T', sub( A ) is overrwritten with the transpose
611* of what would normally be generated,
612* AFORM = 'C', sub( A ) is overwritten with the conjugate
613* transpose of what would normally be genera-
614* ted.
615* AFORM = 'N', a random submatrix is generated.
616*
617* DIAG (global input) CHARACTER*1
618* On entry, DIAG specifies if the generated submatrix is diago-
619* nally dominant or not as follows:
620* DIAG = 'D' : sub( A ) is diagonally dominant,
621* DIAG = 'N' : sub( A ) is not diagonally dominant.
622*
623* OFFA (global input) INTEGER
624* On entry, OFFA specifies the offdiagonal of the underlying
625* matrix A(1:DESCA(M_),1:DESCA(N_)) of interest when the subma-
626* trix is symmetric, Hermitian or diagonally dominant. OFFA = 0
627* specifies the main diagonal, OFFA > 0 specifies a subdiago-
628* nal, and OFFA < 0 specifies a superdiagonal (see further de-
629* tails).
630*
631* M (global input) INTEGER
632* On entry, M specifies the global number of matrix rows of the
633* submatrix sub( A ) to be generated. M must be at least zero.
634*
635* N (global input) INTEGER
636* On entry, N specifies the global number of matrix columns of
637* the submatrix sub( A ) to be generated. N must be at least
638* zero.
639*
640* IA (global input) INTEGER
641* On entry, IA specifies A's global row index, which points to
642* the beginning of the submatrix sub( A ).
643*
644* JA (global input) INTEGER
645* On entry, JA specifies A's global column index, which points
646* to the beginning of the submatrix sub( A ).
647*
648* DESCA (global and local input) INTEGER array
649* On entry, DESCA is an integer array of dimension DLEN_. This
650* is the array descriptor for the matrix A.
651*
652* IASEED (global input) INTEGER
653* On entry, IASEED specifies the seed number to generate the
654* matrix A. IASEED must be at least zero.
655*
656* A (local output) DOUBLE PRECISION array
657* On entry, A is an array of dimension (LLD_A, Ka), where Ka is
658* at least Lc( 1, JA+N-1 ). On exit, this array contains the
659* local entries of the randomly generated submatrix sub( A ).
660*
661* LDA (local input) INTEGER
662* On entry, LDA specifies the local leading dimension of the
663* array A. When INPLACE is .FALSE., LDA is usually DESCA(LLD_).
664* This restriction is however not enforced, and this subroutine
665* requires only that LDA >= MAX( 1, Mp ) where
666*
667* Mp = PB_NUMROC( M, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW ).
668*
669* PB_NUMROC is a ScaLAPACK tool function; MYROW, MYCOL, NPROW
670* and NPCOL can be determined by calling the BLACS subroutine
671* BLACS_GRIDINFO.
672*
673* Further Details
674* ===============
675*
676* OFFD is tied to the matrix described by DESCA, as opposed to the
677* piece that is currently (re)generated. This is a global information
678* independent from the distribution parameters. Below are examples of
679* the meaning of OFFD for a global 7 by 5 matrix:
680*
681* ---------------------------------------------------------------------
682* OFFD | 0 -1 -2 -3 -4 0 -1 -2 -3 -4 0 -1 -2 -3 -4
683* -------|-------------------------------------------------------------
684* | | OFFD=-1 | OFFD=0 OFFD=2
685* | V V
686* 0 | . d . . . -> d . . . . . . . . .
687* 1 | . . d . . . d . . . . . . . .
688* 2 | . . . d . . . d . . -> d . . . .
689* 3 | . . . . d . . . d . . d . . .
690* 4 | . . . . . . . . . d . . d . .
691* 5 | . . . . . . . . . . . . . d .
692* 6 | . . . . . . . . . . . . . . d
693* ---------------------------------------------------------------------
694*
695* -- Written on April 1, 1998 by
696* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
697*
698* =====================================================================
699*
700* .. Parameters ..
701 INTEGER BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
702 $ DTYPE_, IMB_, INB_, LLD_, MB_, M_, NB_, N_,
703 $ RSRC_
704 PARAMETER ( BLOCK_CYCLIC_2D_INB = 2, DLEN_ = 11,
705 $ DTYPE_ = 1, CTXT_ = 2, M_ = 3, N_ = 4,
706 $ IMB_ = 5, INB_ = 6, MB_ = 7, NB_ = 8,
707 $ RSRC_ = 9, CSRC_ = 10, LLD_ = 11 )
708 INTEGER JMP_1, JMP_COL, JMP_IMBV, JMP_INBV, JMP_LEN,
709 $ JMP_MB, JMP_NB, JMP_NPIMBLOC, JMP_NPMB,
710 $ JMP_NQINBLOC, JMP_NQNB, JMP_ROW
711 PARAMETER ( JMP_1 = 1, JMP_ROW = 2, JMP_COL = 3,
712 $ JMP_MB = 4, JMP_IMBV = 5, JMP_NPMB = 6,
713 $ JMP_NPIMBLOC = 7, JMP_NB = 8, JMP_INBV = 9,
714 $ JMP_NQNB = 10, JMP_NQINBLOC = 11,
715 $ JMP_LEN = 11 )
716* ..
717* .. Local Scalars ..
718 LOGICAL DIAGDO, SYMM, HERM, NOTRAN
719 INTEGER CSRC, I, IACOL, IAROW, ICTXT, IIA, ILOCBLK,
720 $ ILOCOFF, ILOW, IMB, IMB1, IMBLOC, IMBVIR, INB,
721 $ INB1, INBLOC, INBVIR, INFO, IOFFDA, ITMP, IUPP,
722 $ IVIR, JJA, JLOCBLK, JLOCOFF, JVIR, LCMT00,
723 $ LMBLOC, LNBLOC, LOW, MAXMN, MB, MBLKS, MP,
724 $ MRCOL, MRROW, MYCDIST, MYCOL, MYRDIST, MYROW,
725 $ NB, NBLKS, NPCOL, NPROW, NQ, NVIR, RSRC, UPP
726 DOUBLE PRECISION ALPHA
727* ..
728* .. Local Arrays ..
729 INTEGER DESCA2( DLEN_ ), IMULADD( 4, JMP_LEN ),
730 $ IRAN( 2 ), JMP( JMP_LEN ), MULADD0( 4 )
731* ..
732* .. External Subroutines ..
733 EXTERNAL BLACS_GRIDINFO, PB_AINFOG2L, PB_BINFO,
734 $ PB_CHKMAT, PB_DESCTRANS, PB_DLAGEN, PB_INITJMP,
735 $ PB_INITMULADD, PB_JUMP, PB_JUMPIT, PB_LOCINFO,
736 $ PB_SETLOCRAN, PB_SETRAN, PDLADOM, PXERBLA
737* ..
738* .. External Functions ..
739 LOGICAL LSAME
740 EXTERNAL LSAME
741* ..
742* .. Intrinsic Functions ..
743 INTRINSIC DBLE, MAX, MIN
744* ..
745* .. Data Statements ..
746 DATA ( MULADD0( I ), I = 1, 4 ) / 20077, 16838,
747 $ 12345, 0 /
748* ..
749* .. Executable Statements ..
750*
751* Convert descriptor
752*
753 CALL PB_DESCTRANS( DESCA, DESCA2 )
754*
755* Test the input arguments
756*
757 ICTXT = DESCA2( CTXT_ )
758 CALL BLACS_GRIDINFO( ICTXT, NPROW, NPCOL, MYROW, MYCOL )
759*
760* Test the input parameters
761*
762 INFO = 0
763.EQ. IF( NPROW-1 ) THEN
764 INFO = -( 1000 + CTXT_ )
765 ELSE
766 SYMM = LSAME( AFORM, 's' )
767 HERM = LSAME( AFORM, 'h' )
768 NOTRAN = LSAME( AFORM, 'n' )
769 DIAGDO = LSAME( DIAG, 'd' )
770.NOT..OR..OR..AND. IF( ( SYMMHERMNOTRAN )
771.NOT. $ ( LSAME( AFORM, 't.AND.' ) )
772.NOT. $ ( LSAME( AFORM, 'c' ) ) ) THEN
773 INFO = -2
774.NOT..AND. ELSE IF( ( DIAGDO )
775.NOT. $ ( LSAME( DIAG, 'n' ) ) ) THEN
776 INFO = -3
777 END IF
778 CALL PB_CHKMAT( ICTXT, M, 5, N, 6, IA, JA, DESCA2, 10, INFO )
779 END IF
780*
781.NE. IF( INFO0 ) THEN
782 CALL PXERBLA( ICTXT, 'pdlagen', -INFO )
783 RETURN
784 END IF
785*
786* Quick return if possible
787*
788.LE..OR..LE. IF( ( M0 )( N0 ) )
789 $ RETURN
790*
791* Start the operations
792*
793 MB = DESCA2( MB_ )
794 NB = DESCA2( NB_ )
795 IMB = DESCA2( IMB_ )
796 INB = DESCA2( INB_ )
797 RSRC = DESCA2( RSRC_ )
798 CSRC = DESCA2( CSRC_ )
799*
800* Figure out local information about the distributed matrix operand
801*
802 CALL PB_AINFOG2L( M, N, IA, JA, DESCA2, NPROW, NPCOL, MYROW,
803 $ MYCOL, IMB1, INB1, MP, NQ, IIA, JJA, IAROW,
804 $ IACOL, MRROW, MRCOL )
805*
806* Decide where the entries shall be stored in memory
807*
808 IF( INPLACE ) THEN
809 IIA = 1
810 JJA = 1
811 END IF
812*
813* Initialize LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC, LMBLOC, LNBLOC,
814* ILOW, LOW, IUPP, and UPP.
815*
816 IOFFDA = JA + OFFA - IA
817 CALL PB_BINFO( IOFFDA, MP, NQ, IMB1, INB1, MB, NB, MRROW,
818 $ MRCOL, LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC,
819 $ LMBLOC, LNBLOC, ILOW, LOW, IUPP, UPP )
820*
821* Initialize ILOCBLK, ILOCOFF, MYRDIST, JLOCBLK, JLOCOFF, MYCDIST
822* This values correspond to the square virtual underlying matrix
823* of size MAX( M_ + MAX( 0, -OFFA ), N_ + MAX( 0, OFFA ) ) used
824* to set up the random sequence. For practical purposes, the size
825* of this virtual matrix is upper bounded by M_ + N_ - 1.
826*
827 ITMP = MAX( 0, -OFFA )
828 IVIR = IA + ITMP
829 IMBVIR = IMB + ITMP
830 NVIR = DESCA2( M_ ) + ITMP
831*
832 CALL PB_LOCINFO( IVIR, IMBVIR, MB, MYROW, RSRC, NPROW, ILOCBLK,
833 $ ILOCOFF, MYRDIST )
834*
835 ITMP = MAX( 0, OFFA )
836 JVIR = JA + ITMP
837 INBVIR = INB + ITMP
838 NVIR = MAX( MAX( NVIR, DESCA2( N_ ) + ITMP ),
839 $ DESCA2( M_ ) + DESCA2( N_ ) - 1 )
840*
841 CALL PB_LOCINFO( JVIR, INBVIR, NB, MYCOL, CSRC, NPCOL, JLOCBLK,
842 $ JLOCOFF, MYCDIST )
843*
844.OR..OR. IF( SYMM HERM NOTRAN ) THEN
845*
846 CALL PB_INITJMP( .TRUE., NVIR, IMBVIR, INBVIR, IMBLOC, INBLOC,
847 $ MB, NB, RSRC, CSRC, NPROW, NPCOL, 1, JMP )
848*
849* Compute constants to jump JMP( * ) numbers in the sequence
850*
851 CALL PB_INITMULADD( MULADD0, JMP, IMULADD )
852*
853* Compute and set the random value corresponding to A( IA, JA )
854*
855 CALL PB_SETLOCRAN( IASEED, ILOCBLK, JLOCBLK, ILOCOFF, JLOCOFF,
856 $ MYRDIST, MYCDIST, NPROW, NPCOL, JMP,
857 $ IMULADD, IRAN )
858*
859 CALL PB_DLAGEN( 'lower', AFORM, A( IIA, JJA ), LDA, LCMT00,
860 $ IRAN, MBLKS, IMBLOC, MB, LMBLOC, NBLKS, INBLOC,
861 $ NB, LNBLOC, JMP, IMULADD )
862*
863 END IF
864*
865.OR..OR..NOT. IF( SYMM HERM ( NOTRAN ) ) THEN
866*
867 CALL PB_INITJMP( .FALSE., NVIR, IMBVIR, INBVIR, IMBLOC, INBLOC,
868 $ MB, NB, RSRC, CSRC, NPROW, NPCOL, 1, JMP )
869*
870* Compute constants to jump JMP( * ) numbers in the sequence
871*
872 CALL PB_INITMULADD( MULADD0, JMP, IMULADD )
873*
874* Compute and set the random value corresponding to A( IA, JA )
875*
876 CALL PB_SETLOCRAN( IASEED, ILOCBLK, JLOCBLK, ILOCOFF, JLOCOFF,
877 $ MYRDIST, MYCDIST, NPROW, NPCOL, JMP,
878 $ IMULADD, IRAN )
879*
880 CALL PB_DLAGEN( 'upper', AFORM, A( IIA, JJA ), LDA, LCMT00,
881 $ IRAN, MBLKS, IMBLOC, MB, LMBLOC, NBLKS, INBLOC,
882 $ NB, LNBLOC, JMP, IMULADD )
883*
884 END IF
885*
886 IF( DIAGDO ) THEN
887*
888 MAXMN = MAX( DESCA2( M_ ), DESCA2( N_ ) )
889 ALPHA = DBLE( MAXMN )
890*
891.GE. IF( IOFFDA0 ) THEN
892 CALL PDLADOM( INPLACE, MIN( MAX( 0, M-IOFFDA ), N ), ALPHA,
893 $ A, MIN( IA+IOFFDA, IA+M-1 ), JA, DESCA )
894 ELSE
895 CALL PDLADOM( INPLACE, MIN( M, MAX( 0, N+IOFFDA ) ), ALPHA,
896 $ A, IA, MIN( JA-IOFFDA, JA+N-1 ), DESCA )
897 END IF
898*
899 END IF
900*
901 RETURN
902*
903* End of PDLAGEN
904*
905 END
906 SUBROUTINE PDLADOM( INPLACE, N, ALPHA, A, IA, JA, DESCA )
907*
908* -- PBLAS test routine (version 2.0) --
909* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
910* and University of California, Berkeley.
911* April 1, 1998
912*
913* .. Scalar Arguments ..
914 LOGICAL INPLACE
915 INTEGER IA, JA, N
916 DOUBLE PRECISION ALPHA
917* ..
918* .. Array Arguments ..
919 INTEGER DESCA( * )
920 DOUBLE PRECISION A( * )
921* ..
922*
923* Purpose
924* =======
925*
926* PDLADOM adds alpha to the diagonal entries of an n by n submatrix
927* sub( A ) denoting A( IA:IA+N-1, JA:JA+N-1 ).
928*
929* Notes
930* =====
931*
932* A description vector is associated with each 2D block-cyclicly dis-
933* tributed matrix. This vector stores the information required to
934* establish the mapping between a matrix entry and its corresponding
935* process and memory location.
936*
937* In the following comments, the character _ should be read as
938* "of the distributed matrix". Let A be a generic term for any 2D
939* block cyclicly distributed matrix. Its description vector is DESCA:
940*
941* NOTATION STORED IN EXPLANATION
942* ---------------- --------------- ------------------------------------
943* DTYPE_A (global) DESCA( DTYPE_ ) The descriptor type.
944* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
945* the NPROW x NPCOL BLACS process grid
946* A is distributed over. The context
947* itself is global, but the handle
948* (the integer value) may vary.
949* M_A (global) DESCA( M_ ) The number of rows in the distribu-
950* ted matrix A, M_A >= 0.
951* N_A (global) DESCA( N_ ) The number of columns in the distri-
952* buted matrix A, N_A >= 0.
953* IMB_A (global) DESCA( IMB_ ) The number of rows of the upper left
954* block of the matrix A, IMB_A > 0.
955* INB_A (global) DESCA( INB_ ) The number of columns of the upper
956* left block of the matrix A,
957* INB_A > 0.
958* MB_A (global) DESCA( MB_ ) The blocking factor used to distri-
959* bute the last M_A-IMB_A rows of A,
960* MB_A > 0.
961* NB_A (global) DESCA( NB_ ) The blocking factor used to distri-
962* bute the last N_A-INB_A columns of
963* A, NB_A > 0.
964* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
965* row of the matrix A is distributed,
966* NPROW > RSRC_A >= 0.
967* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
968* first column of A is distributed.
969* NPCOL > CSRC_A >= 0.
970* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
971* array storing the local blocks of
972* the distributed matrix A,
973* IF( Lc( 1, N_A ) > 0 )
974* LLD_A >= MAX( 1, Lr( 1, M_A ) )
975* ELSE
976* LLD_A >= 1.
977*
978* Let K be the number of rows of a matrix A starting at the global in-
979* dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
980* that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
981* receive if these K rows were distributed over NPROW processes. If K
982* is the number of columns of a matrix A starting at the global index
983* JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
984* lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
985* these K columns were distributed over NPCOL processes.
986*
987* The values of Lr() and Lc() may be determined via a call to the func-
988* tion PB_NUMROC:
989* Lr( IA, K ) = PB_NUMROC( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
990* Lc( JA, K ) = PB_NUMROC( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
991*
992* Arguments
993* =========
994*
995* INPLACE (global input) LOGICAL
996* On entry, INPLACE specifies if the matrix should be generated
997* in place or not. If INPLACE is .TRUE., the local random array
998* to be generated will start in memory at the local memory lo-
999* cation A( 1, 1 ), otherwise it will start at the local posi-
1000* tion induced by IA and JA.
1001*
1002* N (global input) INTEGER
1003* On entry, N specifies the global order of the submatrix
1004* sub( A ) to be modified. N must be at least zero.
1005*
1006* ALPHA (global input) DOUBLE PRECISION
1007* On entry, ALPHA specifies the scalar alpha.
1008*
1009* A (local input/local output) DOUBLE PRECISION array
1010* On entry, A is an array of dimension (LLD_A, Ka), where Ka is
1011* at least Lc( 1, JA+N-1 ). Before entry, this array contains
1012* the local entries of the matrix A. On exit, the local entries
1013* of this array corresponding to the main diagonal of sub( A )
1014* have been updated.
1015*
1016* IA (global input) INTEGER
1017* On entry, IA specifies A's global row index, which points to
1018* the beginning of the submatrix sub( A ).
1019*
1020* JA (global input) INTEGER
1021* On entry, JA specifies A's global column index, which points
1022* to the beginning of the submatrix sub( A ).
1023*
1024* DESCA (global and local input) INTEGER array
1025* On entry, DESCA is an integer array of dimension DLEN_. This
1026* is the array descriptor for the matrix A.
1027*
1028* -- Written on April 1, 1998 by
1029* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
1030*
1031* =====================================================================
1032*
1033* .. Parameters ..
1034 INTEGER BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
1035 $ DTYPE_, IMB_, INB_, LLD_, MB_, M_, NB_, N_,
1036 $ RSRC_
1037 PARAMETER ( BLOCK_CYCLIC_2D_INB = 2, DLEN_ = 11,
1038 $ DTYPE_ = 1, CTXT_ = 2, M_ = 3, N_ = 4,
1039 $ IMB_ = 5, INB_ = 6, MB_ = 7, NB_ = 8,
1040 $ RSRC_ = 9, CSRC_ = 10, LLD_ = 11 )
1041* ..
1042* .. Local Scalars ..
1043 LOGICAL GODOWN, GOLEFT
1044 INTEGER I, IACOL, IAROW, ICTXT, IIA, IJOFFA, ILOW,
1045 $ IMB1, IMBLOC, INB1, INBLOC, IOFFA, IOFFD, IUPP,
1046 $ JJA, JOFFA, JOFFD, LCMT, LCMT00, LDA, LDAP1,
1047 $ LMBLOC, LNBLOC, LOW, MB, MBLKD, MBLKS, MBLOC,
1048 $ MRCOL, MRROW, MYCOL, MYROW, NB, NBLKD, NBLKS,
1049 $ NBLOC, NP, NPCOL, NPROW, NQ, PMB, QNB, UPP
1050 DOUBLE PRECISION ATMP
1051* ..
1052* .. Local Scalars ..
1053 INTEGER DESCA2( DLEN_ )
1054* ..
1055* .. External Subroutines ..
1056 EXTERNAL BLACS_GRIDINFO, PB_AINFOG2L, PB_BINFO,
1057 $ PB_DESCTRANS
1058* ..
1059* .. Intrinsic Functions ..
1060 INTRINSIC ABS, MAX, MIN
1061* ..
1062* .. Executable Statements ..
1063*
1064* Convert descriptor
1065*
1066 CALL PB_DESCTRANS( DESCA, DESCA2 )
1067*
1068* Get grid parameters
1069*
1070 ICTXT = DESCA2( CTXT_ )
1071 CALL BLACS_GRIDINFO( ICTXT, NPROW, NPCOL, MYROW, MYCOL )
1072*
1073.EQ. IF( N0 )
1074 $ RETURN
1075*
1076 CALL PB_AINFOG2L( N, N, IA, JA, DESCA2, NPROW, NPCOL, MYROW,
1077 $ MYCOL, IMB1, INB1, NP, NQ, IIA, JJA, IAROW,
1078 $ IACOL, MRROW, MRCOL )
1079*
1080* Decide where the entries shall be stored in memory
1081*
1082 IF( INPLACE ) THEN
1083 IIA = 1
1084 JJA = 1
1085 END IF
1086*
1087* Initialize LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC, LMBLOC, LNBLOC,
1088* ILOW, LOW, IUPP, and UPP.
1089*
1090 MB = DESCA2( MB_ )
1091 NB = DESCA2( NB_ )
1092*
1093 CALL PB_BINFO( 0, NP, NQ, IMB1, INB1, MB, NB, MRROW, MRCOL,
1094 $ LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC, LMBLOC,
1095 $ LNBLOC, ILOW, LOW, IUPP, UPP )
1096*
1097 IOFFA = IIA - 1
1098 JOFFA = JJA - 1
1099 LDA = DESCA2( LLD_ )
1100 LDAP1 = LDA + 1
1101*
1102.LT. IF( DESCA2( RSRC_ )0 ) THEN
1103 PMB = MB
1104 ELSE
1105 PMB = NPROW * MB
1106 END IF
1107.LT. IF( DESCA2( CSRC_ )0 ) THEN
1108 QNB = NB
1109 ELSE
1110 QNB = NPCOL * NB
1111 END IF
1112*
1113* Handle the first block of rows or columns separately, and update
1114* LCMT00, MBLKS and NBLKS.
1115*
1116.GT. GODOWN = ( LCMT00IUPP )
1117.LT. GOLEFT = ( LCMT00ILOW )
1118*
1119.NOT..AND..NOT. IF( GODOWN GOLEFT ) THEN
1120*
1121* LCMT00 >= ILOW && LCMT00 <= IUPP
1122*
1123.GE. IF( LCMT000 ) THEN
1124 IJOFFA = IOFFA+LCMT00 + ( JOFFA - 1 ) * LDA
1125 DO 10 I = 1, MIN( INBLOC, MAX( 0, IMBLOC - LCMT00 ) )
1126 ATMP = A( IJOFFA + I*LDAP1 )
1127 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1128 10 CONTINUE
1129 ELSE
1130 IJOFFA = IOFFA + ( JOFFA - LCMT00 - 1 ) * LDA
1131 DO 20 I = 1, MIN( IMBLOC, MAX( 0, INBLOC + LCMT00 ) )
1132 ATMP = A( IJOFFA + I*LDAP1 )
1133 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1134 20 CONTINUE
1135 END IF
1136.LT. GOLEFT = ( ( LCMT00 - ( IUPP - UPP + PMB ) )ILOW )
1137.NOT. GODOWN = GOLEFT
1138*
1139 END IF
1140*
1141 IF( GODOWN ) THEN
1142*
1143 LCMT00 = LCMT00 - ( IUPP - UPP + PMB )
1144 MBLKS = MBLKS - 1
1145 IOFFA = IOFFA + IMBLOC
1146*
1147 30 CONTINUE
1148.GT..AND..GT. IF( MBLKS0 LCMT00UPP ) THEN
1149 LCMT00 = LCMT00 - PMB
1150 MBLKS = MBLKS - 1
1151 IOFFA = IOFFA + MB
1152 GO TO 30
1153 END IF
1154*
1155 LCMT = LCMT00
1156 MBLKD = MBLKS
1157 IOFFD = IOFFA
1158*
1159 MBLOC = MB
1160 40 CONTINUE
1161.GT..AND..GE. IF( MBLKD0 LCMTILOW ) THEN
1162.EQ. IF( MBLKD1 )
1163 $ MBLOC = LMBLOC
1164.GE. IF( LCMT0 ) THEN
1165 IJOFFA = IOFFD + LCMT + ( JOFFA - 1 ) * LDA
1166 DO 50 I = 1, MIN( INBLOC, MAX( 0, MBLOC - LCMT ) )
1167 ATMP = A( IJOFFA + I*LDAP1 )
1168 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1169 50 CONTINUE
1170 ELSE
1171 IJOFFA = IOFFD + ( JOFFA - LCMT - 1 ) * LDA
1172 DO 60 I = 1, MIN( MBLOC, MAX( 0, INBLOC + LCMT ) )
1173 ATMP = A( IJOFFA + I*LDAP1 )
1174 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1175 60 CONTINUE
1176 END IF
1177 LCMT00 = LCMT
1178 LCMT = LCMT - PMB
1179 MBLKS = MBLKD
1180 MBLKD = MBLKD - 1
1181 IOFFA = IOFFD
1182 IOFFD = IOFFD + MBLOC
1183 GO TO 40
1184 END IF
1185*
1186 LCMT00 = LCMT00 + LOW - ILOW + QNB
1187 NBLKS = NBLKS - 1
1188 JOFFA = JOFFA + INBLOC
1189*
1190 ELSE IF( GOLEFT ) THEN
1191*
1192 LCMT00 = LCMT00 + LOW - ILOW + QNB
1193 NBLKS = NBLKS - 1
1194 JOFFA = JOFFA + INBLOC
1195*
1196 70 CONTINUE
1197.GT..AND..LT. IF( NBLKS0 LCMT00LOW ) THEN
1198 LCMT00 = LCMT00 + QNB
1199 NBLKS = NBLKS - 1
1200 JOFFA = JOFFA + NB
1201 GO TO 70
1202 END IF
1203*
1204 LCMT = LCMT00
1205 NBLKD = NBLKS
1206 JOFFD = JOFFA
1207*
1208 NBLOC = NB
1209 80 CONTINUE
1210.GT..AND..LE. IF( NBLKD0 LCMTIUPP ) THEN
1211.EQ. IF( NBLKD1 )
1212 $ NBLOC = LNBLOC
1213.GE. IF( LCMT0 ) THEN
1214 IJOFFA = IOFFA + LCMT + ( JOFFD - 1 ) * LDA
1215 DO 90 I = 1, MIN( NBLOC, MAX( 0, IMBLOC - LCMT ) )
1216 ATMP = A( IJOFFA + I*LDAP1 )
1217 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1218 90 CONTINUE
1219 ELSE
1220 IJOFFA = IOFFA + ( JOFFD - LCMT - 1 ) * LDA
1221 DO 100 I = 1, MIN( IMBLOC, MAX( 0, NBLOC + LCMT ) )
1222 ATMP = A( IJOFFA + I*LDAP1 )
1223 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1224 100 CONTINUE
1225 END IF
1226 LCMT00 = LCMT
1227 LCMT = LCMT + QNB
1228 NBLKS = NBLKD
1229 NBLKD = NBLKD - 1
1230 JOFFA = JOFFD
1231 JOFFD = JOFFD + NBLOC
1232 GO TO 80
1233 END IF
1234*
1235 LCMT00 = LCMT00 - ( IUPP - UPP + PMB )
1236 MBLKS = MBLKS - 1
1237 IOFFA = IOFFA + IMBLOC
1238*
1239 END IF
1240*
1241 NBLOC = NB
1242 110 CONTINUE
1243.GT. IF( NBLKS0 ) THEN
1244.EQ. IF( NBLKS1 )
1245 $ NBLOC = LNBLOC
1246 120 CONTINUE
1247.GT..AND..GT. IF( MBLKS0 LCMT00UPP ) THEN
1248 LCMT00 = LCMT00 - PMB
1249 MBLKS = MBLKS - 1
1250 IOFFA = IOFFA + MB
1251 GO TO 120
1252 END IF
1253*
1254 LCMT = LCMT00
1255 MBLKD = MBLKS
1256 IOFFD = IOFFA
1257*
1258 MBLOC = MB
1259 130 CONTINUE
1260.GT..AND..GE. IF( MBLKD0 LCMTLOW ) THEN
1261.EQ. IF( MBLKD1 )
1262 $ MBLOC = LMBLOC
1263.GE. IF( LCMT0 ) THEN
1264 IJOFFA = IOFFD + LCMT + ( JOFFA - 1 ) * LDA
1265 DO 140 I = 1, MIN( NBLOC, MAX( 0, MBLOC - LCMT ) )
1266 ATMP = A( IJOFFA + I*LDAP1 )
1267 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1268 140 CONTINUE
1269 ELSE
1270 IJOFFA = IOFFD + ( JOFFA - LCMT - 1 ) * LDA
1271 DO 150 I = 1, MIN( MBLOC, MAX( 0, NBLOC + LCMT ) )
1272 ATMP = A( IJOFFA + I*LDAP1 )
1273 A( IJOFFA + I*LDAP1 ) = ABS( ATMP ) + ALPHA
1274 150 CONTINUE
1275 END IF
1276 LCMT00 = LCMT
1277 LCMT = LCMT - PMB
1278 MBLKS = MBLKD
1279 MBLKD = MBLKD - 1
1280 IOFFA = IOFFD
1281 IOFFD = IOFFD + MBLOC
1282 GO TO 130
1283 END IF
1284*
1285 LCMT00 = LCMT00 + QNB
1286 NBLKS = NBLKS - 1
1287 JOFFA = JOFFA + NBLOC
1288 GO TO 110
1289*
1290 END IF
1291*
1292 RETURN
1293*
1294* End of PDLADOM
1295*
1296 END
1297 SUBROUTINE PB_DLASCAL( UPLO, M, N, IOFFD, ALPHA, A, LDA )
1298*
1299* -- PBLAS test routine (version 2.0) --
1300* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
1301* and University of California, Berkeley.
1302* April 1, 1998
1303*
1304* .. Scalar Arguments ..
1305 CHARACTER*1 UPLO
1306 INTEGER IOFFD, LDA, M, N
1307 DOUBLE PRECISION ALPHA
1308* ..
1309* .. Array Arguments ..
1310 DOUBLE PRECISION A( LDA, * )
1311* ..
1312*
1313* Purpose
1314* =======
1315*
1316* PB_DLASCAL scales a two-dimensional array A by the scalar alpha.
1317*
1318* Arguments
1319* =========
1320*
1321* UPLO (input) CHARACTER*1
1322* On entry, UPLO specifies which trapezoidal part of the ar-
1323* ray A is to be scaled as follows:
1324* = 'L' or 'l': the lower trapezoid of A is scaled,
1325* = 'U' or 'u': the upper trapezoid of A is scaled,
1326* = 'D' or 'd': diagonal specified by IOFFD is scaled,
1327* Otherwise: all of the array A is scaled.
1328*
1329* M (input) INTEGER
1330* On entry, M specifies the number of rows of the array A. M
1331* must be at least zero.
1332*
1333* N (input) INTEGER
1334* On entry, N specifies the number of columns of the array A.
1335* N must be at least zero.
1336*
1337* IOFFD (input) INTEGER
1338* On entry, IOFFD specifies the position of the offdiagonal de-
1339* limiting the upper and lower trapezoidal part of A as follows
1340* (see the notes below):
1341*
1342* IOFFD = 0 specifies the main diagonal A( i, i ),
1343* with i = 1 ... MIN( M, N ),
1344* IOFFD > 0 specifies the subdiagonal A( i+IOFFD, i ),
1345* with i = 1 ... MIN( M-IOFFD, N ),
1346* IOFFD < 0 specifies the superdiagonal A( i, i-IOFFD ),
1347* with i = 1 ... MIN( M, N+IOFFD ).
1348*
1349* ALPHA (input) DOUBLE PRECISION
1350* On entry, ALPHA specifies the scalar alpha.
1351*
1352* A (input/output) DOUBLE PRECISION array
1353* On entry, A is an array of dimension (LDA,N). Before entry
1354* with UPLO = 'U' or 'u', the leading m by n part of the array
1355* A must contain the upper trapezoidal part of the matrix as
1356* specified by IOFFD to be scaled, and the strictly lower tra-
1357* pezoidal part of A is not referenced; When UPLO = 'L' or 'l',
1358* the leading m by n part of the array A must contain the lower
1359* trapezoidal part of the matrix as specified by IOFFD to be
1360* scaled, and the strictly upper trapezoidal part of A is not
1361* referenced. On exit, the entries of the trapezoid part of A
1362* determined by UPLO and IOFFD are scaled.
1363*
1364* LDA (input) INTEGER
1365* On entry, LDA specifies the leading dimension of the array A.
1366* LDA must be at least max( 1, M ).
1367*
1368* Notes
1369* =====
1370* N N
1371* ---------------------------- -----------
1372* | d | | |
1373* M | d 'U' | | 'U' |
1374* | 'L' 'D' | |d |
1375* | d | M | d |
1376* ---------------------------- | 'D' |
1377* | d |
1378* IOFFD < 0 | 'L' d |
1379* | d|
1380* N | |
1381* ----------- -----------
1382* | d 'U'|
1383* | d | IOFFD > 0
1384* M | 'D' |
1385* | d| N
1386* | 'L' | ----------------------------
1387* | | | 'U' |
1388* | | |d |
1389* | | | 'D' |
1390* | | | d |
1391* | | |'L' d |
1392* ----------- ----------------------------
1393*
1394* -- Written on April 1, 1998 by
1395* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
1396*
1397* =====================================================================
1398*
1399* .. Local Scalars ..
1400 INTEGER I, J, JTMP, MN
1401* ..
1402* .. External Functions ..
1403 LOGICAL LSAME
1404 EXTERNAL LSAME
1405* ..
1406* .. Intrinsic Functions ..
1407 INTRINSIC MAX, MIN
1408* ..
1409* .. Executable Statements ..
1410*
1411* Quick return if possible
1412*
1413.LE..OR..LE. IF( M0 N0 )
1414 $ RETURN
1415*
1416* Start the operations
1417*
1418 IF( LSAME( UPLO, 'l' ) ) THEN
1419*
1420* Scales the lower triangular part of the array by ALPHA.
1421*
1422 MN = MAX( 0, -IOFFD )
1423 DO 20 J = 1, MIN( MN, N )
1424 DO 10 I = 1, M
1425 A( I, J ) = ALPHA * A( I, J )
1426 10 CONTINUE
1427 20 CONTINUE
1428 DO 40 J = MN + 1, MIN( M - IOFFD, N )
1429 DO 30 I = J + IOFFD, M
1430 A( I, J ) = ALPHA * A( I, J )
1431 30 CONTINUE
1432 40 CONTINUE
1433*
1434 ELSE IF( LSAME( UPLO, 'u' ) ) THEN
1435*
1436* Scales the upper triangular part of the array by ALPHA.
1437*
1438 MN = MIN( M - IOFFD, N )
1439 DO 60 J = MAX( 0, -IOFFD ) + 1, MN
1440 DO 50 I = 1, J + IOFFD
1441 A( I, J ) = ALPHA * A( I, J )
1442 50 CONTINUE
1443 60 CONTINUE
1444 DO 80 J = MAX( 0, MN ) + 1, N
1445 DO 70 I = 1, M
1446 A( I, J ) = ALPHA * A( I, J )
1447 70 CONTINUE
1448 80 CONTINUE
1449*
1450 ELSE IF( LSAME( UPLO, 'd' ) ) THEN
1451*
1452* Scales the diagonal entries by ALPHA.
1453*
1454 DO 90 J = MAX( 0, -IOFFD ) + 1, MIN( M - IOFFD, N )
1455 JTMP = J + IOFFD
1456 A( JTMP, J ) = ALPHA * A( JTMP, J )
1457 90 CONTINUE
1458*
1459 ELSE
1460*
1461* Scales the entire array by ALPHA.
1462*
1463 DO 110 J = 1, N
1464 DO 100 I = 1, M
1465 A( I, J ) = ALPHA * A( I, J )
1466 100 CONTINUE
1467 110 CONTINUE
1468*
1469 END IF
1470*
1471 RETURN
1472*
1473* End of PB_DLASCAL
1474*
1475 END
1476 SUBROUTINE PB_DLAGEN( UPLO, AFORM, A, LDA, LCMT00, IRAN, MBLKS,
1477 $ IMBLOC, MB, LMBLOC, NBLKS, INBLOC, NB,
1478 $ LNBLOC, JMP, IMULADD )
1479*
1480* -- PBLAS test routine (version 2.0) --
1481* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
1482* and University of California, Berkeley.
1483* April 1, 1998
1484*
1485* .. Scalar Arguments ..
1486 CHARACTER*1 UPLO, AFORM
1487 INTEGER IMBLOC, INBLOC, LCMT00, LDA, LMBLOC, LNBLOC,
1488 $ MB, MBLKS, NB, NBLKS
1489* ..
1490* .. Array Arguments ..
1491 INTEGER IMULADD( 4, * ), IRAN( * ), JMP( * )
1492 DOUBLE PRECISION A( LDA, * )
1493* ..
1494*
1495* Purpose
1496* =======
1497*
1498* PB_DLAGEN locally initializes an array A.
1499*
1500* Arguments
1501* =========
1502*
1503* UPLO (global input) CHARACTER*1
1504* On entry, UPLO specifies whether the lower (UPLO='L') trape-
1505* zoidal part or the upper (UPLO='U') trapezoidal part is to be
1506* generated when the matrix to be generated is symmetric or
1507* Hermitian. For all the other values of AFORM, the value of
1508* this input argument is ignored.
1509*
1510* AFORM (global input) CHARACTER*1
1511* On entry, AFORM specifies the type of submatrix to be genera-
1512* ted as follows:
1513* AFORM = 'S', sub( A ) is a symmetric matrix,
1514* AFORM = 'H', sub( A ) is a Hermitian matrix,
1515* AFORM = 'T', sub( A ) is overrwritten with the transpose
1516* of what would normally be generated,
1517* AFORM = 'C', sub( A ) is overwritten with the conjugate
1518* transpose of what would normally be genera-
1519* ted.
1520* AFORM = 'N', a random submatrix is generated.
1521*
1522* A (local output) DOUBLE PRECISION array
1523* On entry, A is an array of dimension (LLD_A, *). On exit,
1524* this array contains the local entries of the randomly genera-
1525* ted submatrix sub( A ).
1526*
1527* LDA (local input) INTEGER
1528* On entry, LDA specifies the local leading dimension of the
1529* array A. LDA must be at least one.
1530*
1531* LCMT00 (global input) INTEGER
1532* On entry, LCMT00 is the LCM value specifying the off-diagonal
1533* of the underlying matrix of interest. LCMT00=0 specifies the
1534* main diagonal, LCMT00 > 0 specifies a subdiagonal, LCMT00 < 0
1535* specifies superdiagonals.
1536*
1537* IRAN (local input) INTEGER array
1538* On entry, IRAN is an array of dimension 2 containing respec-
1539* tively the 16-lower and 16-higher bits of the encoding of the
1540* entry of the random sequence corresponding locally to the
1541* first local array entry to generate. Usually, this array is
1542* computed by PB_SETLOCRAN.
1543*
1544* MBLKS (local input) INTEGER
1545* On entry, MBLKS specifies the local number of blocks of rows.
1546* MBLKS is at least zero.
1547*
1548* IMBLOC (local input) INTEGER
1549* On entry, IMBLOC specifies the number of rows (size) of the
1550* local uppest blocks. IMBLOC is at least zero.
1551*
1552* MB (global input) INTEGER
1553* On entry, MB specifies the blocking factor used to partition
1554* the rows of the matrix. MB must be at least one.
1555*
1556* LMBLOC (local input) INTEGER
1557* On entry, LMBLOC specifies the number of rows (size) of the
1558* local lowest blocks. LMBLOC is at least zero.
1559*
1560* NBLKS (local input) INTEGER
1561* On entry, NBLKS specifies the local number of blocks of co-
1562* lumns. NBLKS is at least zero.
1563*
1564* INBLOC (local input) INTEGER
1565* On entry, INBLOC specifies the number of columns (size) of
1566* the local leftmost blocks. INBLOC is at least zero.
1567*
1568* NB (global input) INTEGER
1569* On entry, NB specifies the blocking factor used to partition
1570* the the columns of the matrix. NB must be at least one.
1571*
1572* LNBLOC (local input) INTEGER
1573* On entry, LNBLOC specifies the number of columns (size) of
1574* the local rightmost blocks. LNBLOC is at least zero.
1575*
1576* JMP (local input) INTEGER array
1577* On entry, JMP is an array of dimension JMP_LEN containing the
1578* different jump values used by the random matrix generator.
1579*
1580* IMULADD (local input) INTEGER array
1581* On entry, IMULADD is an array of dimension (4, JMP_LEN). The
1582* jth column of this array contains the encoded initial cons-
1583* tants a_j and c_j to jump from X( n ) to X( n + JMP( j ) )
1584* (= a_j * X( n ) + c_j) in the random sequence. IMULADD(1:2,j)
1585* contains respectively the 16-lower and 16-higher bits of the
1586* constant a_j, and IMULADD(3:4,j) contains the 16-lower and
1587* 16-higher bits of the constant c_j.
1588*
1589* -- Written on April 1, 1998 by
1590* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
1591*
1592* =====================================================================
1593*
1594* .. Parameters ..
1595 INTEGER JMP_1, JMP_COL, JMP_IMBV, JMP_INBV, JMP_LEN,
1596 $ JMP_MB, JMP_NB, JMP_NPIMBLOC, JMP_NPMB,
1597 $ JMP_NQINBLOC, JMP_NQNB, JMP_ROW
1598 PARAMETER ( JMP_1 = 1, JMP_ROW = 2, JMP_COL = 3,
1599 $ JMP_MB = 4, JMP_IMBV = 5, JMP_NPMB = 6,
1600 $ JMP_NPIMBLOC = 7, JMP_NB = 8, JMP_INBV = 9,
1601 $ JMP_NQNB = 10, JMP_NQINBLOC = 11,
1602 $ JMP_LEN = 11 )
1603* ..
1604* .. Local Scalars ..
1605 INTEGER I, IB, IBLK, II, IK, ITMP, JB, JBLK, JJ, JK,
1606 $ JTMP, LCMTC, LCMTR, LOW, MNB, UPP
1607 DOUBLE PRECISION DUMMY
1608* ..
1609* .. Local Arrays ..
1610 INTEGER IB0( 2 ), IB1( 2 ), IB2( 2 ), IB3( 2 )
1611* ..
1612* .. External Subroutines ..
1613 EXTERNAL PB_JUMPIT
1614* ..
1615* .. External Functions ..
1616 LOGICAL LSAME
1617 DOUBLE PRECISION PB_DRAND
1618 EXTERNAL LSAME, PB_DRAND
1619* ..
1620* .. Intrinsic Functions ..
1621 INTRINSIC MAX, MIN
1622* ..
1623* .. Executable Statements ..
1624*
1625 DO 10 I = 1, 2
1626 IB1( I ) = IRAN( I )
1627 IB2( I ) = IRAN( I )
1628 IB3( I ) = IRAN( I )
1629 10 CONTINUE
1630*
1631 IF( LSAME( AFORM, 'n' ) ) THEN
1632*
1633* Generate random matrix
1634*
1635 JJ = 1
1636*
1637 DO 50 JBLK = 1, NBLKS
1638*
1639.EQ. IF( JBLK1 ) THEN
1640 JB = INBLOC
1641.EQ. ELSE IF( JBLKNBLKS ) THEN
1642 JB = LNBLOC
1643 ELSE
1644 JB = NB
1645 END IF
1646*
1647 DO 40 JK = JJ, JJ + JB - 1
1648*
1649 II = 1
1650*
1651 DO 30 IBLK = 1, MBLKS
1652*
1653.EQ. IF( IBLK1 ) THEN
1654 IB = IMBLOC
1655.EQ. ELSE IF( IBLKMBLKS ) THEN
1656 IB = LMBLOC
1657 ELSE
1658 IB = MB
1659 END IF
1660*
1661* Blocks are IB by JB
1662*
1663 DO 20 IK = II, II + IB - 1
1664 A( IK, JK ) = PB_DRAND( 0 )
1665 20 CONTINUE
1666*
1667 II = II + IB
1668*
1669.EQ. IF( IBLK1 ) THEN
1670*
1671* Jump IMBLOC + ( NPROW - 1 ) * MB rows
1672*
1673 CALL PB_JUMPIT( IMULADD( 1, JMP_NPIMBLOC ), IB1,
1674 $ IB0 )
1675*
1676 ELSE
1677*
1678* Jump NPROW * MB rows
1679*
1680 CALL PB_JUMPIT( IMULADD( 1, JMP_NPMB ), IB1, IB0 )
1681*
1682 END IF
1683*
1684 IB1( 1 ) = IB0( 1 )
1685 IB1( 2 ) = IB0( 2 )
1686*
1687 30 CONTINUE
1688*
1689* Jump one column
1690*
1691 CALL PB_JUMPIT( IMULADD( 1, JMP_COL ), IB2, IB0 )
1692*
1693 IB1( 1 ) = IB0( 1 )
1694 IB1( 2 ) = IB0( 2 )
1695 IB2( 1 ) = IB0( 1 )
1696 IB2( 2 ) = IB0( 2 )
1697*
1698 40 CONTINUE
1699*
1700 JJ = JJ + JB
1701*
1702.EQ. IF( JBLK1 ) THEN
1703*
1704* Jump INBLOC + ( NPCOL - 1 ) * NB columns
1705*
1706 CALL PB_JUMPIT( IMULADD( 1, JMP_NQINBLOC ), IB3, IB0 )
1707*
1708 ELSE
1709*
1710* Jump NPCOL * NB columns
1711*
1712 CALL PB_JUMPIT( IMULADD( 1, JMP_NQNB ), IB3, IB0 )
1713*
1714 END IF
1715*
1716 IB1( 1 ) = IB0( 1 )
1717 IB1( 2 ) = IB0( 2 )
1718 IB2( 1 ) = IB0( 1 )
1719 IB2( 2 ) = IB0( 2 )
1720 IB3( 1 ) = IB0( 1 )
1721 IB3( 2 ) = IB0( 2 )
1722*
1723 50 CONTINUE
1724*
1725 ELSE IF( LSAME( AFORM, 't.OR.' ) LSAME( AFORM, 'c' ) ) THEN
1726*
1727* Generate the transpose of the matrix that would be normally
1728* generated.
1729*
1730 II = 1
1731*
1732 DO 90 IBLK = 1, MBLKS
1733*
1734.EQ. IF( IBLK1 ) THEN
1735 IB = IMBLOC
1736.EQ. ELSE IF( IBLKMBLKS ) THEN
1737 IB = LMBLOC
1738 ELSE
1739 IB = MB
1740 END IF
1741*
1742 DO 80 IK = II, II + IB - 1
1743*
1744 JJ = 1
1745*
1746 DO 70 JBLK = 1, NBLKS
1747*
1748.EQ. IF( JBLK1 ) THEN
1749 JB = INBLOC
1750.EQ. ELSE IF( JBLKNBLKS ) THEN
1751 JB = LNBLOC
1752 ELSE
1753 JB = NB
1754 END IF
1755*
1756* Blocks are IB by JB
1757*
1758 DO 60 JK = JJ, JJ + JB - 1
1759 A( IK, JK ) = PB_DRAND( 0 )
1760 60 CONTINUE
1761*
1762 JJ = JJ + JB
1763*
1764.EQ. IF( JBLK1 ) THEN
1765*
1766* Jump INBLOC + ( NPCOL - 1 ) * NB columns
1767*
1768 CALL PB_JUMPIT( IMULADD( 1, JMP_NQINBLOC ), IB1,
1769 $ IB0 )
1770*
1771 ELSE
1772*
1773* Jump NPCOL * NB columns
1774*
1775 CALL PB_JUMPIT( IMULADD( 1, JMP_NQNB ), IB1, IB0 )
1776*
1777 END IF
1778*
1779 IB1( 1 ) = IB0( 1 )
1780 IB1( 2 ) = IB0( 2 )
1781*
1782 70 CONTINUE
1783*
1784* Jump one row
1785*
1786 CALL PB_JUMPIT( IMULADD( 1, JMP_ROW ), IB2, IB0 )
1787*
1788 IB1( 1 ) = IB0( 1 )
1789 IB1( 2 ) = IB0( 2 )
1790 IB2( 1 ) = IB0( 1 )
1791 IB2( 2 ) = IB0( 2 )
1792*
1793 80 CONTINUE
1794*
1795 II = II + IB
1796*
1797.EQ. IF( IBLK1 ) THEN
1798*
1799* Jump IMBLOC + ( NPROW - 1 ) * MB rows
1800*
1801 CALL PB_JUMPIT( IMULADD( 1, JMP_NPIMBLOC ), IB3, IB0 )
1802*
1803 ELSE
1804*
1805* Jump NPROW * MB rows
1806*
1807 CALL PB_JUMPIT( IMULADD( 1, JMP_NPMB ), IB3, IB0 )
1808*
1809 END IF
1810*
1811 IB1( 1 ) = IB0( 1 )
1812 IB1( 2 ) = IB0( 2 )
1813 IB2( 1 ) = IB0( 1 )
1814 IB2( 2 ) = IB0( 2 )
1815 IB3( 1 ) = IB0( 1 )
1816 IB3( 2 ) = IB0( 2 )
1817*
1818 90 CONTINUE
1819*
1820 ELSE IF( ( LSAME( AFORM, 's.OR.' ) )( LSAME( AFORM, 'h' ) ) ) THEN
1821*
1822* Generate a symmetric matrix
1823*
1824 IF( LSAME( UPLO, 'l' ) ) THEN
1825*
1826* generate lower trapezoidal part
1827*
1828 JJ = 1
1829 LCMTC = LCMT00
1830*
1831 DO 170 JBLK = 1, NBLKS
1832*
1833.EQ. IF( JBLK1 ) THEN
1834 JB = INBLOC
1835 LOW = 1 - INBLOC
1836.EQ. ELSE IF( JBLKNBLKS ) THEN
1837 JB = LNBLOC
1838 LOW = 1 - NB
1839 ELSE
1840 JB = NB
1841 LOW = 1 - NB
1842 END IF
1843*
1844 DO 160 JK = JJ, JJ + JB - 1
1845*
1846 II = 1
1847 LCMTR = LCMTC
1848*
1849 DO 150 IBLK = 1, MBLKS
1850*
1851.EQ. IF( IBLK1 ) THEN
1852 IB = IMBLOC
1853 UPP = IMBLOC - 1
1854.EQ. ELSE IF( IBLKMBLKS ) THEN
1855 IB = LMBLOC
1856 UPP = MB - 1
1857 ELSE
1858 IB = MB
1859 UPP = MB - 1
1860 END IF
1861*
1862* Blocks are IB by JB
1863*
1864.GT. IF( LCMTRUPP ) THEN
1865*
1866 DO 100 IK = II, II + IB - 1
1867 DUMMY = PB_DRAND( 0 )
1868 100 CONTINUE
1869*
1870.GE. ELSE IF( LCMTRLOW ) THEN
1871*
1872 JTMP = JK - JJ + 1
1873 MNB = MAX( 0, -LCMTR )
1874*
1875.LE. IF( JTMPMIN( MNB, JB ) ) THEN
1876*
1877 DO 110 IK = II, II + IB - 1
1878 A( IK, JK ) = PB_DRAND( 0 )
1879 110 CONTINUE
1880*
1881.GE..AND. ELSE IF( ( JTMP( MNB + 1 ) )
1882.LE. $ ( JTMPMIN( IB-LCMTR, JB ) ) ) THEN
1883*
1884 ITMP = II + JTMP + LCMTR - 1
1885*
1886 DO 120 IK = II, ITMP - 1
1887 DUMMY = PB_DRAND( 0 )
1888 120 CONTINUE
1889*
1890 DO 130 IK = ITMP, II + IB - 1
1891 A( IK, JK ) = PB_DRAND( 0 )
1892 130 CONTINUE
1893*
1894 END IF
1895*
1896 ELSE
1897*
1898 DO 140 IK = II, II + IB - 1
1899 A( IK, JK ) = PB_DRAND( 0 )
1900 140 CONTINUE
1901*
1902 END IF
1903*
1904 II = II + IB
1905*
1906.EQ. IF( IBLK1 ) THEN
1907*
1908* Jump IMBLOC + ( NPROW - 1 ) * MB rows
1909*
1910 LCMTR = LCMTR - JMP( JMP_NPIMBLOC )
1911 CALL PB_JUMPIT( IMULADD( 1, JMP_NPIMBLOC ), IB1,
1912 $ IB0 )
1913*
1914 ELSE
1915*
1916* Jump NPROW * MB rows
1917*
1918 LCMTR = LCMTR - JMP( JMP_NPMB )
1919 CALL PB_JUMPIT( IMULADD( 1, JMP_NPMB ), IB1,
1920 $ IB0 )
1921*
1922 END IF
1923*
1924 IB1( 1 ) = IB0( 1 )
1925 IB1( 2 ) = IB0( 2 )
1926*
1927 150 CONTINUE
1928*
1929* Jump one column
1930*
1931 CALL PB_JUMPIT( IMULADD( 1, JMP_COL ), IB2, IB0 )
1932*
1933 IB1( 1 ) = IB0( 1 )
1934 IB1( 2 ) = IB0( 2 )
1935 IB2( 1 ) = IB0( 1 )
1936 IB2( 2 ) = IB0( 2 )
1937*
1938 160 CONTINUE
1939*
1940 JJ = JJ + JB
1941*
1942.EQ. IF( JBLK1 ) THEN
1943*
1944* Jump INBLOC + ( NPCOL - 1 ) * NB columns
1945*
1946 LCMTC = LCMTC + JMP( JMP_NQINBLOC )
1947 CALL PB_JUMPIT( IMULADD( 1, JMP_NQINBLOC ), IB3, IB0 )
1948*
1949 ELSE
1950*
1951* Jump NPCOL * NB columns
1952*
1953 LCMTC = LCMTC + JMP( JMP_NQNB )
1954 CALL PB_JUMPIT( IMULADD( 1, JMP_NQNB ), IB3, IB0 )
1955*
1956 END IF
1957*
1958 IB1( 1 ) = IB0( 1 )
1959 IB1( 2 ) = IB0( 2 )
1960 IB2( 1 ) = IB0( 1 )
1961 IB2( 2 ) = IB0( 2 )
1962 IB3( 1 ) = IB0( 1 )
1963 IB3( 2 ) = IB0( 2 )
1964*
1965 170 CONTINUE
1966*
1967 ELSE
1968*
1969* generate upper trapezoidal part
1970*
1971 II = 1
1972 LCMTR = LCMT00
1973*
1974 DO 250 IBLK = 1, MBLKS
1975*
1976.EQ. IF( IBLK1 ) THEN
1977 IB = IMBLOC
1978 UPP = IMBLOC - 1
1979.EQ. ELSE IF( IBLKMBLKS ) THEN
1980 IB = LMBLOC
1981 UPP = MB - 1
1982 ELSE
1983 IB = MB
1984 UPP = MB - 1
1985 END IF
1986*
1987 DO 240 IK = II, II + IB - 1
1988*
1989 JJ = 1
1990 LCMTC = LCMTR
1991*
1992 DO 230 JBLK = 1, NBLKS
1993*
1994.EQ. IF( JBLK1 ) THEN
1995 JB = INBLOC
1996 LOW = 1 - INBLOC
1997.EQ. ELSE IF( JBLKNBLKS ) THEN
1998 JB = LNBLOC
1999 LOW = 1 - NB
2000 ELSE
2001 JB = NB
2002 LOW = 1 - NB
2003 END IF
2004*
2005* Blocks are IB by JB
2006*
2007.LT. IF( LCMTCLOW ) THEN
2008*
2009 DO 180 JK = JJ, JJ + JB - 1
2010 DUMMY = PB_DRAND( 0 )
2011 180 CONTINUE
2012*
2013.LE. ELSE IF( LCMTCUPP ) THEN
2014*
2015 ITMP = IK - II + 1
2016 MNB = MAX( 0, LCMTC )
2017*
2018.LE. IF( ITMPMIN( MNB, IB ) ) THEN
2019*
2020 DO 190 JK = JJ, JJ + JB - 1
2021 A( IK, JK ) = PB_DRAND( 0 )
2022 190 CONTINUE
2023*
2024.GE..AND. ELSE IF( ( ITMP( MNB + 1 ) )
2025.LE. $ ( ITMPMIN( JB+LCMTC, IB ) ) ) THEN
2026*
2027 JTMP = JJ + ITMP - LCMTC - 1
2028*
2029 DO 200 JK = JJ, JTMP - 1
2030 DUMMY = PB_DRAND( 0 )
2031 200 CONTINUE
2032*
2033 DO 210 JK = JTMP, JJ + JB - 1
2034 A( IK, JK ) = PB_DRAND( 0 )
2035 210 CONTINUE
2036*
2037 END IF
2038*
2039 ELSE
2040*
2041 DO 220 JK = JJ, JJ + JB - 1
2042 A( IK, JK ) = PB_DRAND( 0 )
2043 220 CONTINUE
2044*
2045 END IF
2046*
2047 JJ = JJ + JB
2048*
2049.EQ. IF( JBLK1 ) THEN
2050*
2051* Jump INBLOC + ( NPCOL - 1 ) * NB columns
2052*
2053 LCMTC = LCMTC + JMP( JMP_NQINBLOC )
2054 CALL PB_JUMPIT( IMULADD( 1, JMP_NQINBLOC ), IB1,
2055 $ IB0 )
2056*
2057 ELSE
2058*
2059* Jump NPCOL * NB columns
2060*
2061 LCMTC = LCMTC + JMP( JMP_NQNB )
2062 CALL PB_JUMPIT( IMULADD( 1, JMP_NQNB ), IB1,
2063 $ IB0 )
2064*
2065 END IF
2066*
2067 IB1( 1 ) = IB0( 1 )
2068 IB1( 2 ) = IB0( 2 )
2069*
2070 230 CONTINUE
2071*
2072* Jump one row
2073*
2074 CALL PB_JUMPIT( IMULADD( 1, JMP_ROW ), IB2, IB0 )
2075*
2076 IB1( 1 ) = IB0( 1 )
2077 IB1( 2 ) = IB0( 2 )
2078 IB2( 1 ) = IB0( 1 )
2079 IB2( 2 ) = IB0( 2 )
2080*
2081 240 CONTINUE
2082*
2083 II = II + IB
2084*
2085.EQ. IF( IBLK1 ) THEN
2086*
2087* Jump IMBLOC + ( NPROW - 1 ) * MB rows
2088*
2089 LCMTR = LCMTR - JMP( JMP_NPIMBLOC )
2090 CALL PB_JUMPIT( IMULADD( 1, JMP_NPIMBLOC ), IB3, IB0 )
2091*
2092 ELSE
2093*
2094* Jump NPROW * MB rows
2095*
2096 LCMTR = LCMTR - JMP( JMP_NPMB )
2097 CALL PB_JUMPIT( IMULADD( 1, JMP_NPMB ), IB3, IB0 )
2098*
2099 END IF
2100*
2101 IB1( 1 ) = IB0( 1 )
2102 IB1( 2 ) = IB0( 2 )
2103 IB2( 1 ) = IB0( 1 )
2104 IB2( 2 ) = IB0( 2 )
2105 IB3( 1 ) = IB0( 1 )
2106 IB3( 2 ) = IB0( 2 )
2107*
2108 250 CONTINUE
2109*
2110 END IF
2111*
2112 END IF
2113*
2114 RETURN
2115*
2116* End of PB_DLAGEN
2117*
2118 END
2119 DOUBLE PRECISION FUNCTION PB_DRAND( IDUMM )
2120*
2121* -- PBLAS test routine (version 2.0) --
2122* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
2123* and University of California, Berkeley.
2124* April 1, 1998
2125*
2126* .. Scalar Arguments ..
2127 INTEGER IDUMM
2128* ..
2129*
2130* Purpose
2131* =======
2132*
2133* PB_DRAND generates the next number in the random sequence. This func-
2134* tion ensures that this number will be in the interval ( -1.0, 1.0 ).
2135*
2136* Arguments
2137* =========
2138*
2139* IDUMM (local input) INTEGER
2140* This argument is ignored, but necessary to a FORTRAN 77 func-
2141* tion.
2142*
2143* Further Details
2144* ===============
2145*
2146* On entry, the array IRAND stored in the common block RANCOM contains
2147* the information (2 integers) required to generate the next number in
2148* the sequence X( n ). This number is computed as
2149*
2150* X( n ) = ( 2^16 * IRAND( 2 ) + IRAND( 1 ) ) / d,
2151*
2152* where the constant d is the largest 32 bit positive integer. The
2153* array IRAND is then updated for the generation of the next number
2154* X( n+1 ) in the random sequence as follows X( n+1 ) = a * X( n ) + c.
2155* The constants a and c should have been preliminarily stored in the
2156* array IACS as 2 pairs of integers. The initial set up of IRAND and
2157* IACS is performed by the routine PB_SETRAN.
2158*
2159* -- Written on April 1, 1998 by
2160* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
2161*
2162* =====================================================================
2163*
2164* .. Parameters ..
2165 DOUBLE PRECISION ONE, TWO
2166 PARAMETER ( ONE = 1.0D+0, TWO = 2.0D+0 )
2167* ..
2168* .. External Functions ..
2169 DOUBLE PRECISION PB_DRAN
2170 EXTERNAL PB_DRAN
2171* ..
2172* .. Executable Statements ..
2173*
2174 PB_DRAND = ONE - TWO * PB_DRAN( IDUMM )
2175*
2176 RETURN
2177*
2178* End of PB_DRAND
2179*
2180 END
2181 DOUBLE PRECISION FUNCTION PB_DRAN( IDUMM )
2182*
2183* -- PBLAS test routine (version 2.0) --
2184* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
2185* and University of California, Berkeley.
2186* April 1, 1998
2187*
2188* .. Scalar Arguments ..
2189 INTEGER IDUMM
2190* ..
2191*
2192* Purpose
2193* =======
2194*
2195* PB_DRAN generates the next number in the random sequence.
2196*
2197* Arguments
2198* =========
2199*
2200* IDUMM (local input) INTEGER
2201* This argument is ignored, but necessary to a FORTRAN 77 func-
2202* tion.
2203*
2204* Further Details
2205* ===============
2206*
2207* On entry, the array IRAND stored in the common block RANCOM contains
2208* the information (2 integers) required to generate the next number in
2209* the sequence X( n ). This number is computed as
2210*
2211* X( n ) = ( 2^16 * IRAND( 2 ) + IRAND( 1 ) ) / d,
2212*
2213* where the constant d is the largest 32 bit positive integer. The
2214* array IRAND is then updated for the generation of the next number
2215* X( n+1 ) in the random sequence as follows X( n+1 ) = a * X( n ) + c.
2216* The constants a and c should have been preliminarily stored in the
2217* array IACS as 2 pairs of integers. The initial set up of IRAND and
2218* IACS is performed by the routine PB_SETRAN.
2219*
2220* -- Written on April 1, 1998 by
2221* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
2222*
2223* =====================================================================
2224*
2225* .. Parameters ..
2226 DOUBLE PRECISION DIVFAC, POW16
2227 PARAMETER ( DIVFAC = 2.147483648D+9,
2228 $ POW16 = 6.5536D+4 )
2229* ..
2230* .. Local Arrays ..
2231 INTEGER J( 2 )
2232* ..
2233* .. External Subroutines ..
2234 EXTERNAL PB_LADD, PB_LMUL
2235* ..
2236* .. Intrinsic Functions ..
2237 INTRINSIC DBLE
2238* ..
2239* .. Common Blocks ..
2240 INTEGER IACS( 4 ), IRAND( 2 )
2241 COMMON /RANCOM/ IRAND, IACS
2242* ..
2243* .. Save Statements ..
2244 SAVE /RANCOM/
2245* ..
2246* .. Executable Statements ..
2247*
2248 PB_DRAN = ( DBLE( IRAND( 1 ) ) + POW16 * DBLE( IRAND( 2 ) ) ) /
2249 $ DIVFAC
2250*
2251 CALL PB_LMUL( IRAND, IACS, J )
2252 CALL PB_LADD( J, IACS( 3 ), IRAND )
2253*
2254 RETURN
2255*
2256* End of PB_DRAN
2257*
2258 END
#define min(a, b)
Definition macros.h:20
subroutine blacs_gridinfo(cntxt, nprow, npcol, myrow, mycol)
Definition mpi.f:754
subroutine pb_ainfog2l(m, n, i, j, desc, nprow, npcol, myrow, mycol, imb1, inb1, mp, nq, ii, jj, prow, pcol, rprow, rpcol)
Definition pblastst.f:2023
subroutine pb_binfo(offd, m, n, imb1, inb1, mb, nb, mrrow, mrcol, lcmt00, mblks, nblks, imbloc, inbloc, lmbloc, lnbloc, ilow, low, iupp, upp)
Definition pblastst.f:3577
subroutine pb_infog2l(i, j, desc, nprow, npcol, myrow, mycol, ii, jj, prow, pcol)
Definition pblastst.f:1673
subroutine pb_desctrans(descin, descout)
Definition pblastst.f:2964
subroutine pdlagen(inplace, aform, diag, offa, m, n, ia, ja, desca, iaseed, a, lda)
Definition pdblastim.f:510
subroutine pdlascal(type, m, n, alpha, a, ia, ja, desca)
Definition pdblastim.f:2
subroutine pb_dlascal(uplo, m, n, ioffd, alpha, a, lda)
Definition pdblastim.f:1298