LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE CLARGE( N, A, LDA, ISEED, WORK, INFO ) 00002 * 00003 * -- LAPACK auxiliary test routine (version 3.1) -- 00004 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00005 * November 2006 00006 * 00007 * .. Scalar Arguments .. 00008 INTEGER INFO, LDA, N 00009 * .. 00010 * .. Array Arguments .. 00011 INTEGER ISEED( 4 ) 00012 COMPLEX A( LDA, * ), WORK( * ) 00013 * .. 00014 * 00015 * Purpose 00016 * ======= 00017 * 00018 * CLARGE pre- and post-multiplies a complex general n by n matrix A 00019 * with a random unitary matrix: A = U*D*U'. 00020 * 00021 * Arguments 00022 * ========= 00023 * 00024 * N (input) INTEGER 00025 * The order of the matrix A. N >= 0. 00026 * 00027 * A (input/output) COMPLEX array, dimension (LDA,N) 00028 * On entry, the original n by n matrix A. 00029 * On exit, A is overwritten by U*A*U' for some random 00030 * unitary matrix U. 00031 * 00032 * LDA (input) INTEGER 00033 * The leading dimension of the array A. LDA >= N. 00034 * 00035 * ISEED (input/output) INTEGER array, dimension (4) 00036 * On entry, the seed of the random number generator; the array 00037 * elements must be between 0 and 4095, and ISEED(4) must be 00038 * odd. 00039 * On exit, the seed is updated. 00040 * 00041 * WORK (workspace) COMPLEX array, dimension (2*N) 00042 * 00043 * INFO (output) INTEGER 00044 * = 0: successful exit 00045 * < 0: if INFO = -i, the i-th argument had an illegal value 00046 * 00047 * ===================================================================== 00048 * 00049 * .. Parameters .. 00050 COMPLEX ZERO, ONE 00051 PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ), 00052 $ ONE = ( 1.0E+0, 0.0E+0 ) ) 00053 * .. 00054 * .. Local Scalars .. 00055 INTEGER I 00056 REAL WN 00057 COMPLEX TAU, WA, WB 00058 * .. 00059 * .. External Subroutines .. 00060 EXTERNAL CGEMV, CGERC, CLARNV, CSCAL, XERBLA 00061 * .. 00062 * .. Intrinsic Functions .. 00063 INTRINSIC ABS, MAX, REAL 00064 * .. 00065 * .. External Functions .. 00066 REAL SCNRM2 00067 EXTERNAL SCNRM2 00068 * .. 00069 * .. Executable Statements .. 00070 * 00071 * Test the input arguments 00072 * 00073 INFO = 0 00074 IF( N.LT.0 ) THEN 00075 INFO = -1 00076 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 00077 INFO = -3 00078 END IF 00079 IF( INFO.LT.0 ) THEN 00080 CALL XERBLA( 'CLARGE', -INFO ) 00081 RETURN 00082 END IF 00083 * 00084 * pre- and post-multiply A by random unitary matrix 00085 * 00086 DO 10 I = N, 1, -1 00087 * 00088 * generate random reflection 00089 * 00090 CALL CLARNV( 3, ISEED, N-I+1, WORK ) 00091 WN = SCNRM2( N-I+1, WORK, 1 ) 00092 WA = ( WN / ABS( WORK( 1 ) ) )*WORK( 1 ) 00093 IF( WN.EQ.ZERO ) THEN 00094 TAU = ZERO 00095 ELSE 00096 WB = WORK( 1 ) + WA 00097 CALL CSCAL( N-I, ONE / WB, WORK( 2 ), 1 ) 00098 WORK( 1 ) = ONE 00099 TAU = REAL( WB / WA ) 00100 END IF 00101 * 00102 * multiply A(i:n,1:n) by random reflection from the left 00103 * 00104 CALL CGEMV( 'Conjugate transpose', N-I+1, N, ONE, A( I, 1 ), 00105 $ LDA, WORK, 1, ZERO, WORK( N+1 ), 1 ) 00106 CALL CGERC( N-I+1, N, -TAU, WORK, 1, WORK( N+1 ), 1, A( I, 1 ), 00107 $ LDA ) 00108 * 00109 * multiply A(1:n,i:n) by random reflection from the right 00110 * 00111 CALL CGEMV( 'No transpose', N, N-I+1, ONE, A( 1, I ), LDA, 00112 $ WORK, 1, ZERO, WORK( N+1 ), 1 ) 00113 CALL CGERC( N, N-I+1, -TAU, WORK( N+1 ), 1, WORK, 1, A( 1, I ), 00114 $ LDA ) 00115 10 CONTINUE 00116 RETURN 00117 * 00118 * End of CLARGE 00119 * 00120 END