 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE ZQRT15( SCALE, RKSEL, M, N, NRHS, A, LDA, B, LDB, S, 00002 $ RANK, NORMA, NORMB, ISEED, WORK, LWORK ) 00003 * 00004 * -- LAPACK test routine (version 3.1) -- 00005 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 INTEGER LDA, LDB, LWORK, M, N, NRHS, RANK, RKSEL, SCALE 00010 DOUBLE PRECISION NORMA, NORMB 00011 * .. 00012 * .. Array Arguments .. 00013 INTEGER ISEED( 4 ) 00014 DOUBLE PRECISION S( * ) 00015 COMPLEX*16 A( LDA, * ), B( LDB, * ), WORK( LWORK ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * ZQRT15 generates a matrix with full or deficient rank and of various 00022 * norms. 00023 * 00024 * Arguments 00025 * ========= 00026 * 00027 * SCALE (input) INTEGER 00028 * SCALE = 1: normally scaled matrix 00029 * SCALE = 2: matrix scaled up 00030 * SCALE = 3: matrix scaled down 00031 * 00032 * RKSEL (input) INTEGER 00033 * RKSEL = 1: full rank matrix 00034 * RKSEL = 2: rank-deficient matrix 00035 * 00036 * M (input) INTEGER 00037 * The number of rows of the matrix A. 00038 * 00039 * N (input) INTEGER 00040 * The number of columns of A. 00041 * 00042 * NRHS (input) INTEGER 00043 * The number of columns of B. 00044 * 00045 * A (output) COMPLEX*16 array, dimension (LDA,N) 00046 * The M-by-N matrix A. 00047 * 00048 * LDA (input) INTEGER 00049 * The leading dimension of the array A. 00050 * 00051 * B (output) COMPLEX*16 array, dimension (LDB, NRHS) 00052 * A matrix that is in the range space of matrix A. 00053 * 00054 * LDB (input) INTEGER 00055 * The leading dimension of the array B. 00056 * 00057 * S (output) DOUBLE PRECISION array, dimension MIN(M,N) 00058 * Singular values of A. 00059 * 00060 * RANK (output) INTEGER 00061 * number of nonzero singular values of A. 00062 * 00063 * NORMA (output) DOUBLE PRECISION 00064 * one-norm norm of A. 00065 * 00066 * NORMB (output) DOUBLE PRECISION 00067 * one-norm norm of B. 00068 * 00069 * ISEED (input/output) integer array, dimension (4) 00070 * seed for random number generator. 00071 * 00072 * WORK (workspace) COMPLEX*16 array, dimension (LWORK) 00073 * 00074 * LWORK (input) INTEGER 00075 * length of work space required. 00076 * LWORK >= MAX(M+MIN(M,N),NRHS*MIN(M,N),2*N+M) 00077 * 00078 * ===================================================================== 00079 * 00080 * .. Parameters .. 00081 DOUBLE PRECISION ZERO, ONE, TWO, SVMIN 00082 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0, 00083 $ SVMIN = 0.1D+0 ) 00084 COMPLEX*16 CZERO, CONE 00085 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ), 00086 $ CONE = ( 1.0D+0, 0.0D+0 ) ) 00087 * .. 00088 * .. Local Scalars .. 00089 INTEGER INFO, J, MN 00090 DOUBLE PRECISION BIGNUM, EPS, SMLNUM, TEMP 00091 * .. 00092 * .. Local Arrays .. 00093 DOUBLE PRECISION DUMMY( 1 ) 00094 * .. 00095 * .. External Functions .. 00096 DOUBLE PRECISION DASUM, DLAMCH, DLARND, DZNRM2, ZLANGE 00097 EXTERNAL DASUM, DLAMCH, DLARND, DZNRM2, ZLANGE 00098 * .. 00099 * .. External Subroutines .. 00100 EXTERNAL DLABAD, DLAORD, DLASCL, XERBLA, ZDSCAL, ZGEMM, 00101 $ ZLARF, ZLARNV, ZLAROR, ZLASCL, ZLASET 00102 * .. 00103 * .. Intrinsic Functions .. 00104 INTRINSIC ABS, DCMPLX, MAX, MIN 00105 * .. 00106 * .. Executable Statements .. 00107 * 00108 MN = MIN( M, N ) 00109 IF( LWORK.LT.MAX( M+MN, MN*NRHS, 2*N+M ) ) THEN 00110 CALL XERBLA( 'ZQRT15', 16 ) 00111 RETURN 00112 END IF 00113 * 00114 SMLNUM = DLAMCH( 'Safe minimum' ) 00115 BIGNUM = ONE / SMLNUM 00116 CALL DLABAD( SMLNUM, BIGNUM ) 00117 EPS = DLAMCH( 'Epsilon' ) 00118 SMLNUM = ( SMLNUM / EPS ) / EPS 00119 BIGNUM = ONE / SMLNUM 00120 * 00121 * Determine rank and (unscaled) singular values 00122 * 00123 IF( RKSEL.EQ.1 ) THEN 00124 RANK = MN 00125 ELSE IF( RKSEL.EQ.2 ) THEN 00126 RANK = ( 3*MN ) / 4 00127 DO 10 J = RANK + 1, MN 00128 S( J ) = ZERO 00129 10 CONTINUE 00130 ELSE 00131 CALL XERBLA( 'ZQRT15', 2 ) 00132 END IF 00133 * 00134 IF( RANK.GT.0 ) THEN 00135 * 00136 * Nontrivial case 00137 * 00138 S( 1 ) = ONE 00139 DO 30 J = 2, RANK 00140 20 CONTINUE 00141 TEMP = DLARND( 1, ISEED ) 00142 IF( TEMP.GT.SVMIN ) THEN 00143 S( J ) = ABS( TEMP ) 00144 ELSE 00145 GO TO 20 00146 END IF 00147 30 CONTINUE 00148 CALL DLAORD( 'Decreasing', RANK, S, 1 ) 00149 * 00150 * Generate 'rank' columns of a random orthogonal matrix in A 00151 * 00152 CALL ZLARNV( 2, ISEED, M, WORK ) 00153 CALL ZDSCAL( M, ONE / DZNRM2( M, WORK, 1 ), WORK, 1 ) 00154 CALL ZLASET( 'Full', M, RANK, CZERO, CONE, A, LDA ) 00155 CALL ZLARF( 'Left', M, RANK, WORK, 1, DCMPLX( TWO ), A, LDA, 00156 $ WORK( M+1 ) ) 00157 * 00158 * workspace used: m+mn 00159 * 00160 * Generate consistent rhs in the range space of A 00161 * 00162 CALL ZLARNV( 2, ISEED, RANK*NRHS, WORK ) 00163 CALL ZGEMM( 'No transpose', 'No transpose', M, NRHS, RANK, 00164 $ CONE, A, LDA, WORK, RANK, CZERO, B, LDB ) 00165 * 00166 * work space used: <= mn *nrhs 00167 * 00168 * generate (unscaled) matrix A 00169 * 00170 DO 40 J = 1, RANK 00171 CALL ZDSCAL( M, S( J ), A( 1, J ), 1 ) 00172 40 CONTINUE 00173 IF( RANK.LT.N ) 00174 $ CALL ZLASET( 'Full', M, N-RANK, CZERO, CZERO, 00175 $ A( 1, RANK+1 ), LDA ) 00176 CALL ZLAROR( 'Right', 'No initialization', M, N, A, LDA, ISEED, 00177 $ WORK, INFO ) 00178 * 00179 ELSE 00180 * 00181 * work space used 2*n+m 00182 * 00183 * Generate null matrix and rhs 00184 * 00185 DO 50 J = 1, MN 00186 S( J ) = ZERO 00187 50 CONTINUE 00188 CALL ZLASET( 'Full', M, N, CZERO, CZERO, A, LDA ) 00189 CALL ZLASET( 'Full', M, NRHS, CZERO, CZERO, B, LDB ) 00190 * 00191 END IF 00192 * 00193 * Scale the matrix 00194 * 00195 IF( SCALE.NE.1 ) THEN 00196 NORMA = ZLANGE( 'Max', M, N, A, LDA, DUMMY ) 00197 IF( NORMA.NE.ZERO ) THEN 00198 IF( SCALE.EQ.2 ) THEN 00199 * 00200 * matrix scaled up 00201 * 00202 CALL ZLASCL( 'General', 0, 0, NORMA, BIGNUM, M, N, A, 00203 $ LDA, INFO ) 00204 CALL DLASCL( 'General', 0, 0, NORMA, BIGNUM, MN, 1, S, 00205 $ MN, INFO ) 00206 CALL ZLASCL( 'General', 0, 0, NORMA, BIGNUM, M, NRHS, B, 00207 $ LDB, INFO ) 00208 ELSE IF( SCALE.EQ.3 ) THEN 00209 * 00210 * matrix scaled down 00211 * 00212 CALL ZLASCL( 'General', 0, 0, NORMA, SMLNUM, M, N, A, 00213 $ LDA, INFO ) 00214 CALL DLASCL( 'General', 0, 0, NORMA, SMLNUM, MN, 1, S, 00215 $ MN, INFO ) 00216 CALL ZLASCL( 'General', 0, 0, NORMA, SMLNUM, M, NRHS, B, 00217 $ LDB, INFO ) 00218 ELSE 00219 CALL XERBLA( 'ZQRT15', 1 ) 00220 RETURN 00221 END IF 00222 END IF 00223 END IF 00224 * 00225 NORMA = DASUM( MN, S, 1 ) 00226 NORMB = ZLANGE( 'One-norm', M, NRHS, B, LDB, DUMMY ) 00227 * 00228 RETURN 00229 * 00230 * End of ZQRT15 00231 * 00232 END
1.7.3 
