ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pzblastim.f
Go to the documentation of this file.
1  SUBROUTINE pzlascal( 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  COMPLEX*16 ALPHA
12 * ..
13 * .. Array Arguments ..
14  INTEGER DESCA( * )
15  COMPLEX*16 A( * )
16 * ..
17 *
18 * Purpose
19 * =======
20 *
21 * PZLASCAL 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) COMPLEX*16
108 * On entry, ALPHA specifies the scalar alpha.
109 *
110 * A (local input/local output) COMPLEX*16 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 ..
159  EXTERNAL blacs_gridinfo, pb_ainfog2l, pb_binfo,
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_zlascal( '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_zlascal( 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_zlascal( '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_zlascal( '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_zlascal( '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_zlascal( 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_zlascal( '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_zlascal( '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_zlascal( '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_zlascal( 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_zlascal( '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_zlascal( '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_zlascal( 'All', tmp1, n1, 0, alpha,
450  $ a( iia+joffa*lda ), lda )
451  iia = iia + tmp1
452  m1 = m1 - tmp1
453  END IF
454 *
455  IF( mblks.LE.0 )
456  $ RETURN
457 *
458  lcmt = lcmt00
459  mblkd = mblks
460  ioffd = ioffa
461 *
462  mbloc = mb
463  70 CONTINUE
464  IF( mblkd.GT.0 .AND. lcmt.GE.low ) THEN
465  IF( mblkd.EQ.1 )
466  $ mbloc = lmbloc
467  CALL pb_zlascal( 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  IF( lower .AND. tmp1.GT.0 )
480  $ CALL pb_zlascal( '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  IF( upper .AND. tmp1.GT.0 .AND. n1.GT.0 )
491  $ CALL pb_zlascal( '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 PZLASCAL
506 *
507  END
508  SUBROUTINE pzlagen( 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  COMPLEX*16 A( LDA, * )
524 * ..
525 *
526 * Purpose
527 * =======
528 *
529 * PZLAGEN 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) COMPLEX*16 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  DOUBLE PRECISION ZERO
717  parameter( zero = 0.0d+0 )
718 * ..
719 * .. Local Scalars ..
720  LOGICAL DIAGDO, SYMM, HERM, NOTRAN
721  INTEGER CSRC, I, IACOL, IAROW, ICTXT, IIA, ILOCBLK,
722  $ ilocoff, ilow, imb, imb1, imbloc, imbvir, inb,
723  $ inb1, inbloc, inbvir, info, ioffda, itmp, iupp,
724  $ ivir, jja, jlocblk, jlocoff, jvir, lcmt00,
725  $ lmbloc, lnbloc, low, maxmn, mb, mblks, mp,
726  $ mrcol, mrrow, mycdist, mycol, myrdist, myrow,
727  $ nb, nblks, npcol, nprow, nq, nvir, rsrc, upp
728  COMPLEX*16 ALPHA
729 * ..
730 * .. Local Arrays ..
731  INTEGER DESCA2( DLEN_ ), IMULADD( 4, JMP_LEN ),
732  $ iran( 2 ), jmp( jmp_len ), muladd0( 4 )
733 * ..
734 * .. External Subroutines ..
735  EXTERNAL blacs_gridinfo, pb_ainfog2l, pb_binfo,
739  $ pzladom
740 * ..
741 * .. External Functions ..
742  LOGICAL LSAME
743  EXTERNAL lsame
744 * ..
745 * .. Intrinsic Functions ..
746  INTRINSIC dble, dcmplx, max, min
747 * ..
748 * .. Data Statements ..
749  DATA ( muladd0( i ), i = 1, 4 ) / 20077, 16838,
750  $ 12345, 0 /
751 * ..
752 * .. Executable Statements ..
753 *
754 * Convert descriptor
755 *
756  CALL pb_desctrans( desca, desca2 )
757 *
758 * Test the input arguments
759 *
760  ictxt = desca2( ctxt_ )
761  CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
762 *
763 * Test the input parameters
764 *
765  info = 0
766  IF( nprow.EQ.-1 ) THEN
767  info = -( 1000 + ctxt_ )
768  ELSE
769  symm = lsame( aform, 'S' )
770  herm = lsame( aform, 'H' )
771  notran = lsame( aform, 'N' )
772  diagdo = lsame( diag, 'D' )
773  IF( .NOT.( symm.OR.herm.OR.notran ) .AND.
774  $ .NOT.( lsame( aform, 'T' ) ) .AND.
775  $ .NOT.( lsame( aform, 'C' ) ) ) THEN
776  info = -2
777  ELSE IF( ( .NOT.diagdo ) .AND.
778  $ ( .NOT.lsame( diag, 'N' ) ) ) THEN
779  info = -3
780  END IF
781  CALL pb_chkmat( ictxt, m, 5, n, 6, ia, ja, desca2, 10, info )
782  END IF
783 *
784  IF( info.NE.0 ) THEN
785  CALL pxerbla( ictxt, 'PZLAGEN', -info )
786  RETURN
787  END IF
788 *
789 * Quick return if possible
790 *
791  IF( ( m.LE.0 ).OR.( n.LE.0 ) )
792  $ RETURN
793 *
794 * Start the operations
795 *
796  mb = desca2( mb_ )
797  nb = desca2( nb_ )
798  imb = desca2( imb_ )
799  inb = desca2( inb_ )
800  rsrc = desca2( rsrc_ )
801  csrc = desca2( csrc_ )
802 *
803 * Figure out local information about the distributed matrix operand
804 *
805  CALL pb_ainfog2l( m, n, ia, ja, desca2, nprow, npcol, myrow,
806  $ mycol, imb1, inb1, mp, nq, iia, jja, iarow,
807  $ iacol, mrrow, mrcol )
808 *
809 * Decide where the entries shall be stored in memory
810 *
811  IF( inplace ) THEN
812  iia = 1
813  jja = 1
814  END IF
815 *
816 * Initialize LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC, LMBLOC, LNBLOC,
817 * ILOW, LOW, IUPP, and UPP.
818 *
819  ioffda = ja + offa - ia
820  CALL pb_binfo( ioffda, mp, nq, imb1, inb1, mb, nb, mrrow,
821  $ mrcol, lcmt00, mblks, nblks, imbloc, inbloc,
822  $ lmbloc, lnbloc, ilow, low, iupp, upp )
823 *
824 * Initialize ILOCBLK, ILOCOFF, MYRDIST, JLOCBLK, JLOCOFF, MYCDIST
825 * This values correspond to the square virtual underlying matrix
826 * of size MAX( M_ + MAX( 0, -OFFA ), N_ + MAX( 0, OFFA ) ) used
827 * to set up the random sequence. For practical purposes, the size
828 * of this virtual matrix is upper bounded by M_ + N_ - 1.
829 *
830  itmp = max( 0, -offa )
831  ivir = ia + itmp
832  imbvir = imb + itmp
833  nvir = desca2( m_ ) + itmp
834 *
835  CALL pb_locinfo( ivir, imbvir, mb, myrow, rsrc, nprow, ilocblk,
836  $ ilocoff, myrdist )
837 *
838  itmp = max( 0, offa )
839  jvir = ja + itmp
840  inbvir = inb + itmp
841  nvir = max( max( nvir, desca2( n_ ) + itmp ),
842  $ desca2( m_ ) + desca2( n_ ) - 1 )
843 *
844  CALL pb_locinfo( jvir, inbvir, nb, mycol, csrc, npcol, jlocblk,
845  $ jlocoff, mycdist )
846 *
847  IF( symm .OR. herm .OR. notran ) THEN
848 *
849  CALL pb_initjmp( .true., nvir, imbvir, inbvir, imbloc, inbloc,
850  $ mb, nb, rsrc, csrc, nprow, npcol, 2, jmp )
851 *
852 * Compute constants to jump JMP( * ) numbers in the sequence
853 *
854  CALL pb_initmuladd( muladd0, jmp, imuladd )
855 *
856 * Compute and set the random value corresponding to A( IA, JA )
857 *
858  CALL pb_setlocran( iaseed, ilocblk, jlocblk, ilocoff, jlocoff,
859  $ myrdist, mycdist, nprow, npcol, jmp,
860  $ imuladd, iran )
861 *
862  CALL pb_zlagen( 'Lower', aform, a( iia, jja ), lda, lcmt00,
863  $ iran, mblks, imbloc, mb, lmbloc, nblks, inbloc,
864  $ nb, lnbloc, jmp, imuladd )
865 *
866  END IF
867 *
868  IF( symm .OR. herm .OR. ( .NOT. notran ) ) THEN
869 *
870  CALL pb_initjmp( .false., nvir, imbvir, inbvir, imbloc, inbloc,
871  $ mb, nb, rsrc, csrc, nprow, npcol, 2, jmp )
872 *
873 * Compute constants to jump JMP( * ) numbers in the sequence
874 *
875  CALL pb_initmuladd( muladd0, jmp, imuladd )
876 *
877 * Compute and set the random value corresponding to A( IA, JA )
878 *
879  CALL pb_setlocran( iaseed, ilocblk, jlocblk, ilocoff, jlocoff,
880  $ myrdist, mycdist, nprow, npcol, jmp,
881  $ imuladd, iran )
882 *
883  CALL pb_zlagen( 'Upper', aform, a( iia, jja ), lda, lcmt00,
884  $ iran, mblks, imbloc, mb, lmbloc, nblks, inbloc,
885  $ nb, lnbloc, jmp, imuladd )
886 *
887  END IF
888 *
889  IF( diagdo ) THEN
890 *
891  maxmn = max( desca2( m_ ), desca2( n_ ) )
892  IF( herm ) THEN
893  alpha = dcmplx( dble( 2 * maxmn ), zero )
894  ELSE
895  alpha = dcmplx( dble( nvir ), dble( maxmn ) )
896  END IF
897 *
898  IF( ioffda.GE.0 ) THEN
899  CALL pzladom( inplace, min( max( 0, m-ioffda ), n ), alpha,
900  $ a, min( ia+ioffda, ia+m-1 ), ja, desca )
901  ELSE
902  CALL pzladom( inplace, min( m, max( 0, n+ioffda ) ), alpha,
903  $ a, ia, min( ja-ioffda, ja+n-1 ), desca )
904  END IF
905 *
906  END IF
907 *
908  RETURN
909 *
910 * End of PZLAGEN
911 *
912  END
913  SUBROUTINE pzladom( INPLACE, N, ALPHA, A, IA, JA, DESCA )
914 *
915 * -- PBLAS test routine (version 2.0) --
916 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
917 * and University of California, Berkeley.
918 * April 1, 1998
919 *
920 * .. Scalar Arguments ..
921  LOGICAL INPLACE
922  INTEGER IA, JA, N
923  COMPLEX*16 ALPHA
924 * ..
925 * .. Array Arguments ..
926  INTEGER DESCA( * )
927  COMPLEX*16 A( * )
928 * ..
929 *
930 * Purpose
931 * =======
932 *
933 * PZLADOM adds alpha to the diagonal entries of an n by n submatrix
934 * sub( A ) denoting A( IA:IA+N-1, JA:JA+N-1 ).
935 *
936 * Notes
937 * =====
938 *
939 * A description vector is associated with each 2D block-cyclicly dis-
940 * tributed matrix. This vector stores the information required to
941 * establish the mapping between a matrix entry and its corresponding
942 * process and memory location.
943 *
944 * In the following comments, the character _ should be read as
945 * "of the distributed matrix". Let A be a generic term for any 2D
946 * block cyclicly distributed matrix. Its description vector is DESCA:
947 *
948 * NOTATION STORED IN EXPLANATION
949 * ---------------- --------------- ------------------------------------
950 * DTYPE_A (global) DESCA( DTYPE_ ) The descriptor type.
951 * CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
952 * the NPROW x NPCOL BLACS process grid
953 * A is distributed over. The context
954 * itself is global, but the handle
955 * (the integer value) may vary.
956 * M_A (global) DESCA( M_ ) The number of rows in the distribu-
957 * ted matrix A, M_A >= 0.
958 * N_A (global) DESCA( N_ ) The number of columns in the distri-
959 * buted matrix A, N_A >= 0.
960 * IMB_A (global) DESCA( IMB_ ) The number of rows of the upper left
961 * block of the matrix A, IMB_A > 0.
962 * INB_A (global) DESCA( INB_ ) The number of columns of the upper
963 * left block of the matrix A,
964 * INB_A > 0.
965 * MB_A (global) DESCA( MB_ ) The blocking factor used to distri-
966 * bute the last M_A-IMB_A rows of A,
967 * MB_A > 0.
968 * NB_A (global) DESCA( NB_ ) The blocking factor used to distri-
969 * bute the last N_A-INB_A columns of
970 * A, NB_A > 0.
971 * RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
972 * row of the matrix A is distributed,
973 * NPROW > RSRC_A >= 0.
974 * CSRC_A (global) DESCA( CSRC_ ) The process column over which the
975 * first column of A is distributed.
976 * NPCOL > CSRC_A >= 0.
977 * LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
978 * array storing the local blocks of
979 * the distributed matrix A,
980 * IF( Lc( 1, N_A ) > 0 )
981 * LLD_A >= MAX( 1, Lr( 1, M_A ) )
982 * ELSE
983 * LLD_A >= 1.
984 *
985 * Let K be the number of rows of a matrix A starting at the global in-
986 * dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
987 * that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
988 * receive if these K rows were distributed over NPROW processes. If K
989 * is the number of columns of a matrix A starting at the global index
990 * JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
991 * lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
992 * these K columns were distributed over NPCOL processes.
993 *
994 * The values of Lr() and Lc() may be determined via a call to the func-
995 * tion PB_NUMROC:
996 * Lr( IA, K ) = PB_NUMROC( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
997 * Lc( JA, K ) = PB_NUMROC( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
998 *
999 * Arguments
1000 * =========
1001 *
1002 * INPLACE (global input) LOGICAL
1003 * On entry, INPLACE specifies if the matrix should be generated
1004 * in place or not. If INPLACE is .TRUE., the local random array
1005 * to be generated will start in memory at the local memory lo-
1006 * cation A( 1, 1 ), otherwise it will start at the local posi-
1007 * tion induced by IA and JA.
1008 *
1009 * N (global input) INTEGER
1010 * On entry, N specifies the global order of the submatrix
1011 * sub( A ) to be modified. N must be at least zero.
1012 *
1013 * ALPHA (global input) COMPLEX*16
1014 * On entry, ALPHA specifies the scalar alpha.
1015 *
1016 * A (local input/local output) COMPLEX*16 array
1017 * On entry, A is an array of dimension (LLD_A, Ka), where Ka is
1018 * at least Lc( 1, JA+N-1 ). Before entry, this array contains
1019 * the local entries of the matrix A. On exit, the local entries
1020 * of this array corresponding to the main diagonal of sub( A )
1021 * have been updated.
1022 *
1023 * IA (global input) INTEGER
1024 * On entry, IA specifies A's global row index, which points to
1025 * the beginning of the submatrix sub( A ).
1026 *
1027 * JA (global input) INTEGER
1028 * On entry, JA specifies A's global column index, which points
1029 * to the beginning of the submatrix sub( A ).
1030 *
1031 * DESCA (global and local input) INTEGER array
1032 * On entry, DESCA is an integer array of dimension DLEN_. This
1033 * is the array descriptor for the matrix A.
1034 *
1035 * -- Written on April 1, 1998 by
1036 * Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
1037 *
1038 * =====================================================================
1039 *
1040 * .. Parameters ..
1041  INTEGER BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
1042  $ dtype_, imb_, inb_, lld_, mb_, m_, nb_, n_,
1043  $ rsrc_
1044  parameter( block_cyclic_2d_inb = 2, dlen_ = 11,
1045  $ dtype_ = 1, ctxt_ = 2, m_ = 3, n_ = 4,
1046  $ imb_ = 5, inb_ = 6, mb_ = 7, nb_ = 8,
1047  $ rsrc_ = 9, csrc_ = 10, lld_ = 11 )
1048 * ..
1049 * .. Local Scalars ..
1050  LOGICAL GODOWN, GOLEFT
1051  INTEGER I, IACOL, IAROW, ICTXT, IIA, IJOFFA, ILOW,
1052  $ imb1, imbloc, inb1, inbloc, ioffa, ioffd, iupp,
1053  $ jja, joffa, joffd, lcmt, lcmt00, lda, ldap1,
1054  $ lmbloc, lnbloc, low, mb, mblkd, mblks, mbloc,
1055  $ mrcol, mrrow, mycol, myrow, nb, nblkd, nblks,
1056  $ nbloc, np, npcol, nprow, nq, pmb, qnb, upp
1057  COMPLEX*16 ATMP
1058 * ..
1059 * .. Local Scalars ..
1060  INTEGER DESCA2( DLEN_ )
1061 * ..
1062 * .. External Subroutines ..
1063  EXTERNAL blacs_gridinfo, pb_ainfog2l, pb_binfo,
1064  $ pb_desctrans
1065 * ..
1066 * .. Intrinsic Functions ..
1067  INTRINSIC abs, dble, dcmplx, dimag, max, min
1068 * ..
1069 * .. Executable Statements ..
1070 *
1071 * Convert descriptor
1072 *
1073  CALL pb_desctrans( desca, desca2 )
1074 *
1075 * Get grid parameters
1076 *
1077  ictxt = desca2( ctxt_ )
1078  CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
1079 *
1080  IF( n.EQ.0 )
1081  $ RETURN
1082 *
1083  CALL pb_ainfog2l( n, n, ia, ja, desca2, nprow, npcol, myrow,
1084  $ mycol, imb1, inb1, np, nq, iia, jja, iarow,
1085  $ iacol, mrrow, mrcol )
1086 *
1087 * Decide where the entries shall be stored in memory
1088 *
1089  IF( inplace ) THEN
1090  iia = 1
1091  jja = 1
1092  END IF
1093 *
1094 * Initialize LCMT00, MBLKS, NBLKS, IMBLOC, INBLOC, LMBLOC, LNBLOC,
1095 * ILOW, LOW, IUPP, and UPP.
1096 *
1097  mb = desca2( mb_ )
1098  nb = desca2( nb_ )
1099 *
1100  CALL pb_binfo( 0, np, nq, imb1, inb1, mb, nb, mrrow, mrcol,
1101  $ lcmt00, mblks, nblks, imbloc, inbloc, lmbloc,
1102  $ lnbloc, ilow, low, iupp, upp )
1103 *
1104  ioffa = iia - 1
1105  joffa = jja - 1
1106  lda = desca2( lld_ )
1107  ldap1 = lda + 1
1108 *
1109  IF( desca2( rsrc_ ).LT.0 ) THEN
1110  pmb = mb
1111  ELSE
1112  pmb = nprow * mb
1113  END IF
1114  IF( desca2( csrc_ ).LT.0 ) THEN
1115  qnb = nb
1116  ELSE
1117  qnb = npcol * nb
1118  END IF
1119 *
1120 * Handle the first block of rows or columns separately, and update
1121 * LCMT00, MBLKS and NBLKS.
1122 *
1123  godown = ( lcmt00.GT.iupp )
1124  goleft = ( lcmt00.LT.ilow )
1125 *
1126  IF( .NOT.godown .AND. .NOT.goleft ) THEN
1127 *
1128 * LCMT00 >= ILOW && LCMT00 <= IUPP
1129 *
1130  IF( lcmt00.GE.0 ) THEN
1131  ijoffa = ioffa+lcmt00 + ( joffa - 1 ) * lda
1132  DO 10 i = 1, min( inbloc, max( 0, imbloc - lcmt00 ) )
1133  atmp = a( ijoffa + i*ldap1 )
1134  a( ijoffa + i*ldap1 ) = alpha +
1135  $ dcmplx( abs( dble( atmp ) ),
1136  $ abs( dimag( atmp ) ) )
1137  10 CONTINUE
1138  ELSE
1139  ijoffa = ioffa + ( joffa - lcmt00 - 1 ) * lda
1140  DO 20 i = 1, min( imbloc, max( 0, inbloc + lcmt00 ) )
1141  atmp = a( ijoffa + i*ldap1 )
1142  a( ijoffa + i*ldap1 ) = alpha +
1143  $ dcmplx( abs( dble( atmp ) ),
1144  $ abs( dimag( atmp ) ) )
1145  20 CONTINUE
1146  END IF
1147  goleft = ( ( lcmt00 - ( iupp - upp + pmb ) ).LT.ilow )
1148  godown = .NOT.goleft
1149 *
1150  END IF
1151 *
1152  IF( godown ) THEN
1153 *
1154  lcmt00 = lcmt00 - ( iupp - upp + pmb )
1155  mblks = mblks - 1
1156  ioffa = ioffa + imbloc
1157 *
1158  30 CONTINUE
1159  IF( mblks.GT.0 .AND. lcmt00.GT.upp ) THEN
1160  lcmt00 = lcmt00 - pmb
1161  mblks = mblks - 1
1162  ioffa = ioffa + mb
1163  GO TO 30
1164  END IF
1165 *
1166  lcmt = lcmt00
1167  mblkd = mblks
1168  ioffd = ioffa
1169 *
1170  mbloc = mb
1171  40 CONTINUE
1172  IF( mblkd.GT.0 .AND. lcmt.GE.ilow ) THEN
1173  IF( mblkd.EQ.1 )
1174  $ mbloc = lmbloc
1175  IF( lcmt.GE.0 ) THEN
1176  ijoffa = ioffd + lcmt + ( joffa - 1 ) * lda
1177  DO 50 i = 1, min( inbloc, max( 0, mbloc - lcmt ) )
1178  atmp = a( ijoffa + i*ldap1 )
1179  a( ijoffa + i*ldap1 ) = alpha +
1180  $ dcmplx( abs( dble( atmp ) ),
1181  $ abs( dimag( atmp ) ) )
1182  50 CONTINUE
1183  ELSE
1184  ijoffa = ioffd + ( joffa - lcmt - 1 ) * lda
1185  DO 60 i = 1, min( mbloc, max( 0, inbloc + lcmt ) )
1186  atmp = a( ijoffa + i*ldap1 )
1187  a( ijoffa + i*ldap1 ) = alpha +
1188  $ dcmplx( abs( dble( atmp ) ),
1189  $ abs( dimag( atmp ) ) )
1190  60 CONTINUE
1191  END IF
1192  lcmt00 = lcmt
1193  lcmt = lcmt - pmb
1194  mblks = mblkd
1195  mblkd = mblkd - 1
1196  ioffa = ioffd
1197  ioffd = ioffd + mbloc
1198  GO TO 40
1199  END IF
1200 *
1201  lcmt00 = lcmt00 + low - ilow + qnb
1202  nblks = nblks - 1
1203  joffa = joffa + inbloc
1204 *
1205  ELSE IF( goleft ) THEN
1206 *
1207  lcmt00 = lcmt00 + low - ilow + qnb
1208  nblks = nblks - 1
1209  joffa = joffa + inbloc
1210 *
1211  70 CONTINUE
1212  IF( nblks.GT.0 .AND. lcmt00.LT.low ) THEN
1213  lcmt00 = lcmt00 + qnb
1214  nblks = nblks - 1
1215  joffa = joffa + nb
1216  GO TO 70
1217  END IF
1218 *
1219  lcmt = lcmt00
1220  nblkd = nblks
1221  joffd = joffa
1222 *
1223  nbloc = nb
1224  80 CONTINUE
1225  IF( nblkd.GT.0 .AND. lcmt.LE.iupp ) THEN
1226  IF( nblkd.EQ.1 )
1227  $ nbloc = lnbloc
1228  IF( lcmt.GE.0 ) THEN
1229  ijoffa = ioffa + lcmt + ( joffd - 1 ) * lda
1230  DO 90 i = 1, min( nbloc, max( 0, imbloc - lcmt ) )
1231  atmp = a( ijoffa + i*ldap1 )
1232  a( ijoffa + i*ldap1 ) = alpha +
1233  $ dcmplx( abs( dble( atmp ) ),
1234  $ abs( dimag( atmp ) ) )
1235  90 CONTINUE
1236  ELSE
1237  ijoffa = ioffa + ( joffd - lcmt - 1 ) * lda
1238  DO 100 i = 1, min( imbloc, max( 0, nbloc + lcmt ) )
1239  atmp = a( ijoffa + i*ldap1 )
1240  a( ijoffa + i*ldap1 ) = alpha +
1241  $ dcmplx( abs( dble( atmp ) ),
1242  $ abs( dimag( atmp ) ) )
1243  100 CONTINUE
1244  END IF
1245  lcmt00 = lcmt
1246  lcmt = lcmt + qnb
1247  nblks = nblkd
1248  nblkd = nblkd - 1
1249  joffa = joffd
1250  joffd = joffd + nbloc
1251  GO TO 80
1252  END IF
1253 *
1254  lcmt00 = lcmt00 - ( iupp - upp + pmb )
1255  mblks = mblks - 1
1256  ioffa = ioffa + imbloc
1257 *
1258  END IF
1259 *
1260  nbloc = nb
1261  110 CONTINUE
1262  IF( nblks.GT.0 ) THEN
1263  IF( nblks.EQ.1 )
1264  $ nbloc = lnbloc
1265  120 CONTINUE
1266  IF( mblks.GT.0 .AND. lcmt00.GT.upp ) THEN
1267  lcmt00 = lcmt00 - pmb
1268  mblks = mblks - 1
1269  ioffa = ioffa + mb
1270  GO TO 120
1271  END IF
1272 *
1273  lcmt = lcmt00
1274  mblkd = mblks
1275  ioffd = ioffa
1276 *
1277  mbloc = mb
1278  130 CONTINUE
1279  IF( mblkd.GT.0 .AND. lcmt.GE.low ) THEN
1280  IF( mblkd.EQ.1 )
1281  $ mbloc = lmbloc
1282  IF( lcmt.GE.0 ) THEN
1283  ijoffa = ioffd + lcmt + ( joffa - 1 ) * lda
1284  DO 140 i = 1, min( nbloc, max( 0, mbloc - lcmt ) )
1285  atmp = a( ijoffa + i*ldap1 )
1286  a( ijoffa + i*ldap1 ) = alpha +
1287  $ dcmplx( abs( dble( atmp ) ),
1288  $ abs( dimag( atmp ) ) )
1289  140 CONTINUE
1290  ELSE
1291  ijoffa = ioffd + ( joffa - lcmt - 1 ) * lda
1292  DO 150 i = 1, min( mbloc, max( 0, nbloc + lcmt ) )
1293  atmp = a( ijoffa + i*ldap1 )
1294  a( ijoffa + i*ldap1 ) = alpha +
1295  $ dcmplx( abs( dble( atmp ) ),
1296  $ abs( dimag( atmp ) ) )
1297  150 CONTINUE
1298  END IF
1299  lcmt00 = lcmt
1300  lcmt = lcmt - pmb
1301  mblks = mblkd
1302  mblkd = mblkd - 1
1303  ioffa = ioffd
1304  ioffd = ioffd + mbloc
1305  GO TO 130
1306  END IF
1307 *
1308  lcmt00 = lcmt00 + qnb
1309  nblks = nblks - 1
1310  joffa = joffa + nbloc
1311  GO TO 110
1312 *
1313  END IF
1314 *
1315  RETURN
1316 *
1317 * End of PZLADOM
1318 *
1319  END
1320  SUBROUTINE pb_zlascal( UPLO, M, N, IOFFD, ALPHA, A, LDA )
1322 * -- PBLAS test routine (version 2.0) --
1323 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
1324 * and University of California, Berkeley.
1325 * April 1, 1998
1326 *
1327 * .. Scalar Arguments ..
1328  CHARACTER*1 UPLO
1329  INTEGER IOFFD, LDA, M, N
1330  COMPLEX*16 ALPHA
1331 * ..
1332 * .. Array Arguments ..
1333  COMPLEX*16 A( LDA, * )
1334 * ..
1335 *
1336 * Purpose
1337 * =======
1338 *
1339 * PB_ZLASCAL scales a two-dimensional array A by the scalar alpha.
1340 *
1341 * Arguments
1342 * =========
1343 *
1344 * UPLO (input) CHARACTER*1
1345 * On entry, UPLO specifies which trapezoidal part of the ar-
1346 * ray A is to be scaled as follows:
1347 * = 'L' or 'l': the lower trapezoid of A is scaled,
1348 * = 'U' or 'u': the upper trapezoid of A is scaled,
1349 * = 'D' or 'd': diagonal specified by IOFFD is scaled,
1350 * Otherwise: all of the array A is scaled.
1351 *
1352 * M (input) INTEGER
1353 * On entry, M specifies the number of rows of the array A. M
1354 * must be at least zero.
1355 *
1356 * N (input) INTEGER
1357 * On entry, N specifies the number of columns of the array A.
1358 * N must be at least zero.
1359 *
1360 * IOFFD (input) INTEGER
1361 * On entry, IOFFD specifies the position of the offdiagonal de-
1362 * limiting the upper and lower trapezoidal part of A as follows
1363 * (see the notes below):
1364 *
1365 * IOFFD = 0 specifies the main diagonal A( i, i ),
1366 * with i = 1 ... MIN( M, N ),
1367 * IOFFD > 0 specifies the subdiagonal A( i+IOFFD, i ),
1368 * with i = 1 ... MIN( M-IOFFD, N ),
1369 * IOFFD < 0 specifies the superdiagonal A( i, i-IOFFD ),
1370 * with i = 1 ... MIN( M, N+IOFFD ).
1371 *
1372 * ALPHA (input) COMPLEX*16
1373 * On entry, ALPHA specifies the scalar alpha.
1374 *
1375 * A (input/output) COMPLEX*16 array
1376 * On entry, A is an array of dimension (LDA,N). Before entry
1377 * with UPLO = 'U' or 'u', the leading m by n part of the array
1378 * A must contain the upper trapezoidal part of the matrix as
1379 * specified by IOFFD to be scaled, and the strictly lower tra-
1380 * pezoidal part of A is not referenced; When UPLO = 'L' or 'l',
1381 * the leading m by n part of the array A must contain the lower
1382 * trapezoidal part of the matrix as specified by IOFFD to be
1383 * scaled, and the strictly upper trapezoidal part of A is not
1384 * referenced. On exit, the entries of the trapezoid part of A
1385 * determined by UPLO and IOFFD are scaled.
1386 *
1387 * LDA (input) INTEGER
1388 * On entry, LDA specifies the leading dimension of the array A.
1389 * LDA must be at least max( 1, M ).
1390 *
1391 * Notes
1392 * =====
1393 * N N
1394 * ---------------------------- -----------
1395 * | d | | |
1396 * M | d 'U' | | 'U' |
1397 * | 'L' 'D' | |d |
1398 * | d | M | d |
1399 * ---------------------------- | 'D' |
1400 * | d |
1401 * IOFFD < 0 | 'L' d |
1402 * | d|
1403 * N | |
1404 * ----------- -----------
1405 * | d 'U'|
1406 * | d | IOFFD > 0
1407 * M | 'D' |
1408 * | d| N
1409 * | 'L' | ----------------------------
1410 * | | | 'U' |
1411 * | | |d |
1412 * | | | 'D' |
1413 * | | | d |
1414 * | | |'L' d |
1415 * ----------- ----------------------------
1416 *
1417 * -- Written on April 1, 1998 by
1418 * Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
1419 *
1420 * =====================================================================
1421 *
1422 * .. Local Scalars ..
1423  INTEGER I, J, JTMP, MN
1424 * ..
1425 * .. External Functions ..
1426  LOGICAL LSAME
1427  EXTERNAL lsame
1428 * ..
1429 * .. Intrinsic Functions ..
1430  INTRINSIC max, min
1431 * ..
1432 * .. Executable Statements ..
1433 *
1434 * Quick return if possible
1435 *
1436  IF( m.LE.0 .OR. n.LE.0 )
1437  $ RETURN
1438 *
1439 * Start the operations
1440 *
1441  IF( lsame( uplo, 'L' ) ) THEN
1442 *
1443 * Scales the lower triangular part of the array by ALPHA.
1444 *
1445  mn = max( 0, -ioffd )
1446  DO 20 j = 1, min( mn, n )
1447  DO 10 i = 1, m
1448  a( i, j ) = alpha * a( i, j )
1449  10 CONTINUE
1450  20 CONTINUE
1451  DO 40 j = mn + 1, min( m - ioffd, n )
1452  DO 30 i = j + ioffd, m
1453  a( i, j ) = alpha * a( i, j )
1454  30 CONTINUE
1455  40 CONTINUE
1456 *
1457  ELSE IF( lsame( uplo, 'U' ) ) THEN
1458 *
1459 * Scales the upper triangular part of the array by ALPHA.
1460 *
1461  mn = min( m - ioffd, n )
1462  DO 60 j = max( 0, -ioffd ) + 1, mn
1463  DO 50 i = 1, j + ioffd
1464  a( i, j ) = alpha * a( i, j )
1465  50 CONTINUE
1466  60 CONTINUE
1467  DO 80 j = max( 0, mn ) + 1, n
1468  DO 70 i = 1, m
1469  a( i, j ) = alpha * a( i, j )
1470  70 CONTINUE
1471  80 CONTINUE
1472 *
1473  ELSE IF( lsame( uplo, 'D' ) ) THEN
1474 *
1475 * Scales the diagonal entries by ALPHA.
1476 *
1477  DO 90 j = max( 0, -ioffd ) + 1, min( m - ioffd, n )
1478  jtmp = j + ioffd
1479  a( jtmp, j ) = alpha * a( jtmp, j )
1480  90 CONTINUE
1481 *
1482  ELSE
1483 *
1484 * Scales the entire array by ALPHA.
1485 *
1486  DO 110 j = 1, n
1487  DO 100 i = 1, m
1488  a( i, j ) = alpha * a( i, j )
1489  100 CONTINUE
1490  110 CONTINUE
1491 *
1492  END IF
1493 *
1494  RETURN
1495 *
1496 * End of PB_ZLASCAL
1497 *
1498  END
1499  SUBROUTINE pb_zlagen( UPLO, AFORM, A, LDA, LCMT00, IRAN, MBLKS,
1500  $ IMBLOC, MB, LMBLOC, NBLKS, INBLOC, NB,
1501  $ LNBLOC, JMP, IMULADD )
1503 * -- PBLAS test routine (version 2.0) --
1504 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
1505 * and University of California, Berkeley.
1506 * April 1, 1998
1507 *
1508 * .. Scalar Arguments ..
1509  CHARACTER*1 UPLO, AFORM
1510  INTEGER IMBLOC, INBLOC, LCMT00, LDA, LMBLOC, LNBLOC,
1511  $ MB, MBLKS, NB, NBLKS
1512 * ..
1513 * .. Array Arguments ..
1514  INTEGER IMULADD( 4, * ), IRAN( * ), JMP( * )
1515  COMPLEX*16 A( LDA, * )
1516 * ..
1517 *
1518 * Purpose
1519 * =======
1520 *
1521 * PB_ZLAGEN locally initializes an array A.
1522 *
1523 * Arguments
1524 * =========
1525 *
1526 * UPLO (global input) CHARACTER*1
1527 * On entry, UPLO specifies whether the lower (UPLO='L') trape-
1528 * zoidal part or the upper (UPLO='U') trapezoidal part is to be
1529 * generated when the matrix to be generated is symmetric or
1530 * Hermitian. For all the other values of AFORM, the value of
1531 * this input argument is ignored.
1532 *
1533 * AFORM (global input) CHARACTER*1
1534 * On entry, AFORM specifies the type of submatrix to be genera-
1535 * ted as follows:
1536 * AFORM = 'S', sub( A ) is a symmetric matrix,
1537 * AFORM = 'H', sub( A ) is a Hermitian matrix,
1538 * AFORM = 'T', sub( A ) is overrwritten with the transpose
1539 * of what would normally be generated,
1540 * AFORM = 'C', sub( A ) is overwritten with the conjugate
1541 * transpose of what would normally be genera-
1542 * ted.
1543 * AFORM = 'N', a random submatrix is generated.
1544 *
1545 * A (local output) COMPLEX*16 array
1546 * On entry, A is an array of dimension (LLD_A, *). On exit,
1547 * this array contains the local entries of the randomly genera-
1548 * ted submatrix sub( A ).
1549 *
1550 * LDA (local input) INTEGER
1551 * On entry, LDA specifies the local leading dimension of the
1552 * array A. LDA must be at least one.
1553 *
1554 * LCMT00 (global input) INTEGER
1555 * On entry, LCMT00 is the LCM value specifying the off-diagonal
1556 * of the underlying matrix of interest. LCMT00=0 specifies the
1557 * main diagonal, LCMT00 > 0 specifies a subdiagonal, LCMT00 < 0
1558 * specifies superdiagonals.
1559 *
1560 * IRAN (local input) INTEGER array
1561 * On entry, IRAN is an array of dimension 2 containing respec-
1562 * tively the 16-lower and 16-higher bits of the encoding of the
1563 * entry of the random sequence corresponding locally to the
1564 * first local array entry to generate. Usually, this array is
1565 * computed by PB_SETLOCRAN.
1566 *
1567 * MBLKS (local input) INTEGER
1568 * On entry, MBLKS specifies the local number of blocks of rows.
1569 * MBLKS is at least zero.
1570 *
1571 * IMBLOC (local input) INTEGER
1572 * On entry, IMBLOC specifies the number of rows (size) of the
1573 * local uppest blocks. IMBLOC is at least zero.
1574 *
1575 * MB (global input) INTEGER
1576 * On entry, MB specifies the blocking factor used to partition
1577 * the rows of the matrix. MB must be at least one.
1578 *
1579 * LMBLOC (local input) INTEGER
1580 * On entry, LMBLOC specifies the number of rows (size) of the
1581 * local lowest blocks. LMBLOC is at least zero.
1582 *
1583 * NBLKS (local input) INTEGER
1584 * On entry, NBLKS specifies the local number of blocks of co-
1585 * lumns. NBLKS is at least zero.
1586 *
1587 * INBLOC (local input) INTEGER
1588 * On entry, INBLOC specifies the number of columns (size) of
1589 * the local leftmost blocks. INBLOC is at least zero.
1590 *
1591 * NB (global input) INTEGER
1592 * On entry, NB specifies the blocking factor used to partition
1593 * the the columns of the matrix. NB must be at least one.
1594 *
1595 * LNBLOC (local input) INTEGER
1596 * On entry, LNBLOC specifies the number of columns (size) of
1597 * the local rightmost blocks. LNBLOC is at least zero.
1598 *
1599 * JMP (local input) INTEGER array
1600 * On entry, JMP is an array of dimension JMP_LEN containing the
1601 * different jump values used by the random matrix generator.
1602 *
1603 * IMULADD (local input) INTEGER array
1604 * On entry, IMULADD is an array of dimension (4, JMP_LEN). The
1605 * jth column of this array contains the encoded initial cons-
1606 * tants a_j and c_j to jump from X( n ) to X( n + JMP( j ) )
1607 * (= a_j * X( n ) + c_j) in the random sequence. IMULADD(1:2,j)
1608 * contains respectively the 16-lower and 16-higher bits of the
1609 * constant a_j, and IMULADD(3:4,j) contains the 16-lower and
1610 * 16-higher bits of the constant c_j.
1611 *
1612 * -- Written on April 1, 1998 by
1613 * Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
1614 *
1615 * =====================================================================
1616 *
1617 * .. Parameters ..
1618  INTEGER JMP_1, JMP_COL, JMP_IMBV, JMP_INBV, JMP_LEN,
1619  $ JMP_MB, JMP_NB, JMP_NPIMBLOC, JMP_NPMB,
1620  $ JMP_NQINBLOC, JMP_NQNB, JMP_ROW
1621  parameter( jmp_1 = 1, jmp_row = 2, jmp_col = 3,
1622  $ jmp_mb = 4, jmp_imbv = 5, jmp_npmb = 6,
1623  $ jmp_npimbloc = 7, jmp_nb = 8, jmp_inbv = 9,
1624  $ jmp_nqnb = 10, jmp_nqinbloc = 11,
1625  $ jmp_len = 11 )
1626  DOUBLE PRECISION ZERO
1627  PARAMETER ( ZERO = 0.0d+0 )
1628 * ..
1629 * .. Local Scalars ..
1630  INTEGER I, IB, IBLK, II, IK, ITMP, JB, JBLK, JJ, JK,
1631  $ JTMP, LCMTC, LCMTR, LOW, MNB, UPP
1632  COMPLEX*16 DUMMY
1633 * ..
1634 * .. Local Arrays ..
1635  INTEGER IB0( 2 ), IB1( 2 ), IB2( 2 ), IB3( 2 )
1636 * ..
1637 * .. External Subroutines ..
1638  EXTERNAL pb_jumpit
1639 * ..
1640 * .. External Functions ..
1641  LOGICAL LSAME
1642  DOUBLE PRECISION PB_DRAND
1643  EXTERNAL LSAME, PB_DRAND
1644 * ..
1645 * .. Intrinsic Functions ..
1646  INTRINSIC dble, dcmplx, max, min
1647 * ..
1648 * .. Executable Statements ..
1649 *
1650  DO 10 i = 1, 2
1651  ib1( i ) = iran( i )
1652  ib2( i ) = iran( i )
1653  ib3( i ) = iran( i )
1654  10 CONTINUE
1655 *
1656  IF( lsame( aform, 'N' ) ) THEN
1657 *
1658 * Generate random matrix
1659 *
1660  jj = 1
1661 *
1662  DO 50 jblk = 1, nblks
1663 *
1664  IF( jblk.EQ.1 ) THEN
1665  jb = inbloc
1666  ELSE IF( jblk.EQ.nblks ) THEN
1667  jb = lnbloc
1668  ELSE
1669  jb = nb
1670  END IF
1671 *
1672  DO 40 jk = jj, jj + jb - 1
1673 *
1674  ii = 1
1675 *
1676  DO 30 iblk = 1, mblks
1677 *
1678  IF( iblk.EQ.1 ) THEN
1679  ib = imbloc
1680  ELSE IF( iblk.EQ.mblks ) THEN
1681  ib = lmbloc
1682  ELSE
1683  ib = mb
1684  END IF
1685 *
1686 * Blocks are IB by JB
1687 *
1688  DO 20 ik = ii, ii + ib - 1
1689  a( ik, jk ) = dcmplx( pb_drand( 0 ),
1690  $ pb_drand( 0 ) )
1691  20 CONTINUE
1692 *
1693  ii = ii + ib
1694 *
1695  IF( iblk.EQ.1 ) THEN
1696 *
1697 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
1698 *
1699  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib1,
1700  $ ib0 )
1701 *
1702  ELSE
1703 *
1704 * Jump NPROW * MB rows
1705 *
1706  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib1, ib0 )
1707 *
1708  END IF
1709 *
1710  ib1( 1 ) = ib0( 1 )
1711  ib1( 2 ) = ib0( 2 )
1712 *
1713  30 CONTINUE
1714 *
1715 * Jump one column
1716 *
1717  CALL pb_jumpit( imuladd( 1, jmp_col ), ib2, ib0 )
1718 *
1719  ib1( 1 ) = ib0( 1 )
1720  ib1( 2 ) = ib0( 2 )
1721  ib2( 1 ) = ib0( 1 )
1722  ib2( 2 ) = ib0( 2 )
1723 *
1724  40 CONTINUE
1725 *
1726  jj = jj + jb
1727 *
1728  IF( jblk.EQ.1 ) THEN
1729 *
1730 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
1731 *
1732  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib3, ib0 )
1733 *
1734  ELSE
1735 *
1736 * Jump NPCOL * NB columns
1737 *
1738  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib3, ib0 )
1739 *
1740  END IF
1741 *
1742  ib1( 1 ) = ib0( 1 )
1743  ib1( 2 ) = ib0( 2 )
1744  ib2( 1 ) = ib0( 1 )
1745  ib2( 2 ) = ib0( 2 )
1746  ib3( 1 ) = ib0( 1 )
1747  ib3( 2 ) = ib0( 2 )
1748 *
1749  50 CONTINUE
1750 *
1751  ELSE IF( lsame( aform, 'T' ) ) THEN
1752 *
1753 * Generate the transpose of the matrix that would be normally
1754 * generated.
1755 *
1756  ii = 1
1757 *
1758  DO 90 iblk = 1, mblks
1759 *
1760  IF( iblk.EQ.1 ) THEN
1761  ib = imbloc
1762  ELSE IF( iblk.EQ.mblks ) THEN
1763  ib = lmbloc
1764  ELSE
1765  ib = mb
1766  END IF
1767 *
1768  DO 80 ik = ii, ii + ib - 1
1769 *
1770  jj = 1
1771 *
1772  DO 70 jblk = 1, nblks
1773 *
1774  IF( jblk.EQ.1 ) THEN
1775  jb = inbloc
1776  ELSE IF( jblk.EQ.nblks ) THEN
1777  jb = lnbloc
1778  ELSE
1779  jb = nb
1780  END IF
1781 *
1782 * Blocks are IB by JB
1783 *
1784  DO 60 jk = jj, jj + jb - 1
1785  a( ik, jk ) = dcmplx( pb_drand( 0 ),
1786  $ pb_drand( 0 ) )
1787  60 CONTINUE
1788 *
1789  jj = jj + jb
1790 *
1791  IF( jblk.EQ.1 ) THEN
1792 *
1793 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
1794 *
1795  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib1,
1796  $ ib0 )
1797 *
1798  ELSE
1799 *
1800 * Jump NPCOL * NB columns
1801 *
1802  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib1, ib0 )
1803 *
1804  END IF
1805 *
1806  ib1( 1 ) = ib0( 1 )
1807  ib1( 2 ) = ib0( 2 )
1808 *
1809  70 CONTINUE
1810 *
1811 * Jump one row
1812 *
1813  CALL pb_jumpit( imuladd( 1, jmp_row ), ib2, ib0 )
1814 *
1815  ib1( 1 ) = ib0( 1 )
1816  ib1( 2 ) = ib0( 2 )
1817  ib2( 1 ) = ib0( 1 )
1818  ib2( 2 ) = ib0( 2 )
1819 *
1820  80 CONTINUE
1821 *
1822  ii = ii + ib
1823 *
1824  IF( iblk.EQ.1 ) THEN
1825 *
1826 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
1827 *
1828  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib3, ib0 )
1829 *
1830  ELSE
1831 *
1832 * Jump NPROW * MB rows
1833 *
1834  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib3, ib0 )
1835 *
1836  END IF
1837 *
1838  ib1( 1 ) = ib0( 1 )
1839  ib1( 2 ) = ib0( 2 )
1840  ib2( 1 ) = ib0( 1 )
1841  ib2( 2 ) = ib0( 2 )
1842  ib3( 1 ) = ib0( 1 )
1843  ib3( 2 ) = ib0( 2 )
1844 *
1845  90 CONTINUE
1846 *
1847  ELSE IF( lsame( aform, 'S' ) ) THEN
1848 *
1849 * Generate a symmetric matrix
1850 *
1851  IF( lsame( uplo, 'L' ) ) THEN
1852 *
1853 * generate lower trapezoidal part
1854 *
1855  jj = 1
1856  lcmtc = lcmt00
1857 *
1858  DO 170 jblk = 1, nblks
1859 *
1860  IF( jblk.EQ.1 ) THEN
1861  jb = inbloc
1862  low = 1 - inbloc
1863  ELSE IF( jblk.EQ.nblks ) THEN
1864  jb = lnbloc
1865  low = 1 - nb
1866  ELSE
1867  jb = nb
1868  low = 1 - nb
1869  END IF
1870 *
1871  DO 160 jk = jj, jj + jb - 1
1872 *
1873  ii = 1
1874  lcmtr = lcmtc
1875 *
1876  DO 150 iblk = 1, mblks
1877 *
1878  IF( iblk.EQ.1 ) THEN
1879  ib = imbloc
1880  upp = imbloc - 1
1881  ELSE IF( iblk.EQ.mblks ) THEN
1882  ib = lmbloc
1883  upp = mb - 1
1884  ELSE
1885  ib = mb
1886  upp = mb - 1
1887  END IF
1888 *
1889 * Blocks are IB by JB
1890 *
1891  IF( lcmtr.GT.upp ) THEN
1892 *
1893  DO 100 ik = ii, ii + ib - 1
1894  dummy = dcmplx( pb_drand( 0 ),
1895  $ pb_drand( 0 ) )
1896  100 CONTINUE
1897 *
1898  ELSE IF( lcmtr.GE.low ) THEN
1899 *
1900  jtmp = jk - jj + 1
1901  mnb = max( 0, -lcmtr )
1902 *
1903  IF( jtmp.LE.min( mnb, jb ) ) THEN
1904 *
1905  DO 110 ik = ii, ii + ib - 1
1906  a( ik, jk ) = dcmplx( pb_drand( 0 ),
1907  $ pb_drand( 0 ) )
1908  110 CONTINUE
1909 *
1910  ELSE IF( ( jtmp.GE.( mnb + 1 ) ) .AND.
1911  $ ( jtmp.LE.min( ib-lcmtr, jb ) ) ) THEN
1912 *
1913  itmp = ii + jtmp + lcmtr - 1
1914 *
1915  DO 120 ik = ii, itmp - 1
1916  dummy = dcmplx( pb_drand( 0 ),
1917  $ pb_drand( 0 ) )
1918  120 CONTINUE
1919 *
1920  DO 130 ik = itmp, ii + ib - 1
1921  a( ik, jk ) = dcmplx( pb_drand( 0 ),
1922  $ pb_drand( 0 ) )
1923  130 CONTINUE
1924 *
1925  END IF
1926 *
1927  ELSE
1928 *
1929  DO 140 ik = ii, ii + ib - 1
1930  a( ik, jk ) = dcmplx( pb_drand( 0 ),
1931  $ pb_drand( 0 ) )
1932  140 CONTINUE
1933 *
1934  END IF
1935 *
1936  ii = ii + ib
1937 *
1938  IF( iblk.EQ.1 ) THEN
1939 *
1940 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
1941 *
1942  lcmtr = lcmtr - jmp( jmp_npimbloc )
1943  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib1,
1944  $ ib0 )
1945 *
1946  ELSE
1947 *
1948 * Jump NPROW * MB rows
1949 *
1950  lcmtr = lcmtr - jmp( jmp_npmb )
1951  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib1,
1952  $ ib0 )
1953 *
1954  END IF
1955 *
1956  ib1( 1 ) = ib0( 1 )
1957  ib1( 2 ) = ib0( 2 )
1958 *
1959  150 CONTINUE
1960 *
1961 * Jump one column
1962 *
1963  CALL pb_jumpit( imuladd( 1, jmp_col ), ib2, ib0 )
1964 *
1965  ib1( 1 ) = ib0( 1 )
1966  ib1( 2 ) = ib0( 2 )
1967  ib2( 1 ) = ib0( 1 )
1968  ib2( 2 ) = ib0( 2 )
1969 *
1970  160 CONTINUE
1971 *
1972  jj = jj + jb
1973 *
1974  IF( jblk.EQ.1 ) THEN
1975 *
1976 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
1977 *
1978  lcmtc = lcmtc + jmp( jmp_nqinbloc )
1979  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib3, ib0 )
1980 *
1981  ELSE
1982 *
1983 * Jump NPCOL * NB columns
1984 *
1985  lcmtc = lcmtc + jmp( jmp_nqnb )
1986  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib3, ib0 )
1987 *
1988  END IF
1989 *
1990  ib1( 1 ) = ib0( 1 )
1991  ib1( 2 ) = ib0( 2 )
1992  ib2( 1 ) = ib0( 1 )
1993  ib2( 2 ) = ib0( 2 )
1994  ib3( 1 ) = ib0( 1 )
1995  ib3( 2 ) = ib0( 2 )
1996 *
1997  170 CONTINUE
1998 *
1999  ELSE
2000 *
2001 * generate upper trapezoidal part
2002 *
2003  ii = 1
2004  lcmtr = lcmt00
2005 *
2006  DO 250 iblk = 1, mblks
2007 *
2008  IF( iblk.EQ.1 ) THEN
2009  ib = imbloc
2010  upp = imbloc - 1
2011  ELSE IF( iblk.EQ.mblks ) THEN
2012  ib = lmbloc
2013  upp = mb - 1
2014  ELSE
2015  ib = mb
2016  upp = mb - 1
2017  END IF
2018 *
2019  DO 240 ik = ii, ii + ib - 1
2020 *
2021  jj = 1
2022  lcmtc = lcmtr
2023 *
2024  DO 230 jblk = 1, nblks
2025 *
2026  IF( jblk.EQ.1 ) THEN
2027  jb = inbloc
2028  low = 1 - inbloc
2029  ELSE IF( jblk.EQ.nblks ) THEN
2030  jb = lnbloc
2031  low = 1 - nb
2032  ELSE
2033  jb = nb
2034  low = 1 - nb
2035  END IF
2036 *
2037 * Blocks are IB by JB
2038 *
2039  IF( lcmtc.LT.low ) THEN
2040 *
2041  DO 180 jk = jj, jj + jb - 1
2042  dummy = dcmplx( pb_drand( 0 ),
2043  $ pb_drand( 0 ) )
2044  180 CONTINUE
2045 *
2046  ELSE IF( lcmtc.LE.upp ) THEN
2047 *
2048  itmp = ik - ii + 1
2049  mnb = max( 0, lcmtc )
2050 *
2051  IF( itmp.LE.min( mnb, ib ) ) THEN
2052 *
2053  DO 190 jk = jj, jj + jb - 1
2054  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2055  $ pb_drand( 0 ) )
2056  190 CONTINUE
2057 *
2058  ELSE IF( ( itmp.GE.( mnb + 1 ) ) .AND.
2059  $ ( itmp.LE.min( jb+lcmtc, ib ) ) ) THEN
2060 *
2061  jtmp = jj + itmp - lcmtc - 1
2062 *
2063  DO 200 jk = jj, jtmp - 1
2064  dummy = dcmplx( pb_drand( 0 ),
2065  $ pb_drand( 0 ) )
2066  200 CONTINUE
2067 *
2068  DO 210 jk = jtmp, jj + jb - 1
2069  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2070  $ pb_drand( 0 ) )
2071  210 CONTINUE
2072 *
2073  END IF
2074 *
2075  ELSE
2076 *
2077  DO 220 jk = jj, jj + jb - 1
2078  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2079  $ pb_drand( 0 ) )
2080  220 CONTINUE
2081 *
2082  END IF
2083 *
2084  jj = jj + jb
2085 *
2086  IF( jblk.EQ.1 ) THEN
2087 *
2088 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
2089 *
2090  lcmtc = lcmtc + jmp( jmp_nqinbloc )
2091  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib1,
2092  $ ib0 )
2093 *
2094  ELSE
2095 *
2096 * Jump NPCOL * NB columns
2097 *
2098  lcmtc = lcmtc + jmp( jmp_nqnb )
2099  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib1,
2100  $ ib0 )
2101 *
2102  END IF
2103 *
2104  ib1( 1 ) = ib0( 1 )
2105  ib1( 2 ) = ib0( 2 )
2106 *
2107  230 CONTINUE
2108 *
2109 * Jump one row
2110 *
2111  CALL pb_jumpit( imuladd( 1, jmp_row ), ib2, ib0 )
2112 *
2113  ib1( 1 ) = ib0( 1 )
2114  ib1( 2 ) = ib0( 2 )
2115  ib2( 1 ) = ib0( 1 )
2116  ib2( 2 ) = ib0( 2 )
2117 *
2118  240 CONTINUE
2119 *
2120  ii = ii + ib
2121 *
2122  IF( iblk.EQ.1 ) THEN
2123 *
2124 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
2125 *
2126  lcmtr = lcmtr - jmp( jmp_npimbloc )
2127  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib3, ib0 )
2128 *
2129  ELSE
2130 *
2131 * Jump NPROW * MB rows
2132 *
2133  lcmtr = lcmtr - jmp( jmp_npmb )
2134  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib3, ib0 )
2135 *
2136  END IF
2137 *
2138  ib1( 1 ) = ib0( 1 )
2139  ib1( 2 ) = ib0( 2 )
2140  ib2( 1 ) = ib0( 1 )
2141  ib2( 2 ) = ib0( 2 )
2142  ib3( 1 ) = ib0( 1 )
2143  ib3( 2 ) = ib0( 2 )
2144 *
2145  250 CONTINUE
2146 *
2147  END IF
2148 *
2149  ELSE IF( lsame( aform, 'C' ) ) THEN
2150 *
2151 * Generate the conjugate transpose of the matrix that would be
2152 * normally generated.
2153 *
2154  ii = 1
2155 *
2156  DO 290 iblk = 1, mblks
2157 *
2158  IF( iblk.EQ.1 ) THEN
2159  ib = imbloc
2160  ELSE IF( iblk.EQ.mblks ) THEN
2161  ib = lmbloc
2162  ELSE
2163  ib = mb
2164  END IF
2165 *
2166  DO 280 ik = ii, ii + ib - 1
2167 *
2168  jj = 1
2169 *
2170  DO 270 jblk = 1, nblks
2171 *
2172  IF( jblk.EQ.1 ) THEN
2173  jb = inbloc
2174  ELSE IF( jblk.EQ.nblks ) THEN
2175  jb = lnbloc
2176  ELSE
2177  jb = nb
2178  END IF
2179 *
2180 * Blocks are IB by JB
2181 *
2182  DO 260 jk = jj, jj + jb - 1
2183  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2184  $ -pb_drand( 0 ) )
2185  260 CONTINUE
2186 *
2187  jj = jj + jb
2188 *
2189  IF( jblk.EQ.1 ) THEN
2190 *
2191 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
2192 *
2193  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib1,
2194  $ ib0 )
2195 *
2196  ELSE
2197 *
2198 * Jump NPCOL * NB columns
2199 *
2200  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib1,
2201  $ ib0 )
2202 *
2203  END IF
2204 *
2205  ib1( 1 ) = ib0( 1 )
2206  ib1( 2 ) = ib0( 2 )
2207 *
2208  270 CONTINUE
2209 *
2210 * Jump one row
2211 *
2212  CALL pb_jumpit( imuladd( 1, jmp_row ), ib2, ib0 )
2213 *
2214  ib1( 1 ) = ib0( 1 )
2215  ib1( 2 ) = ib0( 2 )
2216  ib2( 1 ) = ib0( 1 )
2217  ib2( 2 ) = ib0( 2 )
2218 *
2219  280 CONTINUE
2220 *
2221  ii = ii + ib
2222 *
2223  IF( iblk.EQ.1 ) THEN
2224 *
2225 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
2226 *
2227  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib3, ib0 )
2228 *
2229  ELSE
2230 *
2231 * Jump NPROW * MB rows
2232 *
2233  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib3, ib0 )
2234 *
2235  END IF
2236 *
2237  ib1( 1 ) = ib0( 1 )
2238  ib1( 2 ) = ib0( 2 )
2239  ib2( 1 ) = ib0( 1 )
2240  ib2( 2 ) = ib0( 2 )
2241  ib3( 1 ) = ib0( 1 )
2242  ib3( 2 ) = ib0( 2 )
2243 *
2244  290 CONTINUE
2245 *
2246  ELSE IF( lsame( aform, 'H' ) ) THEN
2247 *
2248 * Generate a Hermitian matrix
2249 *
2250  IF( lsame( uplo, 'L' ) ) THEN
2251 *
2252 * generate lower trapezoidal part
2253 *
2254  jj = 1
2255  lcmtc = lcmt00
2256 *
2257  DO 370 jblk = 1, nblks
2258 *
2259  IF( jblk.EQ.1 ) THEN
2260  jb = inbloc
2261  low = 1 - inbloc
2262  ELSE IF( jblk.EQ.nblks ) THEN
2263  jb = lnbloc
2264  low = 1 - nb
2265  ELSE
2266  jb = nb
2267  low = 1 - nb
2268  END IF
2269 *
2270  DO 360 jk = jj, jj + jb - 1
2271 *
2272  ii = 1
2273  lcmtr = lcmtc
2274 *
2275  DO 350 iblk = 1, mblks
2276 *
2277  IF( iblk.EQ.1 ) THEN
2278  ib = imbloc
2279  upp = imbloc - 1
2280  ELSE IF( iblk.EQ.mblks ) THEN
2281  ib = lmbloc
2282  upp = mb - 1
2283  ELSE
2284  ib = mb
2285  upp = mb - 1
2286  END IF
2287 *
2288 * Blocks are IB by JB
2289 *
2290  IF( lcmtr.GT.upp ) THEN
2291 *
2292  DO 300 ik = ii, ii + ib - 1
2293  dummy = dcmplx( pb_drand( 0 ),
2294  $ pb_drand( 0 ) )
2295  300 CONTINUE
2296 *
2297  ELSE IF( lcmtr.GE.low ) THEN
2298 *
2299  jtmp = jk - jj + 1
2300  mnb = max( 0, -lcmtr )
2301 *
2302  IF( jtmp.LE.min( mnb, jb ) ) THEN
2303 *
2304  DO 310 ik = ii, ii + ib - 1
2305  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2306  $ pb_drand( 0 ) )
2307  310 CONTINUE
2308 *
2309  ELSE IF( ( jtmp.GE.( mnb + 1 ) ) .AND.
2310  $ ( jtmp.LE.min( ib-lcmtr, jb ) ) ) THEN
2311 *
2312  itmp = ii + jtmp + lcmtr - 1
2313 *
2314  DO 320 ik = ii, itmp - 1
2315  dummy = dcmplx( pb_drand( 0 ),
2316  $ pb_drand( 0 ) )
2317  320 CONTINUE
2318 *
2319  IF( itmp.LE.( ii + ib - 1 ) ) THEN
2320  dummy = dcmplx( pb_drand( 0 ),
2321  $ -pb_drand( 0 ) )
2322  a( itmp, jk ) = dcmplx( dble( dummy ),
2323  $ zero )
2324  END IF
2325 *
2326  DO 330 ik = itmp + 1, ii + ib - 1
2327  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2328  $ pb_drand( 0 ) )
2329  330 CONTINUE
2330 *
2331  END IF
2332 *
2333  ELSE
2334 *
2335  DO 340 ik = ii, ii + ib - 1
2336  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2337  $ pb_drand( 0 ) )
2338  340 CONTINUE
2339 *
2340  END IF
2341 *
2342  ii = ii + ib
2343 *
2344  IF( iblk.EQ.1 ) THEN
2345 *
2346 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
2347 *
2348  lcmtr = lcmtr - jmp( jmp_npimbloc )
2349  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib1,
2350  $ ib0 )
2351 *
2352  ELSE
2353 *
2354 * Jump NPROW * MB rows
2355 *
2356  lcmtr = lcmtr - jmp( jmp_npmb )
2357  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib1,
2358  $ ib0 )
2359 *
2360  END IF
2361 *
2362  ib1( 1 ) = ib0( 1 )
2363  ib1( 2 ) = ib0( 2 )
2364 *
2365  350 CONTINUE
2366 *
2367 * Jump one column
2368 *
2369  CALL pb_jumpit( imuladd( 1, jmp_col ), ib2, ib0 )
2370 *
2371  ib1( 1 ) = ib0( 1 )
2372  ib1( 2 ) = ib0( 2 )
2373  ib2( 1 ) = ib0( 1 )
2374  ib2( 2 ) = ib0( 2 )
2375 *
2376  360 CONTINUE
2377 *
2378  jj = jj + jb
2379 *
2380  IF( jblk.EQ.1 ) THEN
2381 *
2382 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
2383 *
2384  lcmtc = lcmtc + jmp( jmp_nqinbloc )
2385  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib3, ib0 )
2386 *
2387  ELSE
2388 *
2389 * Jump NPCOL * NB columns
2390 *
2391  lcmtc = lcmtc + jmp( jmp_nqnb )
2392  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib3, ib0 )
2393 *
2394  END IF
2395 *
2396  ib1( 1 ) = ib0( 1 )
2397  ib1( 2 ) = ib0( 2 )
2398  ib2( 1 ) = ib0( 1 )
2399  ib2( 2 ) = ib0( 2 )
2400  ib3( 1 ) = ib0( 1 )
2401  ib3( 2 ) = ib0( 2 )
2402 *
2403  370 CONTINUE
2404 *
2405  ELSE
2406 *
2407 * generate upper trapezoidal part
2408 *
2409  ii = 1
2410  lcmtr = lcmt00
2411 *
2412  DO 450 iblk = 1, mblks
2413 *
2414  IF( iblk.EQ.1 ) THEN
2415  ib = imbloc
2416  upp = imbloc - 1
2417  ELSE IF( iblk.EQ.mblks ) THEN
2418  ib = lmbloc
2419  upp = mb - 1
2420  ELSE
2421  ib = mb
2422  upp = mb - 1
2423  END IF
2424 *
2425  DO 440 ik = ii, ii + ib - 1
2426 *
2427  jj = 1
2428  lcmtc = lcmtr
2429 *
2430  DO 430 jblk = 1, nblks
2431 *
2432  IF( jblk.EQ.1 ) THEN
2433  jb = inbloc
2434  low = 1 - inbloc
2435  ELSE IF( jblk.EQ.nblks ) THEN
2436  jb = lnbloc
2437  low = 1 - nb
2438  ELSE
2439  jb = nb
2440  low = 1 - nb
2441  END IF
2442 *
2443 * Blocks are IB by JB
2444 *
2445  IF( lcmtc.LT.low ) THEN
2446 *
2447  DO 380 jk = jj, jj + jb - 1
2448  dummy = dcmplx( pb_drand( 0 ),
2449  $ -pb_drand( 0 ) )
2450  380 CONTINUE
2451 *
2452  ELSE IF( lcmtc.LE.upp ) THEN
2453 *
2454  itmp = ik - ii + 1
2455  mnb = max( 0, lcmtc )
2456 *
2457  IF( itmp.LE.min( mnb, ib ) ) THEN
2458 *
2459  DO 390 jk = jj, jj + jb - 1
2460  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2461  $ -pb_drand( 0 ) )
2462  390 CONTINUE
2463 *
2464  ELSE IF( ( itmp.GE.( mnb + 1 ) ) .AND.
2465  $ ( itmp.LE.min( jb+lcmtc, ib ) ) ) THEN
2466 *
2467  jtmp = jj + itmp - lcmtc - 1
2468 *
2469  DO 400 jk = jj, jtmp - 1
2470  dummy = dcmplx( pb_drand( 0 ),
2471  $ -pb_drand( 0 ) )
2472  400 CONTINUE
2473 *
2474  IF( jtmp.LE.( jj + jb - 1 ) ) THEN
2475  dummy = dcmplx( pb_drand( 0 ),
2476  $ -pb_drand( 0 ) )
2477  a( ik, jtmp ) = dcmplx( dble( dummy ),
2478  $ zero )
2479  END IF
2480 *
2481  DO 410 jk = jtmp + 1, jj + jb - 1
2482  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2483  $ -pb_drand( 0 ) )
2484  410 CONTINUE
2485 *
2486  END IF
2487 *
2488  ELSE
2489 *
2490  DO 420 jk = jj, jj + jb - 1
2491  a( ik, jk ) = dcmplx( pb_drand( 0 ),
2492  $ -pb_drand( 0 ) )
2493  420 CONTINUE
2494 *
2495  END IF
2496 *
2497  jj = jj + jb
2498 *
2499  IF( jblk.EQ.1 ) THEN
2500 *
2501 * Jump INBLOC + ( NPCOL - 1 ) * NB columns
2502 *
2503  lcmtc = lcmtc + jmp( jmp_nqinbloc )
2504  CALL pb_jumpit( imuladd( 1, jmp_nqinbloc ), ib1,
2505  $ ib0 )
2506 *
2507  ELSE
2508 *
2509 * Jump NPCOL * NB columns
2510 *
2511  lcmtc = lcmtc + jmp( jmp_nqnb )
2512  CALL pb_jumpit( imuladd( 1, jmp_nqnb ), ib1,
2513  $ ib0 )
2514 *
2515  END IF
2516 *
2517  ib1( 1 ) = ib0( 1 )
2518  ib1( 2 ) = ib0( 2 )
2519 *
2520  430 CONTINUE
2521 *
2522 * Jump one row
2523 *
2524  CALL pb_jumpit( imuladd( 1, jmp_row ), ib2, ib0 )
2525 *
2526  ib1( 1 ) = ib0( 1 )
2527  ib1( 2 ) = ib0( 2 )
2528  ib2( 1 ) = ib0( 1 )
2529  ib2( 2 ) = ib0( 2 )
2530 *
2531  440 CONTINUE
2532 *
2533  ii = ii + ib
2534 *
2535  IF( iblk.EQ.1 ) THEN
2536 *
2537 * Jump IMBLOC + ( NPROW - 1 ) * MB rows
2538 *
2539  lcmtr = lcmtr - jmp( jmp_npimbloc )
2540  CALL pb_jumpit( imuladd( 1, jmp_npimbloc ), ib3, ib0 )
2541 *
2542  ELSE
2543 *
2544 * Jump NPROW * MB rows
2545 *
2546  lcmtr = lcmtr - jmp( jmp_npmb )
2547  CALL pb_jumpit( imuladd( 1, jmp_npmb ), ib3, ib0 )
2548 *
2549  END IF
2550 *
2551  ib1( 1 ) = ib0( 1 )
2552  ib1( 2 ) = ib0( 2 )
2553  ib2( 1 ) = ib0( 1 )
2554  ib2( 2 ) = ib0( 2 )
2555  ib3( 1 ) = ib0( 1 )
2556  ib3( 2 ) = ib0( 2 )
2557 *
2558  450 CONTINUE
2559 *
2560  END IF
2561 *
2562  END IF
2563 *
2564  RETURN
2565 *
2566 * End of PB_ZLAGEN
2567 *
2568  END
2569  DOUBLE PRECISION FUNCTION pb_drand( IDUMM )
2571 * -- PBLAS test routine (version 2.0) --
2572 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
2573 * and University of California, Berkeley.
2574 * April 1, 1998
2575 *
2576 * .. Scalar Arguments ..
2577  INTEGER idumm
2578 * ..
2579 *
2580 * Purpose
2581 * =======
2582 *
2583 * PB_DRAND generates the next number in the random sequence. This func-
2584 * tion ensures that this number will be in the interval ( -1.0, 1.0 ).
2585 *
2586 * Arguments
2587 * =========
2588 *
2589 * IDUMM (local input) INTEGER
2590 * This argument is ignored, but necessary to a FORTRAN 77 func-
2591 * tion.
2592 *
2593 * Further Details
2594 * ===============
2595 *
2596 * On entry, the array IRAND stored in the common block RANCOM contains
2597 * the information (2 integers) required to generate the next number in
2598 * the sequence X( n ). This number is computed as
2599 *
2600 * X( n ) = ( 2^16 * IRAND( 2 ) + IRAND( 1 ) ) / d,
2601 *
2602 * where the constant d is the largest 32 bit positive integer. The
2603 * array IRAND is then updated for the generation of the next number
2604 * X( n+1 ) in the random sequence as follows X( n+1 ) = a * X( n ) + c.
2605 * The constants a and c should have been preliminarily stored in the
2606 * array IACS as 2 pairs of integers. The initial set up of IRAND and
2607 * IACS is performed by the routine PB_SETRAN.
2608 *
2609 * -- Written on April 1, 1998 by
2610 * Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
2611 *
2612 * =====================================================================
2613 *
2614 * .. Parameters ..
2615  DOUBLE PRECISION one, two
2616  PARAMETER ( one = 1.0d+0, two = 2.0d+0 )
2617 * ..
2618 * .. External Functions ..
2619  DOUBLE PRECISION pb_dran
2620  EXTERNAL pb_dran
2621 * ..
2622 * .. Executable Statements ..
2623 *
2624  pb_drand = one - two * pb_dran( idumm )
2625 *
2626  RETURN
2627 *
2628 * End of PB_DRAND
2629 *
2630  END
2631  DOUBLE PRECISION FUNCTION pb_dran( IDUMM )
2633 * -- PBLAS test routine (version 2.0) --
2634 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
2635 * and University of California, Berkeley.
2636 * April 1, 1998
2637 *
2638 * .. Scalar Arguments ..
2639  INTEGER idumm
2640 * ..
2641 *
2642 * Purpose
2643 * =======
2644 *
2645 * PB_DRAN generates the next number in the random sequence.
2646 *
2647 * Arguments
2648 * =========
2649 *
2650 * IDUMM (local input) INTEGER
2651 * This argument is ignored, but necessary to a FORTRAN 77 func-
2652 * tion.
2653 *
2654 * Further Details
2655 * ===============
2656 *
2657 * On entry, the array IRAND stored in the common block RANCOM contains
2658 * the information (2 integers) required to generate the next number in
2659 * the sequence X( n ). This number is computed as
2660 *
2661 * X( n ) = ( 2^16 * IRAND( 2 ) + IRAND( 1 ) ) / d,
2662 *
2663 * where the constant d is the largest 32 bit positive integer. The
2664 * array IRAND is then updated for the generation of the next number
2665 * X( n+1 ) in the random sequence as follows X( n+1 ) = a * X( n ) + c.
2666 * The constants a and c should have been preliminarily stored in the
2667 * array IACS as 2 pairs of integers. The initial set up of IRAND and
2668 * IACS is performed by the routine PB_SETRAN.
2669 *
2670 * -- Written on April 1, 1998 by
2671 * Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
2672 *
2673 * =====================================================================
2674 *
2675 * .. Parameters ..
2676  DOUBLE PRECISION divfac, pow16
2677  PARAMETER ( divfac = 2.147483648d+9,
2678  $ pow16 = 6.5536d+4 )
2679 * ..
2680 * .. Local Arrays ..
2681  INTEGER j( 2 )
2682 * ..
2683 * .. External Subroutines ..
2684  EXTERNAL pb_ladd, pb_lmul
2685 * ..
2686 * .. Intrinsic Functions ..
2687  INTRINSIC dble
2688 * ..
2689 * .. Common Blocks ..
2690  INTEGER iacs( 4 ), irand( 2 )
2691  common /rancom/ irand, iacs
2692 * ..
2693 * .. Save Statements ..
2694  SAVE /rancom/
2695 * ..
2696 * .. Executable Statements ..
2697 *
2698  pb_dran = ( dble( irand( 1 ) ) + pow16 * dble( irand( 2 ) ) ) /
2699  $ divfac
2700 *
2701  CALL pb_lmul( irand, iacs, j )
2702  CALL pb_ladd( j, iacs( 3 ), irand )
2703 *
2704  RETURN
2705 *
2706 * End of PB_DRAN
2707 *
2708  END
pb_ladd
subroutine pb_ladd(J, K, I)
Definition: pblastst.f:4480
max
#define max(A, B)
Definition: pcgemr.c:180
pb_setlocran
subroutine pb_setlocran(SEED, ILOCBLK, JLOCBLK, ILOCOFF, JLOCOFF, MYRDIST, MYCDIST, NPROW, NPCOL, JMP, IMULADD, IRAN)
Definition: pblastst.f:4302
pb_setran
subroutine pb_setran(IRAN, IAC)
Definition: pblastst.f:4759
pb_zlascal
subroutine pb_zlascal(UPLO, M, N, IOFFD, ALPHA, A, LDA)
Definition: pzblastst.f:10246
pb_lmul
subroutine pb_lmul(K, J, I)
Definition: pblastst.f:4559
pb_desctrans
subroutine pb_desctrans(DESCIN, DESCOUT)
Definition: pblastst.f:2964
pzlascal
subroutine pzlascal(TYPE, M, N, ALPHA, A, IA, JA, DESCA)
Definition: pzblastst.f:7984
pb_dran
double precision function pb_dran(IDUMM)
Definition: pdblastst.f:10440
pb_jumpit
subroutine pb_jumpit(MULADD, IRANN, IRANM)
Definition: pblastst.f:4822
pb_infog2l
subroutine pb_infog2l(I, J, DESC, NPROW, NPCOL, MYROW, MYCOL, II, JJ, PROW, PCOL)
Definition: pblastst.f:1673
pb_drand
double precision function pb_drand(IDUMM)
Definition: pdblastst.f:10378
pb_initmuladd
subroutine pb_initmuladd(MULADD0, JMP, IMULADD)
Definition: pblastst.f:4196
pb_ainfog2l
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
pb_zlagen
subroutine pb_zlagen(UPLO, AFORM, A, LDA, LCMT00, IRAN, MBLKS, IMBLOC, MB, LMBLOC, NBLKS, INBLOC, NB, LNBLOC, JMP, IMULADD)
Definition: pzblastst.f:10427
pb_initjmp
subroutine pb_initjmp(COLMAJ, NVIR, IMBVIR, INBVIR, IMBLOC, INBLOC, MB, NB, RSRC, CSRC, NPROW, NPCOL, STRIDE, JMP)
Definition: pblastst.f:4045
pzladom
subroutine pzladom(INPLACE, N, ALPHA, A, IA, JA, DESCA)
Definition: pzblastst.f:8896
pzlagen
subroutine pzlagen(INPLACE, AFORM, DIAG, OFFA, M, N, IA, JA, DESCA, IASEED, A, LDA)
Definition: pzblastst.f:8492
pb_jump
subroutine pb_jump(K, MULADD, IRANN, IRANM, IMA)
Definition: pblastst.f:4648
pxerbla
subroutine pxerbla(ICTXT, SRNAME, INFO)
Definition: pxerbla.f:2
pb_binfo
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
pb_chkmat
subroutine pb_chkmat(ICTXT, M, MPOS0, N, NPOS0, IA, JA, DESCA, DPOS0, INFO)
Definition: pblastst.f:2742
pb_locinfo
subroutine pb_locinfo(I, INB, NB, MYROC, SRCPROC, NPROCS, ILOCBLK, ILOCOFF, MYDIST)
Definition: pblastst.f:3910
min
#define min(A, B)
Definition: pcgemr.c:181