1*> \brief \b ZQRT12 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* DOUBLE PRECISION FUNCTION ZQRT12( M, N, A, LDA, S, WORK, LWORK, 12* RWORK ) 13* 14* .. Scalar Arguments .. 15* INTEGER LDA, LWORK, M, N 16* .. 17* .. Array Arguments .. 18* DOUBLE PRECISION RWORK( * ), S( * ) 19* COMPLEX*16 A( LDA, * ), WORK( LWORK ) 20* .. 21* 22* 23*> \par Purpose: 24* ============= 25*> 26*> \verbatim 27*> 28*> ZQRT12 computes the singular values `svlues' of the upper trapezoid 29*> of A(1:M,1:N) and returns the ratio 30*> 31*> || s - svlues||/(||svlues||*eps*max(M,N)) 32*> \endverbatim 33* 34* Arguments: 35* ========== 36* 37*> \param[in] M 38*> \verbatim 39*> M is INTEGER 40*> The number of rows of the matrix A. 41*> \endverbatim 42*> 43*> \param[in] N 44*> \verbatim 45*> N is INTEGER 46*> The number of columns of the matrix A. 47*> \endverbatim 48*> 49*> \param[in] A 50*> \verbatim 51*> A is COMPLEX*16 array, dimension (LDA,N) 52*> The M-by-N matrix A. Only the upper trapezoid is referenced. 53*> \endverbatim 54*> 55*> \param[in] LDA 56*> \verbatim 57*> LDA is INTEGER 58*> The leading dimension of the array A. 59*> \endverbatim 60*> 61*> \param[in] S 62*> \verbatim 63*> S is DOUBLE PRECISION array, dimension (min(M,N)) 64*> The singular values of the matrix A. 65*> \endverbatim 66*> 67*> \param[out] WORK 68*> \verbatim 69*> WORK is COMPLEX*16 array, dimension (LWORK) 70*> \endverbatim 71*> 72*> \param[in] LWORK 73*> \verbatim 74*> LWORK is INTEGER 75*> The length of the array WORK. LWORK >= M*N + 2*min(M,N) + 76*> max(M,N). 77*> \endverbatim 78*> 79*> \param[out] RWORK 80*> \verbatim 81*> RWORK is DOUBLE PRECISION array, dimension (2*min(M,N)) 82*> \endverbatim 83* 84* Authors: 85* ======== 86* 87*> \author Univ. of Tennessee 88*> \author Univ. of California Berkeley 89*> \author Univ. of Colorado Denver 90*> \author NAG Ltd. 91* 92*> \ingroup complex16_lin 93* 94* ===================================================================== 95 DOUBLE PRECISION FUNCTION ZQRT12( M, N, A, LDA, S, WORK, LWORK, 96 $ RWORK ) 97* 98* -- LAPACK test routine -- 99* -- LAPACK is a software package provided by Univ. of Tennessee, -- 100* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 101* 102* .. Scalar Arguments .. 103 INTEGER LDA, LWORK, M, N 104* .. 105* .. Array Arguments .. 106 DOUBLE PRECISION RWORK( * ), S( * ) 107 COMPLEX*16 A( LDA, * ), WORK( LWORK ) 108* .. 109* 110* ===================================================================== 111* 112* .. Parameters .. 113 DOUBLE PRECISION ZERO, ONE 114 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) 115* .. 116* .. Local Scalars .. 117 INTEGER I, INFO, ISCL, J, MN 118 DOUBLE PRECISION ANRM, BIGNUM, NRMSVL, SMLNUM 119* .. 120* .. Local Arrays .. 121 DOUBLE PRECISION DUMMY( 1 ) 122* .. 123* .. External Functions .. 124 DOUBLE PRECISION DASUM, DLAMCH, DNRM2, ZLANGE 125 EXTERNAL DASUM, DLAMCH, DNRM2, ZLANGE 126* .. 127* .. External Subroutines .. 128 EXTERNAL DAXPY, DBDSQR, DLABAD, DLASCL, XERBLA, ZGEBD2, 129 $ ZLASCL, ZLASET 130* .. 131* .. Intrinsic Functions .. 132 INTRINSIC DBLE, DCMPLX, MAX, MIN 133* .. 134* .. Executable Statements .. 135* 136 ZQRT12 = ZERO 137* 138* Test that enough workspace is supplied 139* 140 IF( LWORK.LT.M*N+2*MIN( M, N )+MAX( M, N ) ) THEN 141 CALL XERBLA( 'ZQRT12', 7 ) 142 RETURN 143 END IF 144* 145* Quick return if possible 146* 147 MN = MIN( M, N ) 148 IF( MN.LE.ZERO ) 149 $ RETURN 150* 151 NRMSVL = DNRM2( MN, S, 1 ) 152* 153* Copy upper triangle of A into work 154* 155 CALL ZLASET( 'Full', M, N, DCMPLX( ZERO ), DCMPLX( ZERO ), WORK, 156 $ M ) 157 DO 20 J = 1, N 158 DO 10 I = 1, MIN( J, M ) 159 WORK( ( J-1 )*M+I ) = A( I, J ) 160 10 CONTINUE 161 20 CONTINUE 162* 163* Get machine parameters 164* 165 SMLNUM = DLAMCH( 'S' ) / DLAMCH( 'P' ) 166 BIGNUM = ONE / SMLNUM 167 CALL DLABAD( SMLNUM, BIGNUM ) 168* 169* Scale work if max entry outside range [SMLNUM,BIGNUM] 170* 171 ANRM = ZLANGE( 'M', M, N, WORK, M, DUMMY ) 172 ISCL = 0 173 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN 174* 175* Scale matrix norm up to SMLNUM 176* 177 CALL ZLASCL( 'G', 0, 0, ANRM, SMLNUM, M, N, WORK, M, INFO ) 178 ISCL = 1 179 ELSE IF( ANRM.GT.BIGNUM ) THEN 180* 181* Scale matrix norm down to BIGNUM 182* 183 CALL ZLASCL( 'G', 0, 0, ANRM, BIGNUM, M, N, WORK, M, INFO ) 184 ISCL = 1 185 END IF 186* 187 IF( ANRM.NE.ZERO ) THEN 188* 189* Compute SVD of work 190* 191 CALL ZGEBD2( M, N, WORK, M, RWORK( 1 ), RWORK( MN+1 ), 192 $ WORK( M*N+1 ), WORK( M*N+MN+1 ), 193 $ WORK( M*N+2*MN+1 ), INFO ) 194 CALL DBDSQR( 'Upper', MN, 0, 0, 0, RWORK( 1 ), RWORK( MN+1 ), 195 $ DUMMY, MN, DUMMY, 1, DUMMY, MN, RWORK( 2*MN+1 ), 196 $ INFO ) 197* 198 IF( ISCL.EQ.1 ) THEN 199 IF( ANRM.GT.BIGNUM ) THEN 200 CALL DLASCL( 'G', 0, 0, BIGNUM, ANRM, MN, 1, RWORK( 1 ), 201 $ MN, INFO ) 202 END IF 203 IF( ANRM.LT.SMLNUM ) THEN 204 CALL DLASCL( 'G', 0, 0, SMLNUM, ANRM, MN, 1, RWORK( 1 ), 205 $ MN, INFO ) 206 END IF 207 END IF 208* 209 ELSE 210* 211 DO 30 I = 1, MN 212 RWORK( I ) = ZERO 213 30 CONTINUE 214 END IF 215* 216* Compare s and singular values of work 217* 218 CALL DAXPY( MN, -ONE, S, 1, RWORK( 1 ), 1 ) 219 ZQRT12 = DASUM( MN, RWORK( 1 ), 1 ) / 220 $ ( DLAMCH( 'Epsilon' )*DBLE( MAX( M, N ) ) ) 221 IF( NRMSVL.NE.ZERO ) 222 $ ZQRT12 = ZQRT12 / NRMSVL 223* 224 RETURN 225* 226* End of ZQRT12 227* 228 END 229