![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLARGE 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE DLARGE( N, A, LDA, ISEED, WORK, INFO ) 00012 * 00013 * .. Scalar Arguments .. 00014 * INTEGER INFO, LDA, N 00015 * .. 00016 * .. Array Arguments .. 00017 * INTEGER ISEED( 4 ) 00018 * DOUBLE PRECISION A( LDA, * ), WORK( * ) 00019 * .. 00020 * 00021 * 00022 *> \par Purpose: 00023 * ============= 00024 *> 00025 *> \verbatim 00026 *> 00027 *> DLARGE pre- and post-multiplies a real general n by n matrix A 00028 *> with a random orthogonal matrix: A = U*D*U'. 00029 *> \endverbatim 00030 * 00031 * Arguments: 00032 * ========== 00033 * 00034 *> \param[in] N 00035 *> \verbatim 00036 *> N is INTEGER 00037 *> The order of the matrix A. N >= 0. 00038 *> \endverbatim 00039 *> 00040 *> \param[in,out] A 00041 *> \verbatim 00042 *> A is DOUBLE PRECISION array, dimension (LDA,N) 00043 *> On entry, the original n by n matrix A. 00044 *> On exit, A is overwritten by U*A*U' for some random 00045 *> orthogonal matrix U. 00046 *> \endverbatim 00047 *> 00048 *> \param[in] LDA 00049 *> \verbatim 00050 *> LDA is INTEGER 00051 *> The leading dimension of the array A. LDA >= N. 00052 *> \endverbatim 00053 *> 00054 *> \param[in,out] ISEED 00055 *> \verbatim 00056 *> ISEED is INTEGER array, dimension (4) 00057 *> On entry, the seed of the random number generator; the array 00058 *> elements must be between 0 and 4095, and ISEED(4) must be 00059 *> odd. 00060 *> On exit, the seed is updated. 00061 *> \endverbatim 00062 *> 00063 *> \param[out] WORK 00064 *> \verbatim 00065 *> WORK is DOUBLE PRECISION array, dimension (2*N) 00066 *> \endverbatim 00067 *> 00068 *> \param[out] INFO 00069 *> \verbatim 00070 *> INFO is INTEGER 00071 *> = 0: successful exit 00072 *> < 0: if INFO = -i, the i-th argument had an illegal value 00073 *> \endverbatim 00074 * 00075 * Authors: 00076 * ======== 00077 * 00078 *> \author Univ. of Tennessee 00079 *> \author Univ. of California Berkeley 00080 *> \author Univ. of Colorado Denver 00081 *> \author NAG Ltd. 00082 * 00083 *> \date November 2011 00084 * 00085 *> \ingroup double_matgen 00086 * 00087 * ===================================================================== 00088 SUBROUTINE DLARGE( N, A, LDA, ISEED, WORK, INFO ) 00089 * 00090 * -- LAPACK auxiliary routine (version 3.4.0) -- 00091 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00092 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00093 * November 2011 00094 * 00095 * .. Scalar Arguments .. 00096 INTEGER INFO, LDA, N 00097 * .. 00098 * .. Array Arguments .. 00099 INTEGER ISEED( 4 ) 00100 DOUBLE PRECISION A( LDA, * ), WORK( * ) 00101 * .. 00102 * 00103 * ===================================================================== 00104 * 00105 * .. Parameters .. 00106 DOUBLE PRECISION ZERO, ONE 00107 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 00108 * .. 00109 * .. Local Scalars .. 00110 INTEGER I 00111 DOUBLE PRECISION TAU, WA, WB, WN 00112 * .. 00113 * .. External Subroutines .. 00114 EXTERNAL DGEMV, DGER, DLARNV, DSCAL, XERBLA 00115 * .. 00116 * .. Intrinsic Functions .. 00117 INTRINSIC MAX, SIGN 00118 * .. 00119 * .. External Functions .. 00120 DOUBLE PRECISION DNRM2 00121 EXTERNAL DNRM2 00122 * .. 00123 * .. Executable Statements .. 00124 * 00125 * Test the input arguments 00126 * 00127 INFO = 0 00128 IF( N.LT.0 ) THEN 00129 INFO = -1 00130 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 00131 INFO = -3 00132 END IF 00133 IF( INFO.LT.0 ) THEN 00134 CALL XERBLA( 'DLARGE', -INFO ) 00135 RETURN 00136 END IF 00137 * 00138 * pre- and post-multiply A by random orthogonal matrix 00139 * 00140 DO 10 I = N, 1, -1 00141 * 00142 * generate random reflection 00143 * 00144 CALL DLARNV( 3, ISEED, N-I+1, WORK ) 00145 WN = DNRM2( N-I+1, WORK, 1 ) 00146 WA = SIGN( WN, WORK( 1 ) ) 00147 IF( WN.EQ.ZERO ) THEN 00148 TAU = ZERO 00149 ELSE 00150 WB = WORK( 1 ) + WA 00151 CALL DSCAL( N-I, ONE / WB, WORK( 2 ), 1 ) 00152 WORK( 1 ) = ONE 00153 TAU = WB / WA 00154 END IF 00155 * 00156 * multiply A(i:n,1:n) by random reflection from the left 00157 * 00158 CALL DGEMV( 'Transpose', N-I+1, N, ONE, A( I, 1 ), LDA, WORK, 00159 $ 1, ZERO, WORK( N+1 ), 1 ) 00160 CALL DGER( N-I+1, N, -TAU, WORK, 1, WORK( N+1 ), 1, A( I, 1 ), 00161 $ LDA ) 00162 * 00163 * multiply A(1:n,i:n) by random reflection from the right 00164 * 00165 CALL DGEMV( 'No transpose', N, N-I+1, ONE, A( 1, I ), LDA, 00166 $ WORK, 1, ZERO, WORK( N+1 ), 1 ) 00167 CALL DGER( N, N-I+1, -TAU, WORK( N+1 ), 1, WORK, 1, A( 1, I ), 00168 $ LDA ) 00169 10 CONTINUE 00170 RETURN 00171 * 00172 * End of DLARGE 00173 * 00174 END