ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pzlatms.f
Go to the documentation of this file.
1 *
2 *
3  SUBROUTINE pzlatms( M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX,
4  $ KL, KU, PACK, A, IA, JA, DESCA, ORDER, WORK,
5  $ LWORK, INFO )
6 *
7 * -- ScaLAPACK routine (version 1.7) --
8 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
9 * and University of California, Berkeley.
10 * May 1, 1997
11 *
12 * .. Scalar Arguments ..
13  CHARACTER DIST, PACK, SYM
14  INTEGER IA, INFO, JA, KL, KU, LWORK, M, MODE, N, ORDER
15  DOUBLE PRECISION COND, DMAX
16 * ..
17 * .. Array Arguments ..
18  INTEGER DESCA( * ), ISEED( 4 )
19  DOUBLE PRECISION D( * )
20  COMPLEX*16 A( * ), WORK( * )
21 * ..
22 *
23 * Purpose
24 * =======
25 *
26 * PZLATMS generates random Hermitian matrices with specified
27 * eigenvalues for testing SCALAPACK programs.
28 *
29 * PZLATMS operates by applying the following sequence of
30 * operations:
31 *
32 * Set the diagonal to D, where D may be input or
33 * computed according to MODE, COND, DMAX, and SYM
34 * as described below.
35 *
36 * Generate a dense M x N matrix by multiplying D on the left
37 * and the right by random unitary matrices, then:
38 *
39 * Reduce the bandwidth according to KL and KU, using
40 * Householder transformations.
41 * ### bandwidth reduction NOT SUPPORTED ###
42 *
43 * Arguments
44 * =========
45 *
46 * M - (global input) INTEGER
47 * The number of rows of A. Not modified.
48 *
49 * N - (global input) INTEGER
50 * The number of columns of A. Not modified.
51 * ### M .ne. N unsupported
52 *
53 * DIST - (global input) CHARACTER*1
54 * On entry, DIST specifies the type of distribution to be used
55 * to generate the random eigen-/singular values.
56 * 'U' => UNIFORM( 0, 1 ) ( 'U' for uniform )
57 * 'S' => UNIFORM( -1, 1 ) ( 'S' for symmetric )
58 * 'N' => NORMAL( 0, 1 ) ( 'N' for normal )
59 * Not modified.
60 *
61 * ISEED - (global input) INTEGER array, dimension ( 4 )
62 * On entry ISEED specifies the seed of the random number
63 * generator. They should lie between 0 and 4095 inclusive,
64 * and ISEED(4) should be odd. The random number generator
65 * uses a linear congruential sequence limited to small
66 * integers, and so should produce machine independent
67 * random numbers. The values of ISEED are changed on
68 * exit, and can be used in the next call to ZLATMS
69 * to continue the same random number sequence.
70 * Changed on exit.
71 *
72 * SYM - (global input) CHARACTER*1
73 * If SYM='S' or 'H', the generated matrix is Hermitian, with
74 * eigenvalues specified by D, COND, MODE, and DMAX; they
75 * may be positive, negative, or zero.
76 * If SYM='P', the generated matrix is Hermitian, with
77 * eigenvalues (= singular values) specified by D, COND,
78 * MODE, and DMAX; they will not be negative.
79 * If SYM='N', the generated matrix is nonsymmetric, with
80 * singular values specified by D, COND, MODE, and DMAX;
81 * they will not be negative.
82 * ### SYM = 'N' NOT SUPPORTED ###
83 * Not modified.
84 *
85 * D - (local input/output) DOUBLE PRECISION array,
86 * dimension ( MIN( M , N ) )
87 * This array is used to specify the singular values or
88 * eigenvalues of A (see SYM, above.) If MODE=0, then D is
89 * assumed to contain the singular/eigenvalues, otherwise
90 * they will be computed according to MODE, COND, and DMAX,
91 * and placed in D.
92 * Modified if MODE is nonzero.
93 *
94 * MODE - (global input) INTEGER
95 * On entry this describes how the singular/eigenvalues are to
96 * be specified:
97 * MODE = 0 means use D as input
98 * MODE = 1 sets D(1)=1 and D(2:N)=1.0/COND
99 * MODE = 2 sets D(1:N-1)=1 and D(N)=1.0/COND
100 * MODE = 3 sets D(I)=COND**(-(I-1)/(N-1))
101 * MODE = 4 sets D(i)=1 - (i-1)/(N-1)*(1 - 1/COND)
102 * MODE = 5 sets D to random numbers in the range
103 * ( 1/COND , 1 ) such that their logarithms
104 * are uniformly distributed.
105 * MODE = 6 set D to random numbers from same distribution
106 * as the rest of the matrix.
107 * MODE < 0 has the same meaning as ABS(MODE), except that
108 * the order of the elements of D is reversed.
109 * Thus if MODE is positive, D has entries ranging from
110 * 1 to 1/COND, if negative, from 1/COND to 1,
111 * If SYM='S' or 'H', and MODE is neither 0, 6, nor -6, then
112 * the elements of D will also be multiplied by a random
113 * sign (i.e., +1 or -1.)
114 * Not modified.
115 *
116 * COND - (global input) DOUBLE PRECISION
117 * On entry, this is used as described under MODE above.
118 * If used, it must be >= 1. Not modified.
119 *
120 * DMAX - (global input) DOUBLE PRECISION
121 * If MODE is neither -6, 0 nor 6, the contents of D, as
122 * computed according to MODE and COND, will be scaled by
123 * DMAX / max(abs(D(i))); thus, the maximum absolute eigen- or
124 * singular value (which is to say the norm) will be abs(DMAX).
125 * Note that DMAX need not be positive: if DMAX is negative
126 * (or zero), D will be scaled by a negative number (or zero).
127 * Not modified.
128 *
129 * KL - (global input) INTEGER
130 * This specifies the lower bandwidth of the matrix. For
131 * example, KL=0 implies upper triangular, KL=1 implies upper
132 * Hessenberg, and KL being at least M-1 means that the matrix
133 * has full lower bandwidth. KL must equal KU if the matrix
134 * is Hermitian.
135 * Not modified.
136 * ### 1 <= KL < N-1 is NOT SUPPORTED ###
137 *
138 * KU - (global input) INTEGER
139 * This specifies the upper bandwidth of the matrix. For
140 * example, KU=0 implies lower triangular, KU=1 implies lower
141 * Hessenberg, and KU being at least N-1 means that the matrix
142 * has full upper bandwidth. KL must equal KU if the matrix
143 * is Hermitian.
144 * Not modified.
145 * ### 1 <= KU < N-1 is NOT SUPPORTED ###
146 *
147 * PACK - (global input) CHARACTER*1
148 * This specifies packing of matrix as follows:
149 * 'N' => no packing
150 * ### PACK must be 'N' all other options NOT SUPPORTED ###
151 *
152 * A - (local output) COMPLEX*16 array
153 * Global dimension (M, N), local dimension (MP, NQ)
154 * On exit A is the desired test matrix.
155 *
156 * IA (global input) INTEGER
157 * A's global row index, which points to the beginning of the
158 * submatrix which is to be operated on.
159 *
160 * JA (global input) INTEGER
161 * A's global column index, which points to the beginning of
162 * the submatrix which is to be operated on.
163 *
164 * DESCA (global and local input) INTEGER array of dimension DLEN_.
165 * The array descriptor for the distributed matrix A.
166 *
167 * ORDER - (input) INTEGER
168 * The number of reflectors used to define the orthogonal
169 * matrix Q. A = Q * D * Q'
170 * Higher ORDER requires more computation and communication.
171 *
172 * WORK - (local input/output) COMPLEX*16 array,
173 * dimension (LWORK)
174 *
175 * LWORK - (local input) INTEGER dimension of WORK
176 * LWORK >= SIZETMS as returned by PZLASIZESEP
177 *
178 * INFO - (global output) INTEGER
179 * Error code. On exit, INFO will be set to one of the
180 * following values:
181 * 0 => normal return
182 * -1 => M negative or unequal to N and SYM='S', 'H', or 'P'
183 * -2 => N negative
184 * -3 => DIST illegal string
185 * -5 => SYM illegal string
186 * -7 => MODE not in range -6 to 6
187 * -8 => COND less than 1.0, and MODE neither -6, 0 nor 6
188 * -10 => KL negative
189 * -11 => KU negative, or SYM='S' or 'H' and KU not equal to KL
190 * -16 => DESCA is inconsistent
191 * -17 => ORDER not in the range 0 to N inclusive
192 * 1 => Error return from DLATM1
193 * 2 => Cannot scale to DMAX (max. sing. value is 0)
194 * 3 => Error return from PZLAGHE
195 *
196 *-----------------------------------------------------------------------
197 *
198 *
199 * .. Parameters ..
200  INTEGER BLOCK_CYCLIC_2D, DLEN_, DTYPE_, CTXT_, M_, N_,
201  $ MB_, NB_, RSRC_, CSRC_, LLD_
202  parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
203  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
204  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
205  DOUBLE PRECISION ZERO, ONE
206  PARAMETER ( ZERO = 0.0d+0, one = 1.0d+0 )
207  COMPLEX*16 ZZERO
208  parameter( zzero = ( 0.0d+0, 0.0d+0 ) )
209 * ..
210 * .. Local Scalars ..
211  INTEGER I, IDIST, IINFO, IPACK, IRSIGN, ISYM, LLB,
212  $ MNMIN, MYCOL, MYROW, NP, NPCOL, NPROW, NQ
213  DOUBLE PRECISION ALPHA, TEMP
214 * ..
215 * .. Local Arrays ..
216  INTEGER IDUM1( 1 ), IDUM2( 1 )
217 * ..
218 * .. External Functions ..
219  LOGICAL LSAME
220  INTEGER NUMROC
221  EXTERNAL lsame, numroc
222 * ..
223 * .. External Subroutines ..
224  EXTERNAL blacs_gridinfo, chk1mat, dlatm1, dscal,
225  $ pchk1mat, pxerbla, pzlaghe, zlaset
226 * ..
227 * .. Intrinsic Functions ..
228  INTRINSIC abs, max, min, mod
229 * ..
230 * .. Executable Statements ..
231 * This is just to keep ftnchek happy
232  IF( block_cyclic_2d*csrc_*ctxt_*dlen_*dtype_*lld_*mb_*m_*nb_*n_*
233  $ rsrc_.LT.0 )RETURN
234 *
235 * 1) Decode and Test the input parameters.
236 * Initialize flags & seed.
237 *
238 *
239  info = 0
240 *
241  CALL blacs_gridinfo( desca( ctxt_ ), nprow, npcol, myrow, mycol )
242  IF( ( myrow.GE.nprow .OR. myrow.LT.0 ) .OR.
243  $ ( mycol.GE.npcol .OR. mycol.LT.0 ) )RETURN
244 *
245  np = numroc( n, desca( mb_ ), myrow, 0, nprow )
246  nq = numroc( n, desca( nb_ ), mycol, 0, npcol )
247 *
248 * Quick return if possible
249 *
250  IF( m.EQ.0 .OR. n.EQ.0 )
251  $ RETURN
252 *
253 * Decode DIST
254 *
255  IF( lsame( dist, 'U' ) ) THEN
256  idist = 1
257  ELSE IF( lsame( dist, 'S' ) ) THEN
258  idist = 2
259  ELSE IF( lsame( dist, 'N' ) ) THEN
260  idist = 3
261  ELSE
262  idist = -1
263  END IF
264 *
265 * Decode SYM
266 *
267  IF( lsame( sym, 'N' ) ) THEN
268  isym = 1
269  irsign = 0
270  ELSE IF( lsame( sym, 'P' ) ) THEN
271  isym = 2
272  irsign = 0
273  ELSE IF( lsame( sym, 'S' ) ) THEN
274  isym = 2
275  irsign = 1
276  ELSE IF( lsame( sym, 'H' ) ) THEN
277  isym = 2
278  irsign = 1
279  ELSE
280  isym = -1
281  END IF
282 *
283 * Decode PACK
284 *
285  IF( lsame( pack, 'N' ) ) THEN
286  ipack = 0
287  ELSE
288  ipack = 1
289  END IF
290 *
291 * Set certain internal parameters
292 *
293  mnmin = min( m, n )
294  llb = min( kl, m-1 )
295 *
296  IF( order.EQ.0 )
297  $ order = n
298 *
299 * Set INFO if an error
300 *
301  IF( nprow.EQ.-1 ) THEN
302  info = -( 1600+ctxt_ )
303  ELSE
304  CALL chk1mat( m, 1, n, 2, ia, ja, desca, 16, info )
305  IF( info.EQ.0 ) THEN
306  IF( m.NE.n .AND. isym.NE.1 ) THEN
307  info = -2
308  ELSE IF( idist.EQ.-1 ) THEN
309  info = -3
310  ELSE IF( isym.EQ.-1 ) THEN
311  info = -5
312  ELSE IF( abs( mode ).GT.6 ) THEN
313  info = -7
314  ELSE IF( ( mode.NE.0 .AND. abs( mode ).NE.6 ) .AND. cond.LT.
315  $ one ) THEN
316  info = -8
317  ELSE IF( kl.LT.0 ) THEN
318  info = -10
319  ELSE IF( ku.LT.0 .OR. ( isym.NE.1 .AND. kl.NE.ku ) ) THEN
320  info = -11
321  ELSE IF( ( order.LT.0 ) .OR. ( order.GT.n ) ) THEN
322  info = -17
323  END IF
324  END IF
325  CALL pchk1mat( m, 1, n, 2, ia, ja, desca, 16, 0, idum1, idum2,
326  $ info )
327  END IF
328 *
329 * Check for unsupported features
330 *
331  IF( isym.NE.2 ) THEN
332  info = -5
333  ELSE IF( ipack.NE.0 ) THEN
334  info = -12
335  ELSE IF( kl.GT.0 .AND. kl.LT.m-1 ) THEN
336  info = -10
337  ELSE IF( ku.GT.0 .AND. ku.LT.n-1 ) THEN
338  info = -11
339  ELSE IF( llb.NE.0 .AND. llb.NE.m-1 ) THEN
340  info = -10
341  END IF
342  IF( info.NE.0 ) THEN
343  CALL pxerbla( desca( ctxt_ ), 'PZLATMS', -info )
344  RETURN
345  END IF
346 *
347 * Initialize random number generator
348 *
349  DO 10 i = 1, 4
350  iseed( i ) = mod( abs( iseed( i ) ), 4096 )
351  10 CONTINUE
352 *
353  IF( mod( iseed( 4 ), 2 ).NE.1 )
354  $ iseed( 4 ) = iseed( 4 ) + 1
355 *
356 * 2) Set up D if indicated.
357 *
358 * Compute D according to COND and MODE
359 *
360  CALL dlatm1( mode, cond, irsign, idist, iseed, d, mnmin, iinfo )
361 *
362  IF( iinfo.NE.0 ) THEN
363  info = 1
364  RETURN
365  END IF
366 *
367 *
368  IF( mode.NE.0 .AND. abs( mode ).NE.6 ) THEN
369 *
370 * Scale by DMAX
371 *
372  temp = abs( d( 1 ) )
373  DO 20 i = 2, mnmin
374  temp = max( temp, abs( d( i ) ) )
375  20 CONTINUE
376 *
377  IF( temp.GT.zero ) THEN
378  alpha = dmax / temp
379  ELSE
380  info = 2
381  RETURN
382  END IF
383 *
384  CALL dscal( mnmin, alpha, d, 1 )
385 *
386  END IF
387 *
388  CALL zlaset( 'A', np, nq, zzero, zzero, a, desca( lld_ ) )
389 *
390 * Hermitian -- A = U D U'
391 *
392  CALL pzlaghe( m, llb, d, a, ia, ja, desca, iseed, order, work,
393  $ lwork, iinfo )
394 *
395  RETURN
396 *
397 * End of PZLATMS
398 *
399  END
max
#define max(A, B)
Definition: pcgemr.c:180
pzlaghe
subroutine pzlaghe(N, K, D, A, IA, JA, DESCA, ISEED, ORDER, WORK, LWORK, INFO)
Definition: pzlagsy.f:5
pchk1mat
subroutine pchk1mat(MA, MAPOS0, NA, NAPOS0, IA, JA, DESCA, DESCAPOS0, NEXTRA, EX, EXPOS, INFO)
Definition: pchkxmat.f:3
dlatm1
subroutine dlatm1(MODE, COND, IRSIGN, IDIST, ISEED, D, N, INFO)
Definition: dlatm1.f:2
chk1mat
subroutine chk1mat(MA, MAPOS0, NA, NAPOS0, IA, JA, DESCA, DESCAPOS0, INFO)
Definition: chk1mat.f:3
pzlatms
subroutine pzlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, IA, JA, DESCA, ORDER, WORK, LWORK, INFO)
Definition: pzlatms.f:6
pxerbla
subroutine pxerbla(ICTXT, SRNAME, INFO)
Definition: pxerbla.f:2
min
#define min(A, B)
Definition: pcgemr.c:181