ScaLAPACK 2.1  2.1 ScaLAPACK: Scalable Linear Algebra PACKage
psmatgen2.f
Go to the documentation of this file.
1  SUBROUTINE psmatgen2( ICTXT, AFORM, DIAG, M, N, MB, NB, A, LDA,
2  \$ IAROW, IACOL, ISEED, IROFF, IRNUM, ICOFF,
3  \$ ICNUM, MYROW, MYCOL, NPROW, NPCOL )
4 *
5 *
6 * Modified version by K. L. Dackland (U added)
7 * Modified version by Peter Poromaa, Heavy DIAG
8 * Modified version by Robert Granat, R(andom) added
9 *
10 * .. Scalar Arguments ..
11  CHARACTER*1 AFORM, DIAG
12  INTEGER IACOL, IAROW, ICNUM, ICOFF, ICTXT, IRNUM,
13  \$ iroff, iseed, lda, m, mb, mycol, myrow, n,
14  \$ nb, npcol, nprow
15 * ..
16 * .. Array Arguments ..
17  REAL A( LDA, * )
18 * ..
19 *
20 * Purpose
21 * =======
22 *
23 * PSMATGEN2 : Parallel Real Double precision MATrix GENerator.
24 * Generate (or regenerate) a distributed matrix A (or sub-matrix of A).
25 *
26 * Arguments
27 * =========
28 *
29 * ICTXT (global input) INTEGER
30 * The BLACS context handle, indicating the global context of
31 * the operation. The context itself is global.
32 *
33 * AFORM (global input) CHARACTER*1
34 * if AFORM = 'U' : A returned is an Upper triangular matrix.
35 * if AFORM = 'S' : A is returned is a symmetric matrix.
36 * if AFORM = 'H' : A is returned is a Hermitian matrix.
37 * if AFORM = 'T' : A is overwritten with the transpose of
38 * what would normally be generated.
39 * if AFORM = 'C' : A is overwritten with the conjugate trans-
40 * pose of what would normally be generated.
41 * if AFORM = 'R' A random matrix is generated.
42 *
43 * DIAG (global input) CHARACTER*1
44 * if DIAG = 'D' : A is diagonally dominant.
45 *
46 * M (global input) INTEGER
47 * The number of rows in the generated distributed matrix.
48 *
49 * N (global input) INTEGER
50 * The number of columns in the generated distributed
51 * matrix.
52 *
53 * MB (global input) INTEGER
54 * The row blocking factor of the distributed matrix A.
55 *
56 * NB (global input) INTEGER
57 * The column blocking factor of the distributed matrix A.
58 *
59 * A (local output) REAL, pointer into the local
60 * memory to an array of dimension ( LDA, * ) containing the
61 * local pieces of the distributed matrix.
62 *
63 * LDA (local input) INTEGER
64 * The leading dimension of the array containing the local
65 * pieces of the distributed matrix A.
66 *
67 * IAROW (global input) INTEGER
68 * The row processor coordinate which holds the first block
69 * of the distributed matrix A.
70 *
71 * IACOL (global input) INTEGER
72 * The column processor coordinate which holds the first
73 * block of the distributed matrix A.
74 *
75 * ISEED (global input) INTEGER
76 * The seed number to generate the distributed matrix A.
77 *
78 * IROFF (local input) INTEGER
79 * The number of local rows of A that have already been
80 * generated. It should be a multiple of MB.
81 *
82 * IRNUM (local input) INTEGER
83 * The number of local rows to be generated.
84 *
85 * ICOFF (local input) INTEGER
86 * The number of local columns of A that have already been
87 * generated. It should be a multiple of NB.
88 *
89 * ICNUM (local input) INTEGER
90 * The number of local columns to be generated.
91 *
92 * MYROW (local input) INTEGER
93 * The row process coordinate of the calling process.
94 *
95 * MYCOL (local input) INTEGER
96 * The column process coordinate of the calling process.
97 *
98 * NPROW (global input) INTEGER
99 * The number of process rows in the grid.
100 *
101 * NPCOL (global input) INTEGER
102 * The number of process columns in the grid.
103 *
104 * Notes
105 * =====
106 *
107 * The code is originally developed by David Walker, ORNL,
108 * and modified by Jaeyoung Choi, ORNL.
109 *
110 * Reference: G. Fox et al.
111 * Section 12.3 of "Solving problems on concurrent processors Vol. I"
112 *
113 * =====================================================================
114 *
115 * .. Parameters ..
117  PARAMETER ( MULT0=20077, mult1=16838, iadd0=12345,
119  REAL ONE, TWO, ZERO
120  PARAMETER ( ONE = 1.0, two = 2.0, zero = 0.0)
121 * ..
122 * .. Local Scalars ..
123  LOGICAL SYMM, HERM, TRAN, UPPR, RANDOM
124  INTEGER I, IC, IK, INFO, IOFFC, IOFFR, IR, J, JK,
125  \$ jump1, jump2, jump3, jump4, jump5, jump6,
126  \$ jump7, maxmn, mend, moff, mp, mrcol, mrrow,
127  \$ nend, noff, npmb, nq, nqnb
128 * ..
129 * .. Local Arrays ..
130  INTEGER IADD(2), IA1(2), IA2(2), IA3(2), IA4(2),
131  \$ IA5(2), IB1(2), IB2(2), IB3(2), IC1(2), IC2(2),
132  \$ ic3(2), ic4(2), ic5(2), iran1(2), iran2(2),
133  \$ iran3(2), iran4(2), itmp1(2), itmp2(2),
134  \$ itmp3(2), jseed(2), mult(2)
135 * ..
136 * .. External Subroutines ..
137  EXTERNAL jumpit, pxerbla, setran, xjumpm
138 * ..
139 * .. Intrinsic Functions ..
140  INTRINSIC abs, max, mod
141 * ..
142 * .. External Functions ..
143  LOGICAL LSAME
144  INTEGER ICEIL, NUMROC
145  REAL PSRAND
146  EXTERNAL iceil, numroc, lsame, psrand
147 * ..
148 * .. Executable Statements ..
149 *
150 * Test the input arguments
151 *
152  mp = numroc( m, mb, myrow, iarow, nprow )
153  nq = numroc( n, nb, mycol, iacol, npcol )
154  symm = lsame( aform, 'S' )
155  uppr = lsame( aform, 'U' )
156  herm = lsame( aform, 'H' )
157  tran = lsame( aform, 'T' )
158  random = lsame( aform, 'R' )
159 *
160  info = 0
161  IF( .NOT.( uppr.OR.symm.OR.herm.OR.tran.OR.random ) .AND.
162  \$ .NOT.lsame( aform, 'C' ) .AND.
163  \$ .NOT.lsame( aform, 'N' ) ) THEN
164  info = 2
165  ELSE IF( .NOT.lsame( diag, 'D' ) .AND.
166  \$ .NOT.lsame( diag, 'N' ) ) THEN
167  info = 3
168  ELSE IF( uppr.OR.symm.OR.herm ) THEN
169  IF( m.NE.n ) THEN
170  info = 5
171  ELSE IF( mb.NE.nb ) THEN
172  info = 7
173  END IF
174  ELSE IF( m.LT.0 ) THEN
175  info = 4
176  ELSE IF( n.LT.0 ) THEN
177  info = 5
178  ELSE IF( mb.LT.1 ) THEN
179  info = 6
180  ELSE IF( nb.LT.1 ) THEN
181  info = 7
182  ELSE IF( lda.LT.0 ) THEN
183  info = 9
184  ELSE IF( ( iarow.LT.0 ).OR.( iarow.GE.nprow ) ) THEN
185  info = 10
186  ELSE IF( ( iacol.LT.0 ).OR.( iacol.GE.npcol ) ) THEN
187  info = 11
188  ELSE IF( mod(iroff,mb).GT.0 ) THEN
189  info = 13
190  ELSE IF( irnum.GT.(mp-iroff) ) THEN
191  info = 14
192  ELSE IF( mod(icoff,nb).GT.0 ) THEN
193  info = 15
194  ELSE IF( icnum.GT.(nq-icoff) ) THEN
195  info = 16
196  ELSE IF( ( myrow.LT.0 ).OR.( myrow.GE.nprow ) ) THEN
197  info = 17
198  ELSE IF( ( mycol.LT.0 ).OR.( mycol.GE.npcol ) ) THEN
199  info = 18
200  END IF
201  IF( info.NE.0 ) THEN
202  CALL pxerbla( ictxt, 'PSMATGEN2', info )
203  RETURN
204  END IF
205  mrrow = mod( nprow+myrow-iarow, nprow )
206  mrcol = mod( npcol+mycol-iacol, npcol )
207  npmb = nprow * mb
208  nqnb = npcol * nb
209  moff = iroff / mb
210  noff = icoff / nb
211  mend = iceil(irnum, mb) + moff
212  nend = iceil(icnum, nb) + noff
213 *
214  mult(1) = mult0
215  mult(2) = mult1
218  jseed(1) = iseed
219  jseed(2) = 0
220 *
221 * Symmetric or Hermitian matrix will be generated.
222 *
223  IF( symm.OR.herm ) THEN
224 *
225 * First, generate the lower triangular part (with diagonal block)
226 *
227  jump1 = 1
228  jump2 = npmb
229  jump3 = m
230  jump4 = nqnb
231  jump5 = nb
232  jump6 = mrcol
233  jump7 = mb*mrrow
234 *
235  CALL xjumpm( jump1, mult, iadd, jseed, iran1, ia1, ic1 )
236  CALL xjumpm( jump2, mult, iadd, iran1, itmp1, ia2, ic2 )
237  CALL xjumpm( jump3, mult, iadd, iran1, itmp1, ia3, ic3 )
238  CALL xjumpm( jump4, ia3, ic3, iran1, itmp1, ia4, ic4 )
239  CALL xjumpm( jump5, ia3, ic3, iran1, itmp1, ia5, ic5 )
240  CALL xjumpm( jump6, ia5, ic5, iran1, itmp3, itmp1, itmp2 )
241  CALL xjumpm( jump7, mult, iadd, itmp3, iran1, itmp1, itmp2 )
242  CALL xjumpm( noff, ia4, ic4, iran1, itmp1, itmp2, itmp3 )
243  CALL xjumpm( moff, ia2, ic2, itmp1, iran1, itmp2, itmp3 )
244  CALL setran( iran1, ia1, ic1 )
245 *
246  DO 10 i = 1, 2
247  ib1(i) = iran1(i)
248  ib2(i) = iran1(i)
249  ib3(i) = iran1(i)
250  10 CONTINUE
251 *
252  jk = 1
253  DO 80 ic = noff+1, nend
254  ioffc = ((ic-1)*npcol+mrcol) * nb
255  DO 70 i = 1, nb
256  IF( jk .GT. icnum ) GO TO 90
257 *
258  ik = 1
259  DO 50 ir = moff+1, mend
260  ioffr = ((ir-1)*nprow+mrrow) * mb
261 *
262  IF( ioffr .GT. ioffc ) THEN
263  DO 20 j = 1, mb
264  IF( ik .GT. irnum ) GO TO 60
265  a(ik,jk) = one - two*psrand(0)
266  ik = ik + 1
267  20 CONTINUE
268 *
269  ELSE IF( ioffc .EQ. ioffr ) THEN
270  ik = ik + i - 1
271  IF( ik .GT. irnum ) GO TO 60
272  DO 30 j = 1, i-1
273  a(ik,jk) = one - two*psrand(0)
274  30 CONTINUE
275  a(ik,jk) = one - two*psrand(0)
276  DO 40 j = 1, mb-i
277  IF( ik+j .GT. irnum ) GO TO 60
278  a(ik+j,jk) = one - two*psrand(0)
279  a(ik,jk+j) = a(ik+j,jk)
280  40 CONTINUE
281  ik = ik + mb - i + 1
282  ELSE
283  ik = ik + mb
284  END IF
285 *
286  CALL jumpit( ia2, ic2, ib1, iran2 )
287  ib1(1) = iran2(1)
288  ib1(2) = iran2(2)
289  50 CONTINUE
290 *
291  60 CONTINUE
292  jk = jk + 1
293  CALL jumpit( ia3, ic3, ib2, iran3 )
294  ib1(1) = iran3(1)
295  ib1(2) = iran3(2)
296  ib2(1) = iran3(1)
297  ib2(2) = iran3(2)
298  70 CONTINUE
299 *
300  CALL jumpit( ia4, ic4, ib3, iran4 )
301  ib1(1) = iran4(1)
302  ib1(2) = iran4(2)
303  ib2(1) = iran4(1)
304  ib2(2) = iran4(2)
305  ib3(1) = iran4(1)
306  ib3(2) = iran4(2)
307  80 CONTINUE
308 *
309 * Next, generate the upper triangular part.
310 *
311  90 CONTINUE
312  mult(1) = mult0
313  mult(2) = mult1
316  jseed(1) = iseed
317  jseed(2) = 0
318 *
319  jump1 = 1
320  jump2 = nqnb
321  jump3 = n
322  jump4 = npmb
323  jump5 = mb
324  jump6 = mrrow
325  jump7 = nb*mrcol
326 *
327  CALL xjumpm( jump1, mult, iadd, jseed, iran1, ia1, ic1 )
328  CALL xjumpm( jump2, mult, iadd, iran1, itmp1, ia2, ic2 )
329  CALL xjumpm( jump3, mult, iadd, iran1, itmp1, ia3, ic3 )
330  CALL xjumpm( jump4, ia3, ic3, iran1, itmp1, ia4, ic4 )
331  CALL xjumpm( jump5, ia3, ic3, iran1, itmp1, ia5, ic5 )
332  CALL xjumpm( jump6, ia5, ic5, iran1, itmp3, itmp1, itmp2 )
333  CALL xjumpm( jump7, mult, iadd, itmp3, iran1, itmp1, itmp2 )
334  CALL xjumpm( moff, ia4, ic4, iran1, itmp1, itmp2, itmp3 )
335  CALL xjumpm( noff, ia2, ic2, itmp1, iran1, itmp2, itmp3 )
336  CALL setran( iran1, ia1, ic1 )
337 *
338  DO 100 i = 1, 2
339  ib1(i) = iran1(i)
340  ib2(i) = iran1(i)
341  ib3(i) = iran1(i)
342  100 CONTINUE
343 *
344  ik = 1
345  DO 150 ir = moff+1, mend
346  ioffr = ((ir-1)*nprow+mrrow) * mb
347  DO 140 j = 1, mb
348  IF( ik .GT. irnum ) GO TO 160
349  jk = 1
350  DO 120 ic = noff+1, nend
351  ioffc = ((ic-1)*npcol+mrcol) * nb
352  IF( ioffc .GT. ioffr ) THEN
353  DO 110 i = 1, nb
354  IF( jk .GT. icnum ) GO TO 130
355  a(ik,jk) = one - two*psrand(0)
356  jk = jk + 1
357  110 CONTINUE
358  ELSE
359  jk = jk + nb
360  END IF
361  CALL jumpit( ia2, ic2, ib1, iran2 )
362  ib1(1) = iran2(1)
363  ib1(2) = iran2(2)
364  120 CONTINUE
365 *
366  130 CONTINUE
367  ik = ik + 1
368  CALL jumpit( ia3, ic3, ib2, iran3 )
369  ib1(1) = iran3(1)
370  ib1(2) = iran3(2)
371  ib2(1) = iran3(1)
372  ib2(2) = iran3(2)
373  140 CONTINUE
374 *
375  CALL jumpit( ia4, ic4, ib3, iran4 )
376  ib1(1) = iran4(1)
377  ib1(2) = iran4(2)
378  ib2(1) = iran4(1)
379  ib2(2) = iran4(2)
380  ib3(1) = iran4(1)
381  ib3(2) = iran4(2)
382  150 CONTINUE
383  160 CONTINUE
384 *
385 * Generate an upper triangular matrix.
386 *
387  ELSE IF ( uppr ) THEN
388  jump1 = 1
389  jump2 = npmb
390  jump3 = m
391  jump4 = nqnb
392  jump5 = nb
393  jump6 = mrcol
394  jump7 = mb*mrrow
395 *
396  CALL xjumpm( jump1, mult, iadd, jseed, iran1, ia1, ic1 )
397  CALL xjumpm( jump2, mult, iadd, iran1, itmp1, ia2, ic2 )
398  CALL xjumpm( jump3, mult, iadd, iran1, itmp1, ia3, ic3 )
399  CALL xjumpm( jump4, ia3, ic3, iran1, itmp1, ia4, ic4 )
400  CALL xjumpm( jump5, ia3, ic3, iran1, itmp1, ia5, ic5 )
401  CALL xjumpm( jump6, ia5, ic5, iran1, itmp3, itmp1, itmp2 )
402  CALL xjumpm( jump7, mult, iadd, itmp3, iran1, itmp1, itmp2 )
403  CALL xjumpm( noff, ia4, ic4, iran1, itmp1, itmp2, itmp3 )
404  CALL xjumpm( moff, ia2, ic2, itmp1, iran1, itmp2, itmp3 )
405  CALL setran( iran1, ia1, ic1 )
406 *
407  DO 1000 i = 1, 2
408  ib1(i) = iran1(i)
409  ib2(i) = iran1(i)
410  ib3(i) = iran1(i)
411  1000 CONTINUE
412 *
413  jk = 1
414  DO 8000 ic = noff+1, nend
415  ioffc = ((ic-1)*npcol+mrcol) * nb
416  DO 7000 i = 1, nb
417  IF( jk .GT. icnum ) GO TO 8000
418 *
419  ik = 1
420  DO 5000 ir = moff+1, mend
421  ioffr = ((ir-1)*nprow+mrrow) * mb
422 *
423  IF( ioffc .EQ. ioffr ) THEN
424  ik = ik + i - 1
425  IF( ik .GT. irnum ) GO TO 6000
426  DO 3000 j = 1, i-1
427  a(ik,jk) = one - two*psrand(0)
428  3000 CONTINUE
429  a(ik,jk) = one - two*psrand(0)
430  DO 4000 j = 1, mb-i
431  IF( ik+j .GT. irnum ) GO TO 6000
432  a(ik,jk+j) = one - two*psrand(0)
433  4000 CONTINUE
434  ik = ik + mb - i + 1
435  ELSE
436  ik = ik + mb
437  END IF
438 *
439  CALL jumpit( ia2, ic2, ib1, iran2 )
440  ib1(1) = iran2(1)
441  ib1(2) = iran2(2)
442  5000 CONTINUE
443 *
444  6000 CONTINUE
445  jk = jk + 1
446  CALL jumpit( ia3, ic3, ib2, iran3 )
447  ib1(1) = iran3(1)
448  ib1(2) = iran3(2)
449  ib2(1) = iran3(1)
450  ib2(2) = iran3(2)
451  7000 CONTINUE
452 *
453  CALL jumpit( ia4, ic4, ib3, iran4 )
454  ib1(1) = iran4(1)
455  ib1(2) = iran4(2)
456  ib2(1) = iran4(1)
457  ib2(2) = iran4(2)
458  ib3(1) = iran4(1)
459  ib3(2) = iran4(2)
460  8000 CONTINUE
461  mult(1) = mult0
462  mult(2) = mult1
465  jseed(1) = iseed
466  jseed(2) = 0
467 *
468  jump1 = 1
469  jump2 = nqnb
470  jump3 = n
471  jump4 = npmb
472  jump5 = mb
473  jump6 = mrrow
474  jump7 = nb*mrcol
475 *
476  CALL xjumpm( jump1, mult, iadd, jseed, iran1, ia1, ic1 )
477  CALL xjumpm( jump2, mult, iadd, iran1, itmp1, ia2, ic2 )
478  CALL xjumpm( jump3, mult, iadd, iran1, itmp1, ia3, ic3 )
479  CALL xjumpm( jump4, ia3, ic3, iran1, itmp1, ia4, ic4 )
480  CALL xjumpm( jump5, ia3, ic3, iran1, itmp1, ia5, ic5 )
481  CALL xjumpm( jump6, ia5, ic5, iran1, itmp3, itmp1, itmp2 )
482  CALL xjumpm( jump7, mult, iadd, itmp3, iran1, itmp1, itmp2 )
483  CALL xjumpm( moff, ia4, ic4, iran1, itmp1, itmp2, itmp3 )
484  CALL xjumpm( noff, ia2, ic2, itmp1, iran1, itmp2, itmp3 )
485  CALL setran( iran1, ia1, ic1 )
486 *
487  DO 1110 i = 1, 2
488  ib1(i) = iran1(i)
489  ib2(i) = iran1(i)
490  ib3(i) = iran1(i)
491  1110 CONTINUE
492 *
493  ik = 1
494  DO 1500 ir = moff+1, mend
495  ioffr = ((ir-1)*nprow+mrrow) * mb
496  DO 1400 j = 1, mb
497  IF( ik .GT. irnum ) GO TO 1600
498  jk = 1
499  DO 1200 ic = noff+1, nend
500  ioffc = ((ic-1)*npcol+mrcol) * nb
501  IF( ioffc .GT. ioffr ) THEN
502  DO 1100 i = 1, nb
503  IF( jk .GT. icnum ) GO TO 1300
504  a(ik,jk) = one - two*psrand(0)
505  jk = jk + 1
506  1100 CONTINUE
507  ELSE
508  jk = jk + nb
509  END IF
510  CALL jumpit( ia2, ic2, ib1, iran2 )
511  ib1(1) = iran2(1)
512  ib1(2) = iran2(2)
513  1200 CONTINUE
514 *
515  1300 CONTINUE
516  ik = ik + 1
517  CALL jumpit( ia3, ic3, ib2, iran3 )
518  ib1(1) = iran3(1)
519  ib1(2) = iran3(2)
520  ib2(1) = iran3(1)
521  ib2(2) = iran3(2)
522  1400 CONTINUE
523 *
524  CALL jumpit( ia4, ic4, ib3, iran4 )
525  ib1(1) = iran4(1)
526  ib1(2) = iran4(2)
527  ib2(1) = iran4(1)
528  ib2(2) = iran4(2)
529  ib3(1) = iran4(1)
530  ib3(2) = iran4(2)
531  1500 CONTINUE
532  1600 CONTINUE
533 *
534 * (Conjugate) Transposed matrix A will be generated.
535 *
536  ELSE IF( tran .OR. lsame( aform, 'C' ) ) THEN
537 *
538  jump1 = 1
539  jump2 = nqnb
540  jump3 = n
541  jump4 = npmb
542  jump5 = mb
543  jump6 = mrrow
544  jump7 = nb*mrcol
545 *
546  CALL xjumpm( jump1, mult, iadd, jseed, iran1, ia1, ic1 )
547  CALL xjumpm( jump2, mult, iadd, iran1, itmp1, ia2, ic2 )
548  CALL xjumpm( jump3, mult, iadd, iran1, itmp1, ia3, ic3 )
549  CALL xjumpm( jump4, ia3, ic3, iran1, itmp1, ia4, ic4 )
550  CALL xjumpm( jump5, ia3, ic3, iran1, itmp1, ia5, ic5 )
551  CALL xjumpm( jump6, ia5, ic5, iran1, itmp3, itmp1, itmp2 )
552  CALL xjumpm( jump7, mult, iadd, itmp3, iran1, itmp1, itmp2 )
553  CALL xjumpm( moff, ia4, ic4, iran1, itmp1, itmp2, itmp3 )
554  CALL xjumpm( noff, ia2, ic2, itmp1, iran1, itmp2, itmp3 )
555  CALL setran( iran1, ia1, ic1 )
556 *
557  DO 170 i = 1, 2
558  ib1(i) = iran1(i)
559  ib2(i) = iran1(i)
560  ib3(i) = iran1(i)
561  170 CONTINUE
562 *
563  ik = 1
564  DO 220 ir = moff+1, mend
565  ioffr = ((ir-1)*nprow+mrrow) * mb
566  DO 210 j = 1, mb
567  IF( ik .GT. irnum ) GO TO 230
568  jk = 1
569  DO 190 ic = noff+1, nend
570  ioffc = ((ic-1)*npcol+mrcol) * nb
571  DO 180 i = 1, nb
572  IF( jk .GT. icnum ) GO TO 200
573  a(ik,jk) = one - two*psrand(0)
574  jk = jk + 1
575  180 CONTINUE
576  CALL jumpit( ia2, ic2, ib1, iran2 )
577  ib1(1) = iran2(1)
578  ib1(2) = iran2(2)
579  190 CONTINUE
580 *
581  200 CONTINUE
582  ik = ik + 1
583  CALL jumpit( ia3, ic3, ib2, iran3 )
584  ib1(1) = iran3(1)
585  ib1(2) = iran3(2)
586  ib2(1) = iran3(1)
587  ib2(2) = iran3(2)
588  210 CONTINUE
589 *
590  CALL jumpit( ia4, ic4, ib3, iran4 )
591  ib1(1) = iran4(1)
592  ib1(2) = iran4(2)
593  ib2(1) = iran4(1)
594  ib2(2) = iran4(2)
595  ib3(1) = iran4(1)
596  ib3(2) = iran4(2)
597  220 CONTINUE
598  230 CONTINUE
599 *
600 * A random matrix is generated.
601 *
602  ELSEIF( random ) THEN
603 *
604  jump1 = 1
605  jump2 = npmb
606  jump3 = m
607  jump4 = nqnb
608  jump5 = nb
609  jump6 = mrcol
610  jump7 = mb*mrrow
611 *
612  CALL xjumpm( jump1, mult, iadd, jseed, iran1, ia1, ic1 )
613  CALL xjumpm( jump2, mult, iadd, iran1, itmp1, ia2, ic2 )
614  CALL xjumpm( jump3, mult, iadd, iran1, itmp1, ia3, ic3 )
615  CALL xjumpm( jump4, ia3, ic3, iran1, itmp1, ia4, ic4 )
616  CALL xjumpm( jump5, ia3, ic3, iran1, itmp1, ia5, ic5 )
617  CALL xjumpm( jump6, ia5, ic5, iran1, itmp3, itmp1, itmp2 )
618  CALL xjumpm( jump7, mult, iadd, itmp3, iran1, itmp1, itmp2 )
619  CALL xjumpm( noff, ia4, ic4, iran1, itmp1, itmp2, itmp3 )
620  CALL xjumpm( moff, ia2, ic2, itmp1, iran1, itmp2, itmp3 )
621  CALL setran( iran1, ia1, ic1 )
622 *
623  DO 240 i = 1, 2
624  ib1(i) = iran1(i)
625  ib2(i) = iran1(i)
626  ib3(i) = iran1(i)
627  240 CONTINUE
628 *
629  jk = 1
630  DO 290 ic = noff+1, nend
631  ioffc = ((ic-1)*npcol+mrcol) * nb
632  DO 280 i = 1, nb
633  IF( jk .GT. icnum ) GO TO 300
634  ik = 1
635  DO 260 ir = moff+1, mend
636  ioffr = ((ir-1)*nprow+mrrow) * mb
637  DO 250 j = 1, mb
638  IF( ik .GT. irnum ) GO TO 270
639  a(ik,jk) = one - two*psrand(0)
640  ik = ik + 1
641  250 CONTINUE
642  CALL jumpit( ia2, ic2, ib1, iran2 )
643  ib1(1) = iran2(1)
644  ib1(2) = iran2(2)
645  260 CONTINUE
646 *
647  270 CONTINUE
648  jk = jk + 1
649  CALL jumpit( ia3, ic3, ib2, iran3 )
650  ib1(1) = iran3(1)
651  ib1(2) = iran3(2)
652  ib2(1) = iran3(1)
653  ib2(2) = iran3(2)
654  280 CONTINUE
655 *
656  CALL jumpit( ia4, ic4, ib3, iran4 )
657  ib1(1) = iran4(1)
658  ib1(2) = iran4(2)
659  ib2(1) = iran4(1)
660  ib2(2) = iran4(2)
661  ib3(1) = iran4(1)
662  ib3(2) = iran4(2)
663  290 CONTINUE
664  300 CONTINUE
665  END IF
666 *
667 * Diagonally dominant matrix will be generated.
668 *
669  IF( lsame( diag, 'D' ) ) THEN
670  IF( mb.NE.nb ) THEN
671  WRITE(*,*) 'Diagonally dominant matrices with rowNB not'//
672  \$ ' equal colNB is not supported!'
673  RETURN
674  END IF
675 *
676  maxmn = max(m, n)
677  jk = 1
678  DO 340 ic = noff+1, nend
679  ioffc = ((ic-1)*npcol+mrcol) * nb
680  ik = 1
681  DO 320 ir = moff+1, mend
682  ioffr = ((ir-1)*nprow+mrrow) * mb
683  IF( ioffc.EQ.ioffr ) THEN
684  DO 310 j = 0, mb-1
685  IF( ik .GT. irnum ) GO TO 330
686  a(ik,jk+j) = abs(a(ik,jk+j)) + maxmn
687  ik = ik + 1
688  310 CONTINUE
689  ELSE
690  ik = ik + mb
691  END IF
692  320 CONTINUE
693  330 CONTINUE
694  jk = jk + nb
695  340 CONTINUE
696  END IF
697 *
698  RETURN
699 *
700 * End of PSMATGEN2
701 *
702  END
max
#define max(A, B)
Definition: pcgemr.c:180
jumpit