LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlatmr.f
Go to the documentation of this file.
1*> \brief \b DLATMR
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 DLATMR( M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX,
12* RSIGN, GRADE, DL, MODEL, CONDL, DR, MODER,
13* CONDR, PIVTNG, IPIVOT, KL, KU, SPARSE, ANORM,
14* PACK, A, LDA, IWORK, INFO )
15*
16* .. Scalar Arguments ..
17* CHARACTER DIST, GRADE, PACK, PIVTNG, RSIGN, SYM
18* INTEGER INFO, KL, KU, LDA, M, MODE, MODEL, MODER, N
19* DOUBLE PRECISION ANORM, COND, CONDL, CONDR, DMAX, SPARSE
20* ..
21* .. Array Arguments ..
22* INTEGER IPIVOT( * ), ISEED( 4 ), IWORK( * )
23* DOUBLE PRECISION A( LDA, * ), D( * ), DL( * ), DR( * )
24* ..
25*
26*
27*> \par Purpose:
28* =============
29*>
30*> \verbatim
31*>
32*> DLATMR generates random matrices of various types for testing
33*> LAPACK programs.
34*>
35*> DLATMR operates by applying the following sequence of
36*> operations:
37*>
38*> Generate a matrix A with random entries of distribution DIST
39*> which is symmetric if SYM='S', and nonsymmetric
40*> if SYM='N'.
41*>
42*> Set the diagonal to D, where D may be input or
43*> computed according to MODE, COND, DMAX and RSIGN
44*> as described below.
45*>
46*> Grade the matrix, if desired, from the left and/or right
47*> as specified by GRADE. The inputs DL, MODEL, CONDL, DR,
48*> MODER and CONDR also determine the grading as described
49*> below.
50*>
51*> Permute, if desired, the rows and/or columns as specified by
52*> PIVTNG and IPIVOT.
53*>
54*> Set random entries to zero, if desired, to get a random sparse
55*> matrix as specified by SPARSE.
56*>
57*> Make A a band matrix, if desired, by zeroing out the matrix
58*> outside a band of lower bandwidth KL and upper bandwidth KU.
59*>
60*> Scale A, if desired, to have maximum entry ANORM.
61*>
62*> Pack the matrix if desired. Options specified by PACK are:
63*> no packing
64*> zero out upper half (if symmetric)
65*> zero out lower half (if symmetric)
66*> store the upper half columnwise (if symmetric or
67*> square upper triangular)
68*> store the lower half columnwise (if symmetric or
69*> square lower triangular)
70*> same as upper half rowwise if symmetric
71*> store the lower triangle in banded format (if symmetric)
72*> store the upper triangle in banded format (if symmetric)
73*> store the entire matrix in banded format
74*>
75*> Note: If two calls to DLATMR differ only in the PACK parameter,
76*> they will generate mathematically equivalent matrices.
77*>
78*> If two calls to DLATMR both have full bandwidth (KL = M-1
79*> and KU = N-1), and differ only in the PIVTNG and PACK
80*> parameters, then the matrices generated will differ only
81*> in the order of the rows and/or columns, and otherwise
82*> contain the same data. This consistency cannot be and
83*> is not maintained with less than full bandwidth.
84*> \endverbatim
85*
86* Arguments:
87* ==========
88*
89*> \param[in] M
90*> \verbatim
91*> M is INTEGER
92*> Number of rows of A. Not modified.
93*> \endverbatim
94*>
95*> \param[in] N
96*> \verbatim
97*> N is INTEGER
98*> Number of columns of A. Not modified.
99*> \endverbatim
100*>
101*> \param[in] DIST
102*> \verbatim
103*> DIST is CHARACTER*1
104*> On entry, DIST specifies the type of distribution to be used
105*> to generate a random matrix .
106*> 'U' => UNIFORM( 0, 1 ) ( 'U' for uniform )
107*> 'S' => UNIFORM( -1, 1 ) ( 'S' for symmetric )
108*> 'N' => NORMAL( 0, 1 ) ( 'N' for normal )
109*> Not modified.
110*> \endverbatim
111*>
112*> \param[in,out] ISEED
113*> \verbatim
114*> ISEED is INTEGER array, dimension (4)
115*> On entry ISEED specifies the seed of the random number
116*> generator. They should lie between 0 and 4095 inclusive,
117*> and ISEED(4) should be odd. The random number generator
118*> uses a linear congruential sequence limited to small
119*> integers, and so should produce machine independent
120*> random numbers. The values of ISEED are changed on
121*> exit, and can be used in the next call to DLATMR
122*> to continue the same random number sequence.
123*> Changed on exit.
124*> \endverbatim
125*>
126*> \param[in] SYM
127*> \verbatim
128*> SYM is CHARACTER*1
129*> If SYM='S' or 'H', generated matrix is symmetric.
130*> If SYM='N', generated matrix is nonsymmetric.
131*> Not modified.
132*> \endverbatim
133*>
134*> \param[in,out] D
135*> \verbatim
136*> D is DOUBLE PRECISION array, dimension (min(M,N))
137*> On entry this array specifies the diagonal entries
138*> of the diagonal of A. D may either be specified
139*> on entry, or set according to MODE and COND as described
140*> below. May be changed on exit if MODE is nonzero.
141*> \endverbatim
142*>
143*> \param[in] MODE
144*> \verbatim
145*> MODE is INTEGER
146*> On entry describes how D is to be used:
147*> MODE = 0 means use D as input
148*> MODE = 1 sets D(1)=1 and D(2:N)=1.0/COND
149*> MODE = 2 sets D(1:N-1)=1 and D(N)=1.0/COND
150*> MODE = 3 sets D(I)=COND**(-(I-1)/(N-1))
151*> MODE = 4 sets D(i)=1 - (i-1)/(N-1)*(1 - 1/COND)
152*> MODE = 5 sets D to random numbers in the range
153*> ( 1/COND , 1 ) such that their logarithms
154*> are uniformly distributed.
155*> MODE = 6 set D to random numbers from same distribution
156*> as the rest of the matrix.
157*> MODE < 0 has the same meaning as ABS(MODE), except that
158*> the order of the elements of D is reversed.
159*> Thus if MODE is positive, D has entries ranging from
160*> 1 to 1/COND, if negative, from 1/COND to 1,
161*> Not modified.
162*> \endverbatim
163*>
164*> \param[in] COND
165*> \verbatim
166*> COND is DOUBLE PRECISION
167*> On entry, used as described under MODE above.
168*> If used, it must be >= 1. Not modified.
169*> \endverbatim
170*>
171*> \param[in] DMAX
172*> \verbatim
173*> DMAX is DOUBLE PRECISION
174*> If MODE neither -6, 0 nor 6, the diagonal is scaled by
175*> DMAX / max(abs(D(i))), so that maximum absolute entry
176*> of diagonal is abs(DMAX). If DMAX is negative (or zero),
177*> diagonal will be scaled by a negative number (or zero).
178*> \endverbatim
179*>
180*> \param[in] RSIGN
181*> \verbatim
182*> RSIGN is CHARACTER*1
183*> If MODE neither -6, 0 nor 6, specifies sign of diagonal
184*> as follows:
185*> 'T' => diagonal entries are multiplied by 1 or -1
186*> with probability .5
187*> 'F' => diagonal unchanged
188*> Not modified.
189*> \endverbatim
190*>
191*> \param[in] GRADE
192*> \verbatim
193*> GRADE is CHARACTER*1
194*> Specifies grading of matrix as follows:
195*> 'N' => no grading
196*> 'L' => matrix premultiplied by diag( DL )
197*> (only if matrix nonsymmetric)
198*> 'R' => matrix postmultiplied by diag( DR )
199*> (only if matrix nonsymmetric)
200*> 'B' => matrix premultiplied by diag( DL ) and
201*> postmultiplied by diag( DR )
202*> (only if matrix nonsymmetric)
203*> 'S' or 'H' => matrix premultiplied by diag( DL ) and
204*> postmultiplied by diag( DL )
205*> ('S' for symmetric, or 'H' for Hermitian)
206*> 'E' => matrix premultiplied by diag( DL ) and
207*> postmultiplied by inv( diag( DL ) )
208*> ( 'E' for eigenvalue invariance)
209*> (only if matrix nonsymmetric)
210*> Note: if GRADE='E', then M must equal N.
211*> Not modified.
212*> \endverbatim
213*>
214*> \param[in,out] DL
215*> \verbatim
216*> DL is DOUBLE PRECISION array, dimension (M)
217*> If MODEL=0, then on entry this array specifies the diagonal
218*> entries of a diagonal matrix used as described under GRADE
219*> above. If MODEL is not zero, then DL will be set according
220*> to MODEL and CONDL, analogous to the way D is set according
221*> to MODE and COND (except there is no DMAX parameter for DL).
222*> If GRADE='E', then DL cannot have zero entries.
223*> Not referenced if GRADE = 'N' or 'R'. Changed on exit.
224*> \endverbatim
225*>
226*> \param[in] MODEL
227*> \verbatim
228*> MODEL is INTEGER
229*> This specifies how the diagonal array DL is to be computed,
230*> just as MODE specifies how D is to be computed.
231*> Not modified.
232*> \endverbatim
233*>
234*> \param[in] CONDL
235*> \verbatim
236*> CONDL is DOUBLE PRECISION
237*> When MODEL is not zero, this specifies the condition number
238*> of the computed DL. Not modified.
239*> \endverbatim
240*>
241*> \param[in,out] DR
242*> \verbatim
243*> DR is DOUBLE PRECISION array, dimension (N)
244*> If MODER=0, then on entry this array specifies the diagonal
245*> entries of a diagonal matrix used as described under GRADE
246*> above. If MODER is not zero, then DR will be set according
247*> to MODER and CONDR, analogous to the way D is set according
248*> to MODE and COND (except there is no DMAX parameter for DR).
249*> Not referenced if GRADE = 'N', 'L', 'H', 'S' or 'E'.
250*> Changed on exit.
251*> \endverbatim
252*>
253*> \param[in] MODER
254*> \verbatim
255*> MODER is INTEGER
256*> This specifies how the diagonal array DR is to be computed,
257*> just as MODE specifies how D is to be computed.
258*> Not modified.
259*> \endverbatim
260*>
261*> \param[in] CONDR
262*> \verbatim
263*> CONDR is DOUBLE PRECISION
264*> When MODER is not zero, this specifies the condition number
265*> of the computed DR. Not modified.
266*> \endverbatim
267*>
268*> \param[in] PIVTNG
269*> \verbatim
270*> PIVTNG is CHARACTER*1
271*> On entry specifies pivoting permutations as follows:
272*> 'N' or ' ' => none.
273*> 'L' => left or row pivoting (matrix must be nonsymmetric).
274*> 'R' => right or column pivoting (matrix must be
275*> nonsymmetric).
276*> 'B' or 'F' => both or full pivoting, i.e., on both sides.
277*> In this case, M must equal N
278*>
279*> If two calls to DLATMR both have full bandwidth (KL = M-1
280*> and KU = N-1), and differ only in the PIVTNG and PACK
281*> parameters, then the matrices generated will differ only
282*> in the order of the rows and/or columns, and otherwise
283*> contain the same data. This consistency cannot be
284*> maintained with less than full bandwidth.
285*> \endverbatim
286*>
287*> \param[in] IPIVOT
288*> \verbatim
289*> IPIVOT is INTEGER array, dimension (N or M)
290*> This array specifies the permutation used. After the
291*> basic matrix is generated, the rows, columns, or both
292*> are permuted. If, say, row pivoting is selected, DLATMR
293*> starts with the *last* row and interchanges the M-th and
294*> IPIVOT(M)-th rows, then moves to the next-to-last row,
295*> interchanging the (M-1)-th and the IPIVOT(M-1)-th rows,
296*> and so on. In terms of "2-cycles", the permutation is
297*> (1 IPIVOT(1)) (2 IPIVOT(2)) ... (M IPIVOT(M))
298*> where the rightmost cycle is applied first. This is the
299*> *inverse* of the effect of pivoting in LINPACK. The idea
300*> is that factoring (with pivoting) an identity matrix
301*> which has been inverse-pivoted in this way should
302*> result in a pivot vector identical to IPIVOT.
303*> Not referenced if PIVTNG = 'N'. Not modified.
304*> \endverbatim
305*>
306*> \param[in] KL
307*> \verbatim
308*> KL is INTEGER
309*> On entry specifies the lower bandwidth of the matrix. For
310*> example, KL=0 implies upper triangular, KL=1 implies upper
311*> Hessenberg, and KL at least M-1 implies the matrix is not
312*> banded. Must equal KU if matrix is symmetric.
313*> Not modified.
314*> \endverbatim
315*>
316*> \param[in] KU
317*> \verbatim
318*> KU is INTEGER
319*> On entry specifies the upper bandwidth of the matrix. For
320*> example, KU=0 implies lower triangular, KU=1 implies lower
321*> Hessenberg, and KU at least N-1 implies the matrix is not
322*> banded. Must equal KL if matrix is symmetric.
323*> Not modified.
324*> \endverbatim
325*>
326*> \param[in] SPARSE
327*> \verbatim
328*> SPARSE is DOUBLE PRECISION
329*> On entry specifies the sparsity of the matrix if a sparse
330*> matrix is to be generated. SPARSE should lie between
331*> 0 and 1. To generate a sparse matrix, for each matrix entry
332*> a uniform ( 0, 1 ) random number x is generated and
333*> compared to SPARSE; if x is larger the matrix entry
334*> is unchanged and if x is smaller the entry is set
335*> to zero. Thus on the average a fraction SPARSE of the
336*> entries will be set to zero.
337*> Not modified.
338*> \endverbatim
339*>
340*> \param[in] ANORM
341*> \verbatim
342*> ANORM is DOUBLE PRECISION
343*> On entry specifies maximum entry of output matrix
344*> (output matrix will by multiplied by a constant so that
345*> its largest absolute entry equal ANORM)
346*> if ANORM is nonnegative. If ANORM is negative no scaling
347*> is done. Not modified.
348*> \endverbatim
349*>
350*> \param[in] PACK
351*> \verbatim
352*> PACK is CHARACTER*1
353*> On entry specifies packing of matrix as follows:
354*> 'N' => no packing
355*> 'U' => zero out all subdiagonal entries (if symmetric)
356*> 'L' => zero out all superdiagonal entries (if symmetric)
357*> 'C' => store the upper triangle columnwise
358*> (only if matrix symmetric or square upper triangular)
359*> 'R' => store the lower triangle columnwise
360*> (only if matrix symmetric or square lower triangular)
361*> (same as upper half rowwise if symmetric)
362*> 'B' => store the lower triangle in band storage scheme
363*> (only if matrix symmetric)
364*> 'Q' => store the upper triangle in band storage scheme
365*> (only if matrix symmetric)
366*> 'Z' => store the entire matrix in band storage scheme
367*> (pivoting can be provided for by using this
368*> option to store A in the trailing rows of
369*> the allocated storage)
370*>
371*> Using these options, the various LAPACK packed and banded
372*> storage schemes can be obtained:
373*> GB - use 'Z'
374*> PB, SB or TB - use 'B' or 'Q'
375*> PP, SP or TP - use 'C' or 'R'
376*>
377*> If two calls to DLATMR differ only in the PACK parameter,
378*> they will generate mathematically equivalent matrices.
379*> Not modified.
380*> \endverbatim
381*>
382*> \param[out] A
383*> \verbatim
384*> A is DOUBLE PRECISION array, dimension (LDA,N)
385*> On exit A is the desired test matrix. Only those
386*> entries of A which are significant on output
387*> will be referenced (even if A is in packed or band
388*> storage format). The 'unoccupied corners' of A in
389*> band format will be zeroed out.
390*> \endverbatim
391*>
392*> \param[in] LDA
393*> \verbatim
394*> LDA is INTEGER
395*> on entry LDA specifies the first dimension of A as
396*> declared in the calling program.
397*> If PACK='N', 'U' or 'L', LDA must be at least max ( 1, M ).
398*> If PACK='C' or 'R', LDA must be at least 1.
399*> If PACK='B', or 'Q', LDA must be MIN ( KU+1, N )
400*> If PACK='Z', LDA must be at least KUU+KLL+1, where
401*> KUU = MIN ( KU, N-1 ) and KLL = MIN ( KL, M-1 )
402*> Not modified.
403*> \endverbatim
404*>
405*> \param[out] IWORK
406*> \verbatim
407*> IWORK is INTEGER array, dimension ( N or M)
408*> Workspace. Not referenced if PIVTNG = 'N'. Changed on exit.
409*> \endverbatim
410*>
411*> \param[out] INFO
412*> \verbatim
413*> INFO is INTEGER
414*> Error parameter on exit:
415*> 0 => normal return
416*> -1 => M negative or unequal to N and SYM='S' or 'H'
417*> -2 => N negative
418*> -3 => DIST illegal string
419*> -5 => SYM illegal string
420*> -7 => MODE not in range -6 to 6
421*> -8 => COND less than 1.0, and MODE neither -6, 0 nor 6
422*> -10 => MODE neither -6, 0 nor 6 and RSIGN illegal string
423*> -11 => GRADE illegal string, or GRADE='E' and
424*> M not equal to N, or GRADE='L', 'R', 'B' or 'E' and
425*> SYM = 'S' or 'H'
426*> -12 => GRADE = 'E' and DL contains zero
427*> -13 => MODEL not in range -6 to 6 and GRADE= 'L', 'B', 'H',
428*> 'S' or 'E'
429*> -14 => CONDL less than 1.0, GRADE='L', 'B', 'H', 'S' or 'E',
430*> and MODEL neither -6, 0 nor 6
431*> -16 => MODER not in range -6 to 6 and GRADE= 'R' or 'B'
432*> -17 => CONDR less than 1.0, GRADE='R' or 'B', and
433*> MODER neither -6, 0 nor 6
434*> -18 => PIVTNG illegal string, or PIVTNG='B' or 'F' and
435*> M not equal to N, or PIVTNG='L' or 'R' and SYM='S'
436*> or 'H'
437*> -19 => IPIVOT contains out of range number and
438*> PIVTNG not equal to 'N'
439*> -20 => KL negative
440*> -21 => KU negative, or SYM='S' or 'H' and KU not equal to KL
441*> -22 => SPARSE not in range 0. to 1.
442*> -24 => PACK illegal string, or PACK='U', 'L', 'B' or 'Q'
443*> and SYM='N', or PACK='C' and SYM='N' and either KL
444*> not equal to 0 or N not equal to M, or PACK='R' and
445*> SYM='N', and either KU not equal to 0 or N not equal
446*> to M
447*> -26 => LDA too small
448*> 1 => Error return from DLATM1 (computing D)
449*> 2 => Cannot scale diagonal to DMAX (max. entry is 0)
450*> 3 => Error return from DLATM1 (computing DL)
451*> 4 => Error return from DLATM1 (computing DR)
452*> 5 => ANORM is positive, but matrix constructed prior to
453*> attempting to scale it to have norm ANORM, is zero
454*> \endverbatim
455*
456* Authors:
457* ========
458*
459*> \author Univ. of Tennessee
460*> \author Univ. of California Berkeley
461*> \author Univ. of Colorado Denver
462*> \author NAG Ltd.
463*
464*> \ingroup double_matgen
465*
466* =====================================================================
467 SUBROUTINE dlatmr( M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX,
468 $ RSIGN, GRADE, DL, MODEL, CONDL, DR, MODER,
469 $ CONDR, PIVTNG, IPIVOT, KL, KU, SPARSE, ANORM,
470 $ PACK, A, LDA, IWORK, INFO )
471*
472* -- LAPACK computational routine --
473* -- LAPACK is a software package provided by Univ. of Tennessee, --
474* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
475*
476* .. Scalar Arguments ..
477 CHARACTER DIST, GRADE, PACK, PIVTNG, RSIGN, SYM
478 INTEGER INFO, KL, KU, LDA, M, MODE, MODEL, MODER, N
479 DOUBLE PRECISION ANORM, COND, CONDL, CONDR, DMAX, SPARSE
480* ..
481* .. Array Arguments ..
482 INTEGER IPIVOT( * ), ISEED( 4 ), IWORK( * )
483 DOUBLE PRECISION A( LDA, * ), D( * ), DL( * ), DR( * )
484* ..
485*
486* =====================================================================
487*
488* .. Parameters ..
489 DOUBLE PRECISION ZERO
490 PARAMETER ( ZERO = 0.0d0 )
491 DOUBLE PRECISION ONE
492 parameter( one = 1.0d0 )
493* ..
494* .. Local Scalars ..
495 LOGICAL BADPVT, DZERO, FULBND
496 INTEGER I, IDIST, IGRADE, IISUB, IPACK, IPVTNG, IRSIGN,
497 $ ISUB, ISYM, J, JJSUB, JSUB, K, KLL, KUU, MNMIN,
498 $ mnsub, mxsub, npvts
499 DOUBLE PRECISION ALPHA, ONORM, TEMP
500* ..
501* .. Local Arrays ..
502 DOUBLE PRECISION TEMPA( 1 )
503* ..
504* .. External Functions ..
505 LOGICAL LSAME
506 DOUBLE PRECISION DLANGB, DLANGE, DLANSB, DLANSP, DLANSY, DLATM2,
507 $ DLATM3
508 EXTERNAL lsame, dlangb, dlange, dlansb, dlansp, dlansy,
509 $ dlatm2, dlatm3
510* ..
511* .. External Subroutines ..
512 EXTERNAL dlatm1, dscal, xerbla
513* ..
514* .. Intrinsic Functions ..
515 INTRINSIC abs, max, min, mod
516* ..
517* .. Executable Statements ..
518*
519* 1) Decode and Test the input parameters.
520* Initialize flags & seed.
521*
522 info = 0
523*
524* Quick return if possible
525*
526 IF( m.EQ.0 .OR. n.EQ.0 )
527 $ RETURN
528*
529* Decode DIST
530*
531 IF( lsame( dist, 'U' ) ) THEN
532 idist = 1
533 ELSE IF( lsame( dist, 'S' ) ) THEN
534 idist = 2
535 ELSE IF( lsame( dist, 'N' ) ) THEN
536 idist = 3
537 ELSE
538 idist = -1
539 END IF
540*
541* Decode SYM
542*
543 IF( lsame( sym, 'S' ) ) THEN
544 isym = 0
545 ELSE IF( lsame( sym, 'N' ) ) THEN
546 isym = 1
547 ELSE IF( lsame( sym, 'H' ) ) THEN
548 isym = 0
549 ELSE
550 isym = -1
551 END IF
552*
553* Decode RSIGN
554*
555 IF( lsame( rsign, 'F' ) ) THEN
556 irsign = 0
557 ELSE IF( lsame( rsign, 'T' ) ) THEN
558 irsign = 1
559 ELSE
560 irsign = -1
561 END IF
562*
563* Decode PIVTNG
564*
565 IF( lsame( pivtng, 'N' ) ) THEN
566 ipvtng = 0
567 ELSE IF( lsame( pivtng, ' ' ) ) THEN
568 ipvtng = 0
569 ELSE IF( lsame( pivtng, 'L' ) ) THEN
570 ipvtng = 1
571 npvts = m
572 ELSE IF( lsame( pivtng, 'R' ) ) THEN
573 ipvtng = 2
574 npvts = n
575 ELSE IF( lsame( pivtng, 'B' ) ) THEN
576 ipvtng = 3
577 npvts = min( n, m )
578 ELSE IF( lsame( pivtng, 'F' ) ) THEN
579 ipvtng = 3
580 npvts = min( n, m )
581 ELSE
582 ipvtng = -1
583 END IF
584*
585* Decode GRADE
586*
587 IF( lsame( grade, 'N' ) ) THEN
588 igrade = 0
589 ELSE IF( lsame( grade, 'L' ) ) THEN
590 igrade = 1
591 ELSE IF( lsame( grade, 'R' ) ) THEN
592 igrade = 2
593 ELSE IF( lsame( grade, 'B' ) ) THEN
594 igrade = 3
595 ELSE IF( lsame( grade, 'E' ) ) THEN
596 igrade = 4
597 ELSE IF( lsame( grade, 'H' ) .OR. lsame( grade, 'S' ) ) THEN
598 igrade = 5
599 ELSE
600 igrade = -1
601 END IF
602*
603* Decode PACK
604*
605 IF( lsame( pack, 'N' ) ) THEN
606 ipack = 0
607 ELSE IF( lsame( pack, 'U' ) ) THEN
608 ipack = 1
609 ELSE IF( lsame( pack, 'L' ) ) THEN
610 ipack = 2
611 ELSE IF( lsame( pack, 'C' ) ) THEN
612 ipack = 3
613 ELSE IF( lsame( pack, 'R' ) ) THEN
614 ipack = 4
615 ELSE IF( lsame( pack, 'B' ) ) THEN
616 ipack = 5
617 ELSE IF( lsame( pack, 'Q' ) ) THEN
618 ipack = 6
619 ELSE IF( lsame( pack, 'Z' ) ) THEN
620 ipack = 7
621 ELSE
622 ipack = -1
623 END IF
624*
625* Set certain internal parameters
626*
627 mnmin = min( m, n )
628 kll = min( kl, m-1 )
629 kuu = min( ku, n-1 )
630*
631* If inv(DL) is used, check to see if DL has a zero entry.
632*
633 dzero = .false.
634 IF( igrade.EQ.4 .AND. model.EQ.0 ) THEN
635 DO 10 i = 1, m
636 IF( dl( i ).EQ.zero )
637 $ dzero = .true.
638 10 CONTINUE
639 END IF
640*
641* Check values in IPIVOT
642*
643 badpvt = .false.
644 IF( ipvtng.GT.0 ) THEN
645 DO 20 j = 1, npvts
646 IF( ipivot( j ).LE.0 .OR. ipivot( j ).GT.npvts )
647 $ badpvt = .true.
648 20 CONTINUE
649 END IF
650*
651* Set INFO if an error
652*
653 IF( m.LT.0 ) THEN
654 info = -1
655 ELSE IF( m.NE.n .AND. isym.EQ.0 ) THEN
656 info = -1
657 ELSE IF( n.LT.0 ) THEN
658 info = -2
659 ELSE IF( idist.EQ.-1 ) THEN
660 info = -3
661 ELSE IF( isym.EQ.-1 ) THEN
662 info = -5
663 ELSE IF( mode.LT.-6 .OR. mode.GT.6 ) THEN
664 info = -7
665 ELSE IF( ( mode.NE.-6 .AND. mode.NE.0 .AND. mode.NE.6 ) .AND.
666 $ cond.LT.one ) THEN
667 info = -8
668 ELSE IF( ( mode.NE.-6 .AND. mode.NE.0 .AND. mode.NE.6 ) .AND.
669 $ irsign.EQ.-1 ) THEN
670 info = -10
671 ELSE IF( igrade.EQ.-1 .OR. ( igrade.EQ.4 .AND. m.NE.n ) .OR.
672 $ ( ( igrade.GE.1 .AND. igrade.LE.4 ) .AND. isym.EQ.0 ) )
673 $ THEN
674 info = -11
675 ELSE IF( igrade.EQ.4 .AND. dzero ) THEN
676 info = -12
677 ELSE IF( ( igrade.EQ.1 .OR. igrade.EQ.3 .OR. igrade.EQ.4 .OR.
678 $ igrade.EQ.5 ) .AND. ( model.LT.-6 .OR. model.GT.6 ) )
679 $ THEN
680 info = -13
681 ELSE IF( ( igrade.EQ.1 .OR. igrade.EQ.3 .OR. igrade.EQ.4 .OR.
682 $ igrade.EQ.5 ) .AND. ( model.NE.-6 .AND. model.NE.0 .AND.
683 $ model.NE.6 ) .AND. condl.LT.one ) THEN
684 info = -14
685 ELSE IF( ( igrade.EQ.2 .OR. igrade.EQ.3 ) .AND.
686 $ ( moder.LT.-6 .OR. moder.GT.6 ) ) THEN
687 info = -16
688 ELSE IF( ( igrade.EQ.2 .OR. igrade.EQ.3 ) .AND.
689 $ ( moder.NE.-6 .AND. moder.NE.0 .AND. moder.NE.6 ) .AND.
690 $ condr.LT.one ) THEN
691 info = -17
692 ELSE IF( ipvtng.EQ.-1 .OR. ( ipvtng.EQ.3 .AND. m.NE.n ) .OR.
693 $ ( ( ipvtng.EQ.1 .OR. ipvtng.EQ.2 ) .AND. isym.EQ.0 ) )
694 $ THEN
695 info = -18
696 ELSE IF( ipvtng.NE.0 .AND. badpvt ) THEN
697 info = -19
698 ELSE IF( kl.LT.0 ) THEN
699 info = -20
700 ELSE IF( ku.LT.0 .OR. ( isym.EQ.0 .AND. kl.NE.ku ) ) THEN
701 info = -21
702 ELSE IF( sparse.LT.zero .OR. sparse.GT.one ) THEN
703 info = -22
704 ELSE IF( ipack.EQ.-1 .OR. ( ( ipack.EQ.1 .OR. ipack.EQ.2 .OR.
705 $ ipack.EQ.5 .OR. ipack.EQ.6 ) .AND. isym.EQ.1 ) .OR.
706 $ ( ipack.EQ.3 .AND. isym.EQ.1 .AND. ( kl.NE.0 .OR. m.NE.
707 $ n ) ) .OR. ( ipack.EQ.4 .AND. isym.EQ.1 .AND. ( ku.NE.
708 $ 0 .OR. m.NE.n ) ) ) THEN
709 info = -24
710 ELSE IF( ( ( ipack.EQ.0 .OR. ipack.EQ.1 .OR. ipack.EQ.2 ) .AND.
711 $ lda.LT.max( 1, m ) ) .OR. ( ( ipack.EQ.3 .OR. ipack.EQ.
712 $ 4 ) .AND. lda.LT.1 ) .OR. ( ( ipack.EQ.5 .OR. ipack.EQ.
713 $ 6 ) .AND. lda.LT.kuu+1 ) .OR.
714 $ ( ipack.EQ.7 .AND. lda.LT.kll+kuu+1 ) ) THEN
715 info = -26
716 END IF
717*
718 IF( info.NE.0 ) THEN
719 CALL xerbla( 'DLATMR', -info )
720 RETURN
721 END IF
722*
723* Decide if we can pivot consistently
724*
725 fulbnd = .false.
726 IF( kuu.EQ.n-1 .AND. kll.EQ.m-1 )
727 $ fulbnd = .true.
728*
729* Initialize random number generator
730*
731 DO 30 i = 1, 4
732 iseed( i ) = mod( abs( iseed( i ) ), 4096 )
733 30 CONTINUE
734*
735 iseed( 4 ) = 2*( iseed( 4 ) / 2 ) + 1
736*
737* 2) Set up D, DL, and DR, if indicated.
738*
739* Compute D according to COND and MODE
740*
741 CALL dlatm1( mode, cond, irsign, idist, iseed, d, mnmin, info )
742 IF( info.NE.0 ) THEN
743 info = 1
744 RETURN
745 END IF
746 IF( mode.NE.0 .AND. mode.NE.-6 .AND. mode.NE.6 ) THEN
747*
748* Scale by DMAX
749*
750 temp = abs( d( 1 ) )
751 DO 40 i = 2, mnmin
752 temp = max( temp, abs( d( i ) ) )
753 40 CONTINUE
754 IF( temp.EQ.zero .AND. dmax.NE.zero ) THEN
755 info = 2
756 RETURN
757 END IF
758 IF( temp.NE.zero ) THEN
759 alpha = dmax / temp
760 ELSE
761 alpha = one
762 END IF
763 DO 50 i = 1, mnmin
764 d( i ) = alpha*d( i )
765 50 CONTINUE
766*
767 END IF
768*
769* Compute DL if grading set
770*
771 IF( igrade.EQ.1 .OR. igrade.EQ.3 .OR. igrade.EQ.4 .OR. igrade.EQ.
772 $ 5 ) THEN
773 CALL dlatm1( model, condl, 0, idist, iseed, dl, m, info )
774 IF( info.NE.0 ) THEN
775 info = 3
776 RETURN
777 END IF
778 END IF
779*
780* Compute DR if grading set
781*
782 IF( igrade.EQ.2 .OR. igrade.EQ.3 ) THEN
783 CALL dlatm1( moder, condr, 0, idist, iseed, dr, n, info )
784 IF( info.NE.0 ) THEN
785 info = 4
786 RETURN
787 END IF
788 END IF
789*
790* 3) Generate IWORK if pivoting
791*
792 IF( ipvtng.GT.0 ) THEN
793 DO 60 i = 1, npvts
794 iwork( i ) = i
795 60 CONTINUE
796 IF( fulbnd ) THEN
797 DO 70 i = 1, npvts
798 k = ipivot( i )
799 j = iwork( i )
800 iwork( i ) = iwork( k )
801 iwork( k ) = j
802 70 CONTINUE
803 ELSE
804 DO 80 i = npvts, 1, -1
805 k = ipivot( i )
806 j = iwork( i )
807 iwork( i ) = iwork( k )
808 iwork( k ) = j
809 80 CONTINUE
810 END IF
811 END IF
812*
813* 4) Generate matrices for each kind of PACKing
814* Always sweep matrix columnwise (if symmetric, upper
815* half only) so that matrix generated does not depend
816* on PACK
817*
818 IF( fulbnd ) THEN
819*
820* Use DLATM3 so matrices generated with differing PIVOTing only
821* differ only in the order of their rows and/or columns.
822*
823 IF( ipack.EQ.0 ) THEN
824 IF( isym.EQ.0 ) THEN
825 DO 100 j = 1, n
826 DO 90 i = 1, j
827 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku,
828 $ idist, iseed, d, igrade, dl, dr, ipvtng,
829 $ iwork, sparse )
830 a( isub, jsub ) = temp
831 a( jsub, isub ) = temp
832 90 CONTINUE
833 100 CONTINUE
834 ELSE IF( isym.EQ.1 ) THEN
835 DO 120 j = 1, n
836 DO 110 i = 1, m
837 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku,
838 $ idist, iseed, d, igrade, dl, dr, ipvtng,
839 $ iwork, sparse )
840 a( isub, jsub ) = temp
841 110 CONTINUE
842 120 CONTINUE
843 END IF
844*
845 ELSE IF( ipack.EQ.1 ) THEN
846*
847 DO 140 j = 1, n
848 DO 130 i = 1, j
849 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku, idist,
850 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
851 $ sparse )
852 mnsub = min( isub, jsub )
853 mxsub = max( isub, jsub )
854 a( mnsub, mxsub ) = temp
855 IF( mnsub.NE.mxsub )
856 $ a( mxsub, mnsub ) = zero
857 130 CONTINUE
858 140 CONTINUE
859*
860 ELSE IF( ipack.EQ.2 ) THEN
861*
862 DO 160 j = 1, n
863 DO 150 i = 1, j
864 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku, idist,
865 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
866 $ sparse )
867 mnsub = min( isub, jsub )
868 mxsub = max( isub, jsub )
869 a( mxsub, mnsub ) = temp
870 IF( mnsub.NE.mxsub )
871 $ a( mnsub, mxsub ) = zero
872 150 CONTINUE
873 160 CONTINUE
874*
875 ELSE IF( ipack.EQ.3 ) THEN
876*
877 DO 180 j = 1, n
878 DO 170 i = 1, j
879 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku, idist,
880 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
881 $ sparse )
882*
883* Compute K = location of (ISUB,JSUB) entry in packed
884* array
885*
886 mnsub = min( isub, jsub )
887 mxsub = max( isub, jsub )
888 k = mxsub*( mxsub-1 ) / 2 + mnsub
889*
890* Convert K to (IISUB,JJSUB) location
891*
892 jjsub = ( k-1 ) / lda + 1
893 iisub = k - lda*( jjsub-1 )
894*
895 a( iisub, jjsub ) = temp
896 170 CONTINUE
897 180 CONTINUE
898*
899 ELSE IF( ipack.EQ.4 ) THEN
900*
901 DO 200 j = 1, n
902 DO 190 i = 1, j
903 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku, idist,
904 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
905 $ sparse )
906*
907* Compute K = location of (I,J) entry in packed array
908*
909 mnsub = min( isub, jsub )
910 mxsub = max( isub, jsub )
911 IF( mnsub.EQ.1 ) THEN
912 k = mxsub
913 ELSE
914 k = n*( n+1 ) / 2 - ( n-mnsub+1 )*( n-mnsub+2 ) /
915 $ 2 + mxsub - mnsub + 1
916 END IF
917*
918* Convert K to (IISUB,JJSUB) location
919*
920 jjsub = ( k-1 ) / lda + 1
921 iisub = k - lda*( jjsub-1 )
922*
923 a( iisub, jjsub ) = temp
924 190 CONTINUE
925 200 CONTINUE
926*
927 ELSE IF( ipack.EQ.5 ) THEN
928*
929 DO 220 j = 1, n
930 DO 210 i = j - kuu, j
931 IF( i.LT.1 ) THEN
932 a( j-i+1, i+n ) = zero
933 ELSE
934 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku,
935 $ idist, iseed, d, igrade, dl, dr, ipvtng,
936 $ iwork, sparse )
937 mnsub = min( isub, jsub )
938 mxsub = max( isub, jsub )
939 a( mxsub-mnsub+1, mnsub ) = temp
940 END IF
941 210 CONTINUE
942 220 CONTINUE
943*
944 ELSE IF( ipack.EQ.6 ) THEN
945*
946 DO 240 j = 1, n
947 DO 230 i = j - kuu, j
948 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku, idist,
949 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
950 $ sparse )
951 mnsub = min( isub, jsub )
952 mxsub = max( isub, jsub )
953 a( mnsub-mxsub+kuu+1, mxsub ) = temp
954 230 CONTINUE
955 240 CONTINUE
956*
957 ELSE IF( ipack.EQ.7 ) THEN
958*
959 IF( isym.EQ.0 ) THEN
960 DO 260 j = 1, n
961 DO 250 i = j - kuu, j
962 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku,
963 $ idist, iseed, d, igrade, dl, dr, ipvtng,
964 $ iwork, sparse )
965 mnsub = min( isub, jsub )
966 mxsub = max( isub, jsub )
967 a( mnsub-mxsub+kuu+1, mxsub ) = temp
968 IF( i.LT.1 )
969 $ a( j-i+1+kuu, i+n ) = zero
970 IF( i.GE.1 .AND. mnsub.NE.mxsub )
971 $ a( mxsub-mnsub+1+kuu, mnsub ) = temp
972 250 CONTINUE
973 260 CONTINUE
974 ELSE IF( isym.EQ.1 ) THEN
975 DO 280 j = 1, n
976 DO 270 i = j - kuu, j + kll
977 temp = dlatm3( m, n, i, j, isub, jsub, kl, ku,
978 $ idist, iseed, d, igrade, dl, dr, ipvtng,
979 $ iwork, sparse )
980 a( isub-jsub+kuu+1, jsub ) = temp
981 270 CONTINUE
982 280 CONTINUE
983 END IF
984*
985 END IF
986*
987 ELSE
988*
989* Use DLATM2
990*
991 IF( ipack.EQ.0 ) THEN
992 IF( isym.EQ.0 ) THEN
993 DO 300 j = 1, n
994 DO 290 i = 1, j
995 a( i, j ) = dlatm2( m, n, i, j, kl, ku, idist,
996 $ iseed, d, igrade, dl, dr, ipvtng,
997 $ iwork, sparse )
998 a( j, i ) = a( i, j )
999 290 CONTINUE
1000 300 CONTINUE
1001 ELSE IF( isym.EQ.1 ) THEN
1002 DO 320 j = 1, n
1003 DO 310 i = 1, m
1004 a( i, j ) = dlatm2( m, n, i, j, kl, ku, idist,
1005 $ iseed, d, igrade, dl, dr, ipvtng,
1006 $ iwork, sparse )
1007 310 CONTINUE
1008 320 CONTINUE
1009 END IF
1010*
1011 ELSE IF( ipack.EQ.1 ) THEN
1012*
1013 DO 340 j = 1, n
1014 DO 330 i = 1, j
1015 a( i, j ) = dlatm2( m, n, i, j, kl, ku, idist, iseed,
1016 $ d, igrade, dl, dr, ipvtng, iwork, sparse )
1017 IF( i.NE.j )
1018 $ a( j, i ) = zero
1019 330 CONTINUE
1020 340 CONTINUE
1021*
1022 ELSE IF( ipack.EQ.2 ) THEN
1023*
1024 DO 360 j = 1, n
1025 DO 350 i = 1, j
1026 a( j, i ) = dlatm2( m, n, i, j, kl, ku, idist, iseed,
1027 $ d, igrade, dl, dr, ipvtng, iwork, sparse )
1028 IF( i.NE.j )
1029 $ a( i, j ) = zero
1030 350 CONTINUE
1031 360 CONTINUE
1032*
1033 ELSE IF( ipack.EQ.3 ) THEN
1034*
1035 isub = 0
1036 jsub = 1
1037 DO 380 j = 1, n
1038 DO 370 i = 1, j
1039 isub = isub + 1
1040 IF( isub.GT.lda ) THEN
1041 isub = 1
1042 jsub = jsub + 1
1043 END IF
1044 a( isub, jsub ) = dlatm2( m, n, i, j, kl, ku, idist,
1045 $ iseed, d, igrade, dl, dr, ipvtng,
1046 $ iwork, sparse )
1047 370 CONTINUE
1048 380 CONTINUE
1049*
1050 ELSE IF( ipack.EQ.4 ) THEN
1051*
1052 IF( isym.EQ.0 ) THEN
1053 DO 400 j = 1, n
1054 DO 390 i = 1, j
1055*
1056* Compute K = location of (I,J) entry in packed array
1057*
1058 IF( i.EQ.1 ) THEN
1059 k = j
1060 ELSE
1061 k = n*( n+1 ) / 2 - ( n-i+1 )*( n-i+2 ) / 2 +
1062 $ j - i + 1
1063 END IF
1064*
1065* Convert K to (ISUB,JSUB) location
1066*
1067 jsub = ( k-1 ) / lda + 1
1068 isub = k - lda*( jsub-1 )
1069*
1070 a( isub, jsub ) = dlatm2( m, n, i, j, kl, ku,
1071 $ idist, iseed, d, igrade, dl, dr,
1072 $ ipvtng, iwork, sparse )
1073 390 CONTINUE
1074 400 CONTINUE
1075 ELSE
1076 isub = 0
1077 jsub = 1
1078 DO 420 j = 1, n
1079 DO 410 i = j, m
1080 isub = isub + 1
1081 IF( isub.GT.lda ) THEN
1082 isub = 1
1083 jsub = jsub + 1
1084 END IF
1085 a( isub, jsub ) = dlatm2( m, n, i, j, kl, ku,
1086 $ idist, iseed, d, igrade, dl, dr,
1087 $ ipvtng, iwork, sparse )
1088 410 CONTINUE
1089 420 CONTINUE
1090 END IF
1091*
1092 ELSE IF( ipack.EQ.5 ) THEN
1093*
1094 DO 440 j = 1, n
1095 DO 430 i = j - kuu, j
1096 IF( i.LT.1 ) THEN
1097 a( j-i+1, i+n ) = zero
1098 ELSE
1099 a( j-i+1, i ) = dlatm2( m, n, i, j, kl, ku, idist,
1100 $ iseed, d, igrade, dl, dr, ipvtng,
1101 $ iwork, sparse )
1102 END IF
1103 430 CONTINUE
1104 440 CONTINUE
1105*
1106 ELSE IF( ipack.EQ.6 ) THEN
1107*
1108 DO 460 j = 1, n
1109 DO 450 i = j - kuu, j
1110 a( i-j+kuu+1, j ) = dlatm2( m, n, i, j, kl, ku, idist,
1111 $ iseed, d, igrade, dl, dr, ipvtng,
1112 $ iwork, sparse )
1113 450 CONTINUE
1114 460 CONTINUE
1115*
1116 ELSE IF( ipack.EQ.7 ) THEN
1117*
1118 IF( isym.EQ.0 ) THEN
1119 DO 480 j = 1, n
1120 DO 470 i = j - kuu, j
1121 a( i-j+kuu+1, j ) = dlatm2( m, n, i, j, kl, ku,
1122 $ idist, iseed, d, igrade, dl,
1123 $ dr, ipvtng, iwork, sparse )
1124 IF( i.LT.1 )
1125 $ a( j-i+1+kuu, i+n ) = zero
1126 IF( i.GE.1 .AND. i.NE.j )
1127 $ a( j-i+1+kuu, i ) = a( i-j+kuu+1, j )
1128 470 CONTINUE
1129 480 CONTINUE
1130 ELSE IF( isym.EQ.1 ) THEN
1131 DO 500 j = 1, n
1132 DO 490 i = j - kuu, j + kll
1133 a( i-j+kuu+1, j ) = dlatm2( m, n, i, j, kl, ku,
1134 $ idist, iseed, d, igrade, dl,
1135 $ dr, ipvtng, iwork, sparse )
1136 490 CONTINUE
1137 500 CONTINUE
1138 END IF
1139*
1140 END IF
1141*
1142 END IF
1143*
1144* 5) Scaling the norm
1145*
1146 IF( ipack.EQ.0 ) THEN
1147 onorm = dlange( 'M', m, n, a, lda, tempa )
1148 ELSE IF( ipack.EQ.1 ) THEN
1149 onorm = dlansy( 'M', 'U', n, a, lda, tempa )
1150 ELSE IF( ipack.EQ.2 ) THEN
1151 onorm = dlansy( 'M', 'L', n, a, lda, tempa )
1152 ELSE IF( ipack.EQ.3 ) THEN
1153 onorm = dlansp( 'M', 'U', n, a, tempa )
1154 ELSE IF( ipack.EQ.4 ) THEN
1155 onorm = dlansp( 'M', 'L', n, a, tempa )
1156 ELSE IF( ipack.EQ.5 ) THEN
1157 onorm = dlansb( 'M', 'L', n, kll, a, lda, tempa )
1158 ELSE IF( ipack.EQ.6 ) THEN
1159 onorm = dlansb( 'M', 'U', n, kuu, a, lda, tempa )
1160 ELSE IF( ipack.EQ.7 ) THEN
1161 onorm = dlangb( 'M', n, kll, kuu, a, lda, tempa )
1162 END IF
1163*
1164 IF( anorm.GE.zero ) THEN
1165*
1166 IF( anorm.GT.zero .AND. onorm.EQ.zero ) THEN
1167*
1168* Desired scaling impossible
1169*
1170 info = 5
1171 RETURN
1172*
1173 ELSE IF( ( anorm.GT.one .AND. onorm.LT.one ) .OR.
1174 $ ( anorm.LT.one .AND. onorm.GT.one ) ) THEN
1175*
1176* Scale carefully to avoid over / underflow
1177*
1178 IF( ipack.LE.2 ) THEN
1179 DO 510 j = 1, n
1180 CALL dscal( m, one / onorm, a( 1, j ), 1 )
1181 CALL dscal( m, anorm, a( 1, j ), 1 )
1182 510 CONTINUE
1183*
1184 ELSE IF( ipack.EQ.3 .OR. ipack.EQ.4 ) THEN
1185*
1186 CALL dscal( n*( n+1 ) / 2, one / onorm, a, 1 )
1187 CALL dscal( n*( n+1 ) / 2, anorm, a, 1 )
1188*
1189 ELSE IF( ipack.GE.5 ) THEN
1190*
1191 DO 520 j = 1, n
1192 CALL dscal( kll+kuu+1, one / onorm, a( 1, j ), 1 )
1193 CALL dscal( kll+kuu+1, anorm, a( 1, j ), 1 )
1194 520 CONTINUE
1195*
1196 END IF
1197*
1198 ELSE
1199*
1200* Scale straightforwardly
1201*
1202 IF( ipack.LE.2 ) THEN
1203 DO 530 j = 1, n
1204 CALL dscal( m, anorm / onorm, a( 1, j ), 1 )
1205 530 CONTINUE
1206*
1207 ELSE IF( ipack.EQ.3 .OR. ipack.EQ.4 ) THEN
1208*
1209 CALL dscal( n*( n+1 ) / 2, anorm / onorm, a, 1 )
1210*
1211 ELSE IF( ipack.GE.5 ) THEN
1212*
1213 DO 540 j = 1, n
1214 CALL dscal( kll+kuu+1, anorm / onorm, a( 1, j ), 1 )
1215 540 CONTINUE
1216 END IF
1217*
1218 END IF
1219*
1220 END IF
1221*
1222* End of DLATMR
1223*
1224 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dlatm1(mode, cond, irsign, idist, iseed, d, n, info)
DLATM1
Definition dlatm1.f:135
subroutine dlatmr(m, n, dist, iseed, sym, d, mode, cond, dmax, rsign, grade, dl, model, condl, dr, moder, condr, pivtng, ipivot, kl, ku, sparse, anorm, pack, a, lda, iwork, info)
DLATMR
Definition dlatmr.f:471
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79