LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
zlatme.f
Go to the documentation of this file.
1 *> \brief \b ZLATME
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE ZLATME( N, DIST, ISEED, D, MODE, COND, DMAX,
12 * RSIGN,
13 * UPPER, SIM, DS, MODES, CONDS, KL, KU, ANORM,
14 * A,
15 * LDA, WORK, INFO )
16 *
17 * .. Scalar Arguments ..
18 * CHARACTER DIST, RSIGN, SIM, UPPER
19 * INTEGER INFO, KL, KU, LDA, MODE, MODES, N
20 * DOUBLE PRECISION ANORM, COND, CONDS
21 * COMPLEX*16 DMAX
22 * ..
23 * .. Array Arguments ..
24 * INTEGER ISEED( 4 )
25 * DOUBLE PRECISION DS( * )
26 * COMPLEX*16 A( LDA, * ), D( * ), WORK( * )
27 * ..
28 *
29 *
30 *> \par Purpose:
31 * =============
32 *>
33 *> \verbatim
34 *>
35 *> ZLATME generates random non-symmetric square matrices with
36 *> specified eigenvalues for testing LAPACK programs.
37 *>
38 *> ZLATME operates by applying the following sequence of
39 *> operations:
40 *>
41 *> 1. Set the diagonal to D, where D may be input or
42 *> computed according to MODE, COND, DMAX, and RSIGN
43 *> as described below.
44 *>
45 *> 2. If UPPER='T', the upper triangle of A is set to random values
46 *> out of distribution DIST.
47 *>
48 *> 3. If SIM='T', A is multiplied on the left by a random matrix
49 *> X, whose singular values are specified by DS, MODES, and
50 *> CONDS, and on the right by X inverse.
51 *>
52 *> 4. If KL < N-1, the lower bandwidth is reduced to KL using
53 *> Householder transformations. If KU < N-1, the upper
54 *> bandwidth is reduced to KU.
55 *>
56 *> 5. If ANORM is not negative, the matrix is scaled to have
57 *> maximum-element-norm ANORM.
58 *>
59 *> (Note: since the matrix cannot be reduced beyond Hessenberg form,
60 *> no packing options are available.)
61 *> \endverbatim
62 *
63 * Arguments:
64 * ==========
65 *
66 *> \param[in] N
67 *> \verbatim
68 *> N is INTEGER
69 *> The number of columns (or rows) of A. Not modified.
70 *> \endverbatim
71 *>
72 *> \param[in] DIST
73 *> \verbatim
74 *> DIST is CHARACTER*1
75 *> On entry, DIST specifies the type of distribution to be used
76 *> to generate the random eigen-/singular values, and on the
77 *> upper triangle (see UPPER).
78 *> 'U' => UNIFORM( 0, 1 ) ( 'U' for uniform )
79 *> 'S' => UNIFORM( -1, 1 ) ( 'S' for symmetric )
80 *> 'N' => NORMAL( 0, 1 ) ( 'N' for normal )
81 *> 'D' => uniform on the complex disc |z| < 1.
82 *> Not modified.
83 *> \endverbatim
84 *>
85 *> \param[in,out] ISEED
86 *> \verbatim
87 *> ISEED is INTEGER array, dimension ( 4 )
88 *> On entry ISEED specifies the seed of the random number
89 *> generator. They should lie between 0 and 4095 inclusive,
90 *> and ISEED(4) should be odd. The random number generator
91 *> uses a linear congruential sequence limited to small
92 *> integers, and so should produce machine independent
93 *> random numbers. The values of ISEED are changed on
94 *> exit, and can be used in the next call to ZLATME
95 *> to continue the same random number sequence.
96 *> Changed on exit.
97 *> \endverbatim
98 *>
99 *> \param[in,out] D
100 *> \verbatim
101 *> D is COMPLEX*16 array, dimension ( N )
102 *> This array is used to specify the eigenvalues of A. If
103 *> MODE=0, then D is assumed to contain the eigenvalues
104 *> otherwise they will be computed according to MODE, COND,
105 *> DMAX, and RSIGN and placed in D.
106 *> Modified if MODE is nonzero.
107 *> \endverbatim
108 *>
109 *> \param[in] MODE
110 *> \verbatim
111 *> MODE is INTEGER
112 *> On entry this describes how the eigenvalues are to
113 *> be specified:
114 *> MODE = 0 means use D as input
115 *> MODE = 1 sets D(1)=1 and D(2:N)=1.0/COND
116 *> MODE = 2 sets D(1:N-1)=1 and D(N)=1.0/COND
117 *> MODE = 3 sets D(I)=COND**(-(I-1)/(N-1))
118 *> MODE = 4 sets D(i)=1 - (i-1)/(N-1)*(1 - 1/COND)
119 *> MODE = 5 sets D to random numbers in the range
120 *> ( 1/COND , 1 ) such that their logarithms
121 *> are uniformly distributed.
122 *> MODE = 6 set D to random numbers from same distribution
123 *> as the rest of the matrix.
124 *> MODE < 0 has the same meaning as ABS(MODE), except that
125 *> the order of the elements of D is reversed.
126 *> Thus if MODE is between 1 and 4, D has entries ranging
127 *> from 1 to 1/COND, if between -1 and -4, D has entries
128 *> ranging from 1/COND to 1,
129 *> Not modified.
130 *> \endverbatim
131 *>
132 *> \param[in] COND
133 *> \verbatim
134 *> COND is DOUBLE PRECISION
135 *> On entry, this is used as described under MODE above.
136 *> If used, it must be >= 1. Not modified.
137 *> \endverbatim
138 *>
139 *> \param[in] DMAX
140 *> \verbatim
141 *> DMAX is COMPLEX*16
142 *> If MODE is neither -6, 0 nor 6, the contents of D, as
143 *> computed according to MODE and COND, will be scaled by
144 *> DMAX / max(abs(D(i))). Note that DMAX need not be
145 *> positive or real: if DMAX is negative or complex (or zero),
146 *> D will be scaled by a negative or complex number (or zero).
147 *> If RSIGN='F' then the largest (absolute) eigenvalue will be
148 *> equal to DMAX.
149 *> Not modified.
150 *> \endverbatim
151 *>
152 *> \param[in] RSIGN
153 *> \verbatim
154 *> RSIGN is CHARACTER*1
155 *> If MODE is not 0, 6, or -6, and RSIGN='T', then the
156 *> elements of D, as computed according to MODE and COND, will
157 *> be multiplied by a random complex number from the unit
158 *> circle |z| = 1. If RSIGN='F', they will not be. RSIGN may
159 *> only have the values 'T' or 'F'.
160 *> Not modified.
161 *> \endverbatim
162 *>
163 *> \param[in] UPPER
164 *> \verbatim
165 *> UPPER is CHARACTER*1
166 *> If UPPER='T', then the elements of A above the diagonal
167 *> will be set to random numbers out of DIST. If UPPER='F',
168 *> they will not. UPPER may only have the values 'T' or 'F'.
169 *> Not modified.
170 *> \endverbatim
171 *>
172 *> \param[in] SIM
173 *> \verbatim
174 *> SIM is CHARACTER*1
175 *> If SIM='T', then A will be operated on by a "similarity
176 *> transform", i.e., multiplied on the left by a matrix X and
177 *> on the right by X inverse. X = U S V, where U and V are
178 *> random unitary matrices and S is a (diagonal) matrix of
179 *> singular values specified by DS, MODES, and CONDS. If
180 *> SIM='F', then A will not be transformed.
181 *> Not modified.
182 *> \endverbatim
183 *>
184 *> \param[in,out] DS
185 *> \verbatim
186 *> DS is DOUBLE PRECISION array, dimension ( N )
187 *> This array is used to specify the singular values of X,
188 *> in the same way that D specifies the eigenvalues of A.
189 *> If MODE=0, the DS contains the singular values, which
190 *> may not be zero.
191 *> Modified if MODE is nonzero.
192 *> \endverbatim
193 *>
194 *> \param[in] MODES
195 *> \verbatim
196 *> MODES is INTEGER
197 *> \endverbatim
198 *>
199 *> \param[in] CONDS
200 *> \verbatim
201 *> CONDS is DOUBLE PRECISION
202 *> Similar to MODE and COND, but for specifying the diagonal
203 *> of S. MODES=-6 and +6 are not allowed (since they would
204 *> result in randomly ill-conditioned eigenvalues.)
205 *> \endverbatim
206 *>
207 *> \param[in] KL
208 *> \verbatim
209 *> KL is INTEGER
210 *> This specifies the lower bandwidth of the matrix. KL=1
211 *> specifies upper Hessenberg form. If KL is at least N-1,
212 *> then A will have full lower bandwidth.
213 *> Not modified.
214 *> \endverbatim
215 *>
216 *> \param[in] KU
217 *> \verbatim
218 *> KU is INTEGER
219 *> This specifies the upper bandwidth of the matrix. KU=1
220 *> specifies lower Hessenberg form. If KU is at least N-1,
221 *> then A will have full upper bandwidth; if KU and KL
222 *> are both at least N-1, then A will be dense. Only one of
223 *> KU and KL may be less than N-1.
224 *> Not modified.
225 *> \endverbatim
226 *>
227 *> \param[in] ANORM
228 *> \verbatim
229 *> ANORM is DOUBLE PRECISION
230 *> If ANORM is not negative, then A will be scaled by a non-
231 *> negative real number to make the maximum-element-norm of A
232 *> to be ANORM.
233 *> Not modified.
234 *> \endverbatim
235 *>
236 *> \param[out] A
237 *> \verbatim
238 *> A is COMPLEX*16 array, dimension ( LDA, N )
239 *> On exit A is the desired test matrix.
240 *> Modified.
241 *> \endverbatim
242 *>
243 *> \param[in] LDA
244 *> \verbatim
245 *> LDA is INTEGER
246 *> LDA specifies the first dimension of A as declared in the
247 *> calling program. LDA must be at least M.
248 *> Not modified.
249 *> \endverbatim
250 *>
251 *> \param[out] WORK
252 *> \verbatim
253 *> WORK is COMPLEX*16 array, dimension ( 3*N )
254 *> Workspace.
255 *> Modified.
256 *> \endverbatim
257 *>
258 *> \param[out] INFO
259 *> \verbatim
260 *> INFO is INTEGER
261 *> Error code. On exit, INFO will be set to one of the
262 *> following values:
263 *> 0 => normal return
264 *> -1 => N negative
265 *> -2 => DIST illegal string
266 *> -5 => MODE not in range -6 to 6
267 *> -6 => COND less than 1.0, and MODE neither -6, 0 nor 6
268 *> -9 => RSIGN is not 'T' or 'F'
269 *> -10 => UPPER is not 'T' or 'F'
270 *> -11 => SIM is not 'T' or 'F'
271 *> -12 => MODES=0 and DS has a zero singular value.
272 *> -13 => MODES is not in the range -5 to 5.
273 *> -14 => MODES is nonzero and CONDS is less than 1.
274 *> -15 => KL is less than 1.
275 *> -16 => KU is less than 1, or KL and KU are both less than
276 *> N-1.
277 *> -19 => LDA is less than M.
278 *> 1 => Error return from ZLATM1 (computing D)
279 *> 2 => Cannot scale to DMAX (max. eigenvalue is 0)
280 *> 3 => Error return from DLATM1 (computing DS)
281 *> 4 => Error return from ZLARGE
282 *> 5 => Zero singular value from DLATM1.
283 *> \endverbatim
284 *
285 * Authors:
286 * ========
287 *
288 *> \author Univ. of Tennessee
289 *> \author Univ. of California Berkeley
290 *> \author Univ. of Colorado Denver
291 *> \author NAG Ltd.
292 *
293 *> \date November 2011
294 *
295 *> \ingroup complex16_matgen
296 *
297 * =====================================================================
298  SUBROUTINE zlatme( N, DIST, ISEED, D, MODE, COND, DMAX,
299  $ rsign,
300  $ upper, sim, ds, modes, conds, kl, ku, anorm,
301  $ a,
302  $ lda, work, info )
303 *
304 * -- LAPACK computational routine (version 3.4.0) --
305 * -- LAPACK is a software package provided by Univ. of Tennessee, --
306 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
307 * November 2011
308 *
309 * .. Scalar Arguments ..
310  CHARACTER dist, rsign, sim, upper
311  INTEGER info, kl, ku, lda, mode, modes, n
312  DOUBLE PRECISION anorm, cond, conds
313  COMPLEX*16 dmax
314 * ..
315 * .. Array Arguments ..
316  INTEGER iseed( 4 )
317  DOUBLE PRECISION ds( * )
318  COMPLEX*16 a( lda, * ), d( * ), work( * )
319 * ..
320 *
321 * =====================================================================
322 *
323 * .. Parameters ..
324  DOUBLE PRECISION zero
325  parameter( zero = 0.0d+0 )
326  DOUBLE PRECISION one
327  parameter( one = 1.0d+0 )
328  COMPLEX*16 czero
329  parameter( czero = ( 0.0d+0, 0.0d+0 ) )
330  COMPLEX*16 cone
331  parameter( cone = ( 1.0d+0, 0.0d+0 ) )
332 * ..
333 * .. Local Scalars ..
334  LOGICAL bads
335  INTEGER i, ic, icols, idist, iinfo, ir, irows, irsign,
336  $ isim, iupper, j, jc, jcr
337  DOUBLE PRECISION ralpha, temp
338  COMPLEX*16 alpha, tau, xnorms
339 * ..
340 * .. Local Arrays ..
341  DOUBLE PRECISION tempa( 1 )
342 * ..
343 * .. External Functions ..
344  LOGICAL lsame
345  DOUBLE PRECISION zlange
346  COMPLEX*16 zlarnd
347  EXTERNAL lsame, zlange, zlarnd
348 * ..
349 * .. External Subroutines ..
350  EXTERNAL dlatm1, xerbla, zcopy, zdscal, zgemv, zgerc,
352  $ zscal
353 * ..
354 * .. Intrinsic Functions ..
355  INTRINSIC abs, dconjg, max, mod
356 * ..
357 * .. Executable Statements ..
358 *
359 * 1) Decode and Test the input parameters.
360 * Initialize flags & seed.
361 *
362  info = 0
363 *
364 * Quick return if possible
365 *
366  IF( n.EQ.0 )
367  $ RETURN
368 *
369 * Decode DIST
370 *
371  IF( lsame( dist, 'U' ) ) THEN
372  idist = 1
373  ELSE IF( lsame( dist, 'S' ) ) THEN
374  idist = 2
375  ELSE IF( lsame( dist, 'N' ) ) THEN
376  idist = 3
377  ELSE IF( lsame( dist, 'D' ) ) THEN
378  idist = 4
379  ELSE
380  idist = -1
381  END IF
382 *
383 * Decode RSIGN
384 *
385  IF( lsame( rsign, 'T' ) ) THEN
386  irsign = 1
387  ELSE IF( lsame( rsign, 'F' ) ) THEN
388  irsign = 0
389  ELSE
390  irsign = -1
391  END IF
392 *
393 * Decode UPPER
394 *
395  IF( lsame( upper, 'T' ) ) THEN
396  iupper = 1
397  ELSE IF( lsame( upper, 'F' ) ) THEN
398  iupper = 0
399  ELSE
400  iupper = -1
401  END IF
402 *
403 * Decode SIM
404 *
405  IF( lsame( sim, 'T' ) ) THEN
406  isim = 1
407  ELSE IF( lsame( sim, 'F' ) ) THEN
408  isim = 0
409  ELSE
410  isim = -1
411  END IF
412 *
413 * Check DS, if MODES=0 and ISIM=1
414 *
415  bads = .false.
416  IF( modes.EQ.0 .AND. isim.EQ.1 ) THEN
417  DO 10 j = 1, n
418  IF( ds( j ).EQ.zero )
419  $ bads = .true.
420  10 CONTINUE
421  END IF
422 *
423 * Set INFO if an error
424 *
425  IF( n.LT.0 ) THEN
426  info = -1
427  ELSE IF( idist.EQ.-1 ) THEN
428  info = -2
429  ELSE IF( abs( mode ).GT.6 ) THEN
430  info = -5
431  ELSE IF( ( mode.NE.0 .AND. abs( mode ).NE.6 ) .AND. cond.LT.one )
432  $ THEN
433  info = -6
434  ELSE IF( irsign.EQ.-1 ) THEN
435  info = -9
436  ELSE IF( iupper.EQ.-1 ) THEN
437  info = -10
438  ELSE IF( isim.EQ.-1 ) THEN
439  info = -11
440  ELSE IF( bads ) THEN
441  info = -12
442  ELSE IF( isim.EQ.1 .AND. abs( modes ).GT.5 ) THEN
443  info = -13
444  ELSE IF( isim.EQ.1 .AND. modes.NE.0 .AND. conds.LT.one ) THEN
445  info = -14
446  ELSE IF( kl.LT.1 ) THEN
447  info = -15
448  ELSE IF( ku.LT.1 .OR. ( ku.LT.n-1 .AND. kl.LT.n-1 ) ) THEN
449  info = -16
450  ELSE IF( lda.LT.max( 1, n ) ) THEN
451  info = -19
452  END IF
453 *
454  IF( info.NE.0 ) THEN
455  CALL xerbla( 'ZLATME', -info )
456  RETURN
457  END IF
458 *
459 * Initialize random number generator
460 *
461  DO 20 i = 1, 4
462  iseed( i ) = mod( abs( iseed( i ) ), 4096 )
463  20 CONTINUE
464 *
465  IF( mod( iseed( 4 ), 2 ).NE.1 )
466  $ iseed( 4 ) = iseed( 4 ) + 1
467 *
468 * 2) Set up diagonal of A
469 *
470 * Compute D according to COND and MODE
471 *
472  CALL zlatm1( mode, cond, irsign, idist, iseed, d, n, iinfo )
473  IF( iinfo.NE.0 ) THEN
474  info = 1
475  RETURN
476  END IF
477  IF( mode.NE.0 .AND. abs( mode ).NE.6 ) THEN
478 *
479 * Scale by DMAX
480 *
481  temp = abs( d( 1 ) )
482  DO 30 i = 2, n
483  temp = max( temp, abs( d( i ) ) )
484  30 CONTINUE
485 *
486  IF( temp.GT.zero ) THEN
487  alpha = dmax / temp
488  ELSE
489  info = 2
490  RETURN
491  END IF
492 *
493  CALL zscal( n, alpha, d, 1 )
494 *
495  END IF
496 *
497  CALL zlaset( 'Full', n, n, czero, czero, a, lda )
498  CALL zcopy( n, d, 1, a, lda+1 )
499 *
500 * 3) If UPPER='T', set upper triangle of A to random numbers.
501 *
502  IF( iupper.NE.0 ) THEN
503  DO 40 jc = 2, n
504  CALL zlarnv( idist, iseed, jc-1, a( 1, jc ) )
505  40 CONTINUE
506  END IF
507 *
508 * 4) If SIM='T', apply similarity transformation.
509 *
510 * -1
511 * Transform is X A X , where X = U S V, thus
512 *
513 * it is U S V A V' (1/S) U'
514 *
515  IF( isim.NE.0 ) THEN
516 *
517 * Compute S (singular values of the eigenvector matrix)
518 * according to CONDS and MODES
519 *
520  CALL dlatm1( modes, conds, 0, 0, iseed, ds, n, iinfo )
521  IF( iinfo.NE.0 ) THEN
522  info = 3
523  RETURN
524  END IF
525 *
526 * Multiply by V and V'
527 *
528  CALL zlarge( n, a, lda, iseed, work, iinfo )
529  IF( iinfo.NE.0 ) THEN
530  info = 4
531  RETURN
532  END IF
533 *
534 * Multiply by S and (1/S)
535 *
536  DO 50 j = 1, n
537  CALL zdscal( n, ds( j ), a( j, 1 ), lda )
538  IF( ds( j ).NE.zero ) THEN
539  CALL zdscal( n, one / ds( j ), a( 1, j ), 1 )
540  ELSE
541  info = 5
542  RETURN
543  END IF
544  50 CONTINUE
545 *
546 * Multiply by U and U'
547 *
548  CALL zlarge( n, a, lda, iseed, work, iinfo )
549  IF( iinfo.NE.0 ) THEN
550  info = 4
551  RETURN
552  END IF
553  END IF
554 *
555 * 5) Reduce the bandwidth.
556 *
557  IF( kl.LT.n-1 ) THEN
558 *
559 * Reduce bandwidth -- kill column
560 *
561  DO 60 jcr = kl + 1, n - 1
562  ic = jcr - kl
563  irows = n + 1 - jcr
564  icols = n + kl - jcr
565 *
566  CALL zcopy( irows, a( jcr, ic ), 1, work, 1 )
567  xnorms = work( 1 )
568  CALL zlarfg( irows, xnorms, work( 2 ), 1, tau )
569  tau = dconjg( tau )
570  work( 1 ) = cone
571  alpha = zlarnd( 5, iseed )
572 *
573  CALL zgemv( 'C', irows, icols, cone, a( jcr, ic+1 ), lda,
574  $ work, 1, czero, work( irows+1 ), 1 )
575  CALL zgerc( irows, icols, -tau, work, 1, work( irows+1 ), 1,
576  $ a( jcr, ic+1 ), lda )
577 *
578  CALL zgemv( 'N', n, irows, cone, a( 1, jcr ), lda, work, 1,
579  $ czero, work( irows+1 ), 1 )
580  CALL zgerc( n, irows, -dconjg( tau ), work( irows+1 ), 1,
581  $ work, 1, a( 1, jcr ), lda )
582 *
583  a( jcr, ic ) = xnorms
584  CALL zlaset( 'Full', irows-1, 1, czero, czero,
585  $ a( jcr+1, ic ), lda )
586 *
587  CALL zscal( icols+1, alpha, a( jcr, ic ), lda )
588  CALL zscal( n, dconjg( alpha ), a( 1, jcr ), 1 )
589  60 CONTINUE
590  ELSE IF( ku.LT.n-1 ) THEN
591 *
592 * Reduce upper bandwidth -- kill a row at a time.
593 *
594  DO 70 jcr = ku + 1, n - 1
595  ir = jcr - ku
596  irows = n + ku - jcr
597  icols = n + 1 - jcr
598 *
599  CALL zcopy( icols, a( ir, jcr ), lda, work, 1 )
600  xnorms = work( 1 )
601  CALL zlarfg( icols, xnorms, work( 2 ), 1, tau )
602  tau = dconjg( tau )
603  work( 1 ) = cone
604  CALL zlacgv( icols-1, work( 2 ), 1 )
605  alpha = zlarnd( 5, iseed )
606 *
607  CALL zgemv( 'N', irows, icols, cone, a( ir+1, jcr ), lda,
608  $ work, 1, czero, work( icols+1 ), 1 )
609  CALL zgerc( irows, icols, -tau, work( icols+1 ), 1, work, 1,
610  $ a( ir+1, jcr ), lda )
611 *
612  CALL zgemv( 'C', icols, n, cone, a( jcr, 1 ), lda, work, 1,
613  $ czero, work( icols+1 ), 1 )
614  CALL zgerc( icols, n, -dconjg( tau ), work, 1,
615  $ work( icols+1 ), 1, a( jcr, 1 ), lda )
616 *
617  a( ir, jcr ) = xnorms
618  CALL zlaset( 'Full', 1, icols-1, czero, czero,
619  $ a( ir, jcr+1 ), lda )
620 *
621  CALL zscal( irows+1, alpha, a( ir, jcr ), 1 )
622  CALL zscal( n, dconjg( alpha ), a( jcr, 1 ), lda )
623  70 CONTINUE
624  END IF
625 *
626 * Scale the matrix to have norm ANORM
627 *
628  IF( anorm.GE.zero ) THEN
629  temp = zlange( 'M', n, n, a, lda, tempa )
630  IF( temp.GT.zero ) THEN
631  ralpha = anorm / temp
632  DO 80 j = 1, n
633  CALL zdscal( n, ralpha, a( 1, j ), 1 )
634  80 CONTINUE
635  END IF
636  END IF
637 *
638  RETURN
639 *
640 * End of ZLATME
641 *
642  END