LAPACK  3.10.1
LAPACK: Linear Algebra PACKage

◆ 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)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine cscal(N, CA, CX, INCX)
CSCAL
Definition: cscal.f:78
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
Definition: cgemv.f:158
subroutine cgerc(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
CGERC
Definition: cgerc.f:130
complex function clarnd(IDIST, ISEED)
CLARND
Definition: clarnd.f:75
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
real(wp) function scnrm2(n, x, incx)
SCNRM2
Definition: scnrm2.f90:90
Here is the call graph for this function:
Here is the caller graph for this function: