LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ claror()

subroutine claror ( character side,
character init,
integer m,
integer n,
complex, dimension( lda, * ) a,
integer lda,
integer, dimension( 4 ) iseed,
complex, dimension( * ) x,
integer info )

CLAROR

Purpose:
!>
!>    CLAROR pre- or post-multiplies an M by N matrix A by a random
!>    unitary matrix U, overwriting A. A may optionally be
!>    initialized to the identity matrix before multiplying by U.
!>    U is generated using the method of G.W. Stewart
!>    ( SIAM J. Numer. Anal. 17, 1980, pp. 403-409 ).
!>    (BLAS-2 version)
!>
Parameters
[in]SIDE
!>          SIDE is CHARACTER*1
!>           SIDE specifies whether A is multiplied on the left or right
!>           by U.
!>       SIDE = 'L'   Multiply A on the left (premultiply) by U
!>       SIDE = 'R'   Multiply A on the right (postmultiply) by UC>       SIDE = 'C'   Multiply A on the left by U and the right by UC>       SIDE = 'T'   Multiply A on the left by U and the right by U'
!>           Not modified.
!>
[in]INIT
!>          INIT is CHARACTER*1
!>           INIT specifies whether or not A should be initialized to
!>           the identity matrix.
!>              INIT = 'I'   Initialize A to (a section of) the
!>                           identity matrix before applying U.
!>              INIT = 'N'   No initialization.  Apply U to the
!>                           input matrix A.
!>
!>           INIT = 'I' may be used to generate square (i.e., unitary)
!>           or rectangular orthogonal matrices (orthogonality being
!>           in the sense of CDOTC):
!>
!>           For square matrices, M=N, and SIDE many be either 'L' or
!>           'R'; the rows will be orthogonal to each other, as will the
!>           columns.
!>           For rectangular matrices where M < N, SIDE = 'R' will
!>           produce a dense matrix whose rows will be orthogonal and
!>           whose columns will not, while SIDE = 'L' will produce a
!>           matrix whose rows will be orthogonal, and whose first M
!>           columns will be orthogonal, the remaining columns being
!>           zero.
!>           For matrices where M > N, just use the previous
!>           explanation, interchanging 'L' and 'R' and  and
!>           .
!>
!>           Not modified.
!>
[in]M
!>          M is INTEGER
!>           Number of rows of A. Not modified.
!>
[in]N
!>          N is INTEGER
!>           Number of columns of A. Not modified.
!>
[in,out]A
!>          A is COMPLEX array, dimension ( LDA, N )
!>           Input and output array. Overwritten by U A ( if SIDE = 'L' )
!>           or by A U ( if SIDE = 'R' )
!>           or by U A U* ( if SIDE = 'C')
!>           or by U A U' ( if SIDE = 'T') on exit.
!>
[in]LDA
!>          LDA is INTEGER
!>           Leading dimension of A. Must be at least MAX ( 1, M ).
!>           Not modified.
!>
[in,out]ISEED
!>          ISEED is INTEGER array, dimension ( 4 )
!>           On entry ISEED specifies the seed of the random number
!>           generator. The array elements should be between 0 and 4095;
!>           if not they will be reduced mod 4096.  Also, ISEED(4) must
!>           be odd.  The random number generator uses a linear
!>           congruential sequence limited to small integers, and so
!>           should produce machine independent random numbers. The
!>           values of ISEED are changed on exit, and can be used in the
!>           next call to CLAROR to continue the same random number
!>           sequence.
!>           Modified.
!>
[out]X
!>          X is COMPLEX array, dimension ( 3*MAX( M, N ) )
!>           Workspace. Of length:
!>               2*M + N if SIDE = 'L',
!>               2*N + M if SIDE = 'R',
!>               3*N     if SIDE = 'C' or 'T'.
!>           Modified.
!>
[out]INFO
!>          INFO is INTEGER
!>           An error flag.  It is set to:
!>            0  if no error.
!>            1  if CLARND returned a bad random number (installation
!>               problem)
!>           -1  if SIDE is not L, R, C, or T.
!>           -3  if M is negative.
!>           -4  if N is negative or if SIDE is C or T and N is not equal
!>               to M.
!>           -6  if LDA is less than M.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 157 of file claror.f.

158*
159* -- LAPACK auxiliary routine --
160* -- LAPACK is a software package provided by Univ. of Tennessee, --
161* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
162*
163* .. Scalar Arguments ..
164 CHARACTER INIT, SIDE
165 INTEGER INFO, LDA, M, N
166* ..
167* .. Array Arguments ..
168 INTEGER ISEED( 4 )
169 COMPLEX A( LDA, * ), X( * )
170* ..
171*
172* =====================================================================
173*
174* .. Parameters ..
175 REAL ZERO, ONE, TOOSML
176 parameter( zero = 0.0e+0, one = 1.0e+0,
177 $ toosml = 1.0e-20 )
178 COMPLEX CZERO, CONE
179 parameter( czero = ( 0.0e+0, 0.0e+0 ),
180 $ cone = ( 1.0e+0, 0.0e+0 ) )
181* ..
182* .. Local Scalars ..
183 INTEGER IROW, ITYPE, IXFRM, J, JCOL, KBEG, NXFRM
184 REAL FACTOR, XABS, XNORM
185 COMPLEX CSIGN, XNORMS
186* ..
187* .. External Functions ..
188 LOGICAL LSAME
189 REAL SCNRM2
190 COMPLEX CLARND
191 EXTERNAL lsame, scnrm2, clarnd
192* ..
193* .. External Subroutines ..
194 EXTERNAL cgemv, cgerc, clacgv, claset, cscal, xerbla
195* ..
196* .. Intrinsic Functions ..
197 INTRINSIC abs, cmplx, conjg
198* ..
199* .. Executable Statements ..
200*
201 info = 0
202 IF( n.EQ.0 .OR. m.EQ.0 )
203 $ RETURN
204*
205 itype = 0
206 IF( lsame( side, 'L' ) ) THEN
207 itype = 1
208 ELSE IF( lsame( side, 'R' ) ) THEN
209 itype = 2
210 ELSE IF( lsame( side, 'C' ) ) THEN
211 itype = 3
212 ELSE IF( lsame( side, 'T' ) ) THEN
213 itype = 4
214 END IF
215*
216* Check for argument errors.
217*
218 IF( itype.EQ.0 ) THEN
219 info = -1
220 ELSE IF( m.LT.0 ) THEN
221 info = -3
222 ELSE IF( n.LT.0 .OR. ( itype.EQ.3 .AND. n.NE.m ) ) THEN
223 info = -4
224 ELSE IF( lda.LT.m ) THEN
225 info = -6
226 END IF
227 IF( info.NE.0 ) THEN
228 CALL xerbla( 'CLAROR', -info )
229 RETURN
230 END IF
231*
232 IF( itype.EQ.1 ) THEN
233 nxfrm = m
234 ELSE
235 nxfrm = n
236 END IF
237*
238* Initialize A to the identity matrix if desired
239*
240 IF( lsame( init, 'I' ) )
241 $ CALL claset( 'Full', m, n, czero, cone, a, lda )
242*
243* If no rotation possible, still multiply by
244* a random complex number from the circle |x| = 1
245*
246* 2) Compute Rotation by computing Householder
247* Transformations H(2), H(3), ..., H(n). Note that the
248* order in which they are computed is irrelevant.
249*
250 DO 40 j = 1, nxfrm
251 x( j ) = czero
252 40 CONTINUE
253*
254 DO 60 ixfrm = 2, nxfrm
255 kbeg = nxfrm - ixfrm + 1
256*
257* Generate independent normal( 0, 1 ) random numbers
258*
259 DO 50 j = kbeg, nxfrm
260 x( j ) = clarnd( 3, iseed )
261 50 CONTINUE
262*
263* Generate a Householder transformation from the random vector X
264*
265 xnorm = scnrm2( ixfrm, x( kbeg ), 1 )
266 xabs = abs( x( kbeg ) )
267 IF( xabs.NE.czero ) THEN
268 csign = x( kbeg ) / xabs
269 ELSE
270 csign = cone
271 END IF
272 xnorms = csign*xnorm
273 x( nxfrm+kbeg ) = -csign
274 factor = xnorm*( xnorm+xabs )
275 IF( abs( factor ).LT.toosml ) THEN
276 info = 1
277 CALL xerbla( 'CLAROR', -info )
278 RETURN
279 ELSE
280 factor = one / factor
281 END IF
282 x( kbeg ) = x( kbeg ) + xnorms
283*
284* Apply Householder transformation to A
285*
286 IF( itype.EQ.1 .OR. itype.EQ.3 .OR. itype.EQ.4 ) THEN
287*
288* Apply H(k) on the left of A
289*
290 CALL cgemv( 'C', ixfrm, n, cone, a( kbeg, 1 ), lda,
291 $ x( kbeg ), 1, czero, x( 2*nxfrm+1 ), 1 )
292 CALL cgerc( ixfrm, n, -cmplx( factor ), x( kbeg ), 1,
293 $ x( 2*nxfrm+1 ), 1, a( kbeg, 1 ), lda )
294*
295 END IF
296*
297 IF( itype.GE.2 .AND. itype.LE.4 ) THEN
298*
299* Apply H(k)* (or H(k)') on the right of A
300*
301 IF( itype.EQ.4 ) THEN
302 CALL clacgv( ixfrm, x( kbeg ), 1 )
303 END IF
304*
305 CALL cgemv( 'N', m, ixfrm, cone, a( 1, kbeg ), lda,
306 $ x( kbeg ), 1, czero, x( 2*nxfrm+1 ), 1 )
307 CALL cgerc( m, ixfrm, -cmplx( factor ), x( 2*nxfrm+1 ), 1,
308 $ x( kbeg ), 1, a( 1, kbeg ), lda )
309*
310 END IF
311 60 CONTINUE
312*
313 x( 1 ) = clarnd( 3, iseed )
314 xabs = abs( x( 1 ) )
315 IF( xabs.NE.zero ) THEN
316 csign = x( 1 ) / xabs
317 ELSE
318 csign = cone
319 END IF
320 x( 2*nxfrm ) = csign
321*
322* Scale the matrix A by D.
323*
324 IF( itype.EQ.1 .OR. itype.EQ.3 .OR. itype.EQ.4 ) THEN
325 DO 70 irow = 1, m
326 CALL cscal( n, conjg( x( nxfrm+irow ) ), a( irow, 1 ), lda )
327 70 CONTINUE
328 END IF
329*
330 IF( itype.EQ.2 .OR. itype.EQ.3 ) THEN
331 DO 80 jcol = 1, n
332 CALL cscal( m, x( nxfrm+jcol ), a( 1, jcol ), 1 )
333 80 CONTINUE
334 END IF
335*
336 IF( itype.EQ.4 ) THEN
337 DO 90 jcol = 1, n
338 CALL cscal( m, conjg( x( nxfrm+jcol ) ), a( 1, jcol ), 1 )
339 90 CONTINUE
340 END IF
341 RETURN
342*
343* End of CLAROR
344*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
clarnd
complex function clarnd(idist, iseed)
CLARND
Definition clarnd.f:75
cgemv
subroutine cgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
CGEMV
Definition cgemv.f:160
cgerc
subroutine cgerc(m, n, alpha, x, incx, y, incy, a, lda)
CGERC
Definition cgerc.f:130
clacgv
subroutine clacgv(n, x, incx)
CLACGV conjugates a complex vector.
Definition clacgv.f:72
claset
subroutine claset(uplo, m, n, alpha, beta, a, lda)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition claset.f:104
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
scnrm2
real(wp) function scnrm2(n, x, incx)
SCNRM2
Definition scnrm2.f90:90
cscal
subroutine cscal(n, ca, cx, incx)
CSCAL
Definition cscal.f:78
Here is the call graph for this function:
Here is the caller graph for this function: