LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlaror.f
Go to the documentation of this file.
1*> \brief \b DLAROR
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 DLAROR( SIDE, INIT, M, N, A, LDA, ISEED, X, INFO )
12*
13* .. Scalar Arguments ..
14* CHARACTER INIT, SIDE
15* INTEGER INFO, LDA, M, N
16* ..
17* .. Array Arguments ..
18* INTEGER ISEED( 4 )
19* DOUBLE PRECISION A( LDA, * ), X( * )
20* ..
21*
22*
23*> \par Purpose:
24* =============
25*>
26*> \verbatim
27*>
28*> DLAROR pre- or post-multiplies an M by N matrix A by a random
29*> orthogonal matrix U, overwriting A. A may optionally be initialized
30*> to the identity matrix before multiplying by U. U is generated using
31*> the method of G.W. Stewart (SIAM J. Numer. Anal. 17, 1980, 403-409).
32*> \endverbatim
33*
34* Arguments:
35* ==========
36*
37*> \param[in] SIDE
38*> \verbatim
39*> SIDE is CHARACTER*1
40*> Specifies whether A is multiplied on the left or right by U.
41*> = 'L': Multiply A on the left (premultiply) by U
42*> = 'R': Multiply A on the right (postmultiply) by U'
43*> = 'C' or 'T': Multiply A on the left by U and the right
44*> by U' (Here, U' means U-transpose.)
45*> \endverbatim
46*>
47*> \param[in] INIT
48*> \verbatim
49*> INIT is CHARACTER*1
50*> Specifies whether or not A should be initialized to the
51*> identity matrix.
52*> = 'I': Initialize A to (a section of) the identity matrix
53*> before applying U.
54*> = 'N': No initialization. Apply U to the input matrix A.
55*>
56*> INIT = 'I' may be used to generate square or rectangular
57*> orthogonal matrices:
58*>
59*> For M = N and SIDE = 'L' or 'R', the rows will be orthogonal
60*> to each other, as will the columns.
61*>
62*> If M < N, SIDE = 'R' produces a dense matrix whose rows are
63*> orthogonal and whose columns are not, while SIDE = 'L'
64*> produces a matrix whose rows are orthogonal, and whose first
65*> M columns are orthogonal, and whose remaining columns are
66*> zero.
67*>
68*> If M > N, SIDE = 'L' produces a dense matrix whose columns
69*> are orthogonal and whose rows are not, while SIDE = 'R'
70*> produces a matrix whose columns are orthogonal, and whose
71*> first M rows are orthogonal, and whose remaining rows are
72*> zero.
73*> \endverbatim
74*>
75*> \param[in] M
76*> \verbatim
77*> M is INTEGER
78*> The number of rows of A.
79*> \endverbatim
80*>
81*> \param[in] N
82*> \verbatim
83*> N is INTEGER
84*> The number of columns of A.
85*> \endverbatim
86*>
87*> \param[in,out] A
88*> \verbatim
89*> A is DOUBLE PRECISION array, dimension (LDA, N)
90*> On entry, the array A.
91*> On exit, overwritten by U A ( if SIDE = 'L' ),
92*> or by A U ( if SIDE = 'R' ),
93*> or by U A U' ( if SIDE = 'C' or 'T').
94*> \endverbatim
95*>
96*> \param[in] LDA
97*> \verbatim
98*> LDA is INTEGER
99*> The leading dimension of the array A. LDA >= max(1,M).
100*> \endverbatim
101*>
102*> \param[in,out] ISEED
103*> \verbatim
104*> ISEED is INTEGER array, dimension (4)
105*> On entry ISEED specifies the seed of the random number
106*> generator. The array elements should be between 0 and 4095;
107*> if not they will be reduced mod 4096. Also, ISEED(4) must
108*> be odd. The random number generator uses a linear
109*> congruential sequence limited to small integers, and so
110*> should produce machine independent random numbers. The
111*> values of ISEED are changed on exit, and can be used in the
112*> next call to DLAROR to continue the same random number
113*> sequence.
114*> \endverbatim
115*>
116*> \param[out] X
117*> \verbatim
118*> X is DOUBLE PRECISION array, dimension (3*MAX( M, N ))
119*> Workspace of length
120*> 2*M + N if SIDE = 'L',
121*> 2*N + M if SIDE = 'R',
122*> 3*N if SIDE = 'C' or 'T'.
123*> \endverbatim
124*>
125*> \param[out] INFO
126*> \verbatim
127*> INFO is INTEGER
128*> An error flag. It is set to:
129*> = 0: normal return
130*> < 0: if INFO = -k, the k-th argument had an illegal value
131*> = 1: if the random numbers generated by DLARND are bad.
132*> \endverbatim
133*
134* Authors:
135* ========
136*
137*> \author Univ. of Tennessee
138*> \author Univ. of California Berkeley
139*> \author Univ. of Colorado Denver
140*> \author NAG Ltd.
141*
142*> \ingroup double_matgen
143*
144* =====================================================================
145 SUBROUTINE dlaror( SIDE, INIT, M, N, A, LDA, ISEED, X, INFO )
146*
147* -- LAPACK auxiliary routine --
148* -- LAPACK is a software package provided by Univ. of Tennessee, --
149* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
150*
151* .. Scalar Arguments ..
152 CHARACTER INIT, SIDE
153 INTEGER INFO, LDA, M, N
154* ..
155* .. Array Arguments ..
156 INTEGER ISEED( 4 )
157 DOUBLE PRECISION A( LDA, * ), X( * )
158* ..
159*
160* =====================================================================
161*
162* .. Parameters ..
163 DOUBLE PRECISION ZERO, ONE, TOOSML
164 parameter( zero = 0.0d+0, one = 1.0d+0,
165 $ toosml = 1.0d-20 )
166* ..
167* .. Local Scalars ..
168 INTEGER IROW, ITYPE, IXFRM, J, JCOL, KBEG, NXFRM
169 DOUBLE PRECISION FACTOR, XNORM, XNORMS
170* ..
171* .. External Functions ..
172 LOGICAL LSAME
173 DOUBLE PRECISION DLARND, DNRM2
174 EXTERNAL lsame, dlarnd, dnrm2
175* ..
176* .. External Subroutines ..
177 EXTERNAL dgemv, dger, dlaset, dscal, xerbla
178* ..
179* .. Intrinsic Functions ..
180 INTRINSIC abs, sign
181* ..
182* .. Executable Statements ..
183*
184 info = 0
185 IF( n.EQ.0 .OR. m.EQ.0 )
186 $ RETURN
187*
188 itype = 0
189 IF( lsame( side, 'L' ) ) THEN
190 itype = 1
191 ELSE IF( lsame( side, 'R' ) ) THEN
192 itype = 2
193 ELSE IF( lsame( side, 'C' ) .OR. lsame( side, 'T' ) ) THEN
194 itype = 3
195 END IF
196*
197* Check for argument errors.
198*
199 IF( itype.EQ.0 ) THEN
200 info = -1
201 ELSE IF( m.LT.0 ) THEN
202 info = -3
203 ELSE IF( n.LT.0 .OR. ( itype.EQ.3 .AND. n.NE.m ) ) THEN
204 info = -4
205 ELSE IF( lda.LT.m ) THEN
206 info = -6
207 END IF
208 IF( info.NE.0 ) THEN
209 CALL xerbla( 'DLAROR', -info )
210 RETURN
211 END IF
212*
213 IF( itype.EQ.1 ) THEN
214 nxfrm = m
215 ELSE
216 nxfrm = n
217 END IF
218*
219* Initialize A to the identity matrix if desired
220*
221 IF( lsame( init, 'I' ) )
222 $ CALL dlaset( 'Full', m, n, zero, one, a, lda )
223*
224* If no rotation possible, multiply by random +/-1
225*
226* Compute rotation by computing Householder transformations
227* H(2), H(3), ..., H(nhouse)
228*
229 DO 10 j = 1, nxfrm
230 x( j ) = zero
231 10 CONTINUE
232*
233 DO 30 ixfrm = 2, nxfrm
234 kbeg = nxfrm - ixfrm + 1
235*
236* Generate independent normal( 0, 1 ) random numbers
237*
238 DO 20 j = kbeg, nxfrm
239 x( j ) = dlarnd( 3, iseed )
240 20 CONTINUE
241*
242* Generate a Householder transformation from the random vector X
243*
244 xnorm = dnrm2( ixfrm, x( kbeg ), 1 )
245 xnorms = sign( xnorm, x( kbeg ) )
246 x( kbeg+nxfrm ) = sign( one, -x( kbeg ) )
247 factor = xnorms*( xnorms+x( kbeg ) )
248 IF( abs( factor ).LT.toosml ) THEN
249 info = 1
250 CALL xerbla( 'DLAROR', info )
251 RETURN
252 ELSE
253 factor = one / factor
254 END IF
255 x( kbeg ) = x( kbeg ) + xnorms
256*
257* Apply Householder transformation to A
258*
259 IF( itype.EQ.1 .OR. itype.EQ.3 ) THEN
260*
261* Apply H(k) from the left.
262*
263 CALL dgemv( 'T', ixfrm, n, one, a( kbeg, 1 ), lda,
264 $ x( kbeg ), 1, zero, x( 2*nxfrm+1 ), 1 )
265 CALL dger( ixfrm, n, -factor, x( kbeg ), 1, x( 2*nxfrm+1 ),
266 $ 1, a( kbeg, 1 ), lda )
267*
268 END IF
269*
270 IF( itype.EQ.2 .OR. itype.EQ.3 ) THEN
271*
272* Apply H(k) from the right.
273*
274 CALL dgemv( 'N', m, ixfrm, one, a( 1, kbeg ), lda,
275 $ x( kbeg ), 1, zero, x( 2*nxfrm+1 ), 1 )
276 CALL dger( m, ixfrm, -factor, x( 2*nxfrm+1 ), 1, x( kbeg ),
277 $ 1, a( 1, kbeg ), lda )
278*
279 END IF
280 30 CONTINUE
281*
282 x( 2*nxfrm ) = sign( one, dlarnd( 3, iseed ) )
283*
284* Scale the matrix A by D.
285*
286 IF( itype.EQ.1 .OR. itype.EQ.3 ) THEN
287 DO 40 irow = 1, m
288 CALL dscal( n, x( nxfrm+irow ), a( irow, 1 ), lda )
289 40 CONTINUE
290 END IF
291*
292 IF( itype.EQ.2 .OR. itype.EQ.3 ) THEN
293 DO 50 jcol = 1, n
294 CALL dscal( m, x( nxfrm+jcol ), a( 1, jcol ), 1 )
295 50 CONTINUE
296 END IF
297 RETURN
298*
299* End of DLAROR
300*
301 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dlaror(side, init, m, n, a, lda, iseed, x, info)
DLAROR
Definition dlaror.f:146
subroutine dgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
DGEMV
Definition dgemv.f:158
subroutine dger(m, n, alpha, x, incx, y, incy, a, lda)
DGER
Definition dger.f:130
subroutine dlaset(uplo, m, n, alpha, beta, a, lda)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition dlaset.f:110
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79