LAPACK 3.12.0
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 "rows" and
           "columns".

           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*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
complex function clarnd(idist, iseed)
CLARND
Definition clarnd.f:75
subroutine cgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
CGEMV
Definition cgemv.f:160
subroutine cgerc(m, n, alpha, x, incx, y, incy, a, lda)
CGERC
Definition cgerc.f:130
subroutine clacgv(n, x, incx)
CLACGV conjugates a complex vector.
Definition clacgv.f:74
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:106
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
real(wp) function scnrm2(n, x, incx)
SCNRM2
Definition scnrm2.f90:90
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: