1*> \brief \b ZQRT05 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE ZQRT05(M,N,L,NB,RESULT) 12* 13* .. Scalar Arguments .. 14* INTEGER LWORK, M, N, L, NB, LDT 15* .. Return values .. 16* DOUBLE PRECISION RESULT(6) 17* 18* 19*> \par Purpose: 20* ============= 21*> 22*> \verbatim 23*> 24*> ZQRT05 tests ZTPQRT and ZTPMQRT. 25*> \endverbatim 26* 27* Arguments: 28* ========== 29* 30*> \param[in] M 31*> \verbatim 32*> M is INTEGER 33*> Number of rows in lower part of the test matrix. 34*> \endverbatim 35*> 36*> \param[in] N 37*> \verbatim 38*> N is INTEGER 39*> Number of columns in test matrix. 40*> \endverbatim 41*> 42*> \param[in] L 43*> \verbatim 44*> L is INTEGER 45*> The number of rows of the upper trapezoidal part the 46*> lower test matrix. 0 <= L <= M. 47*> \endverbatim 48*> 49*> \param[in] NB 50*> \verbatim 51*> NB is INTEGER 52*> Block size of test matrix. NB <= N. 53*> \endverbatim 54*> 55*> \param[out] RESULT 56*> \verbatim 57*> RESULT is DOUBLE PRECISION array, dimension (6) 58*> Results of each of the six tests below. 59*> 60*> RESULT(1) = | A - Q R | 61*> RESULT(2) = | I - Q^H Q | 62*> RESULT(3) = | Q C - Q C | 63*> RESULT(4) = | Q^H C - Q^H C | 64*> RESULT(5) = | C Q - C Q | 65*> RESULT(6) = | C Q^H - C Q^H | 66*> \endverbatim 67* 68* Authors: 69* ======== 70* 71*> \author Univ. of Tennessee 72*> \author Univ. of California Berkeley 73*> \author Univ. of Colorado Denver 74*> \author NAG Ltd. 75* 76*> \date April 2012 77* 78*> \ingroup complex16_lin 79* 80* ===================================================================== 81 SUBROUTINE ZQRT05(M,N,L,NB,RESULT) 82 IMPLICIT NONE 83* 84* -- LAPACK test routine (version 3.4.1) -- 85* -- LAPACK is a software package provided by Univ. of Tennessee, -- 86* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 87* April 2012 88* 89* .. Scalar Arguments .. 90 INTEGER LWORK, M, N, L, NB, LDT 91* .. Return values .. 92 DOUBLE PRECISION RESULT(6) 93* 94* ===================================================================== 95* 96* .. 97* .. Local allocatable arrays 98 COMPLEX*16, ALLOCATABLE :: AF(:,:), Q(:,:), 99 $ R(:,:), RWORK(:), WORK( : ), T(:,:), 100 $ CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:) 101* 102* .. Parameters .. 103 DOUBLE PRECISION ZERO 104 COMPLEX*16 ONE, CZERO 105 PARAMETER( ZERO = 0.0, ONE = (1.0,0.0), CZERO=(0.0,0.0) ) 106* .. 107* .. Local Scalars .. 108 INTEGER INFO, J, K, M2, NP1 109 DOUBLE PRECISION ANORM, EPS, RESID, CNORM, DNORM 110* .. 111* .. Local Arrays .. 112 INTEGER ISEED( 4 ) 113* .. 114* .. External Functions .. 115 DOUBLE PRECISION DLAMCH 116 DOUBLE PRECISION ZLANGE, ZLANSY 117 LOGICAL LSAME 118 EXTERNAL DLAMCH, ZLANGE, ZLANSY, LSAME 119* .. 120* .. Data statements .. 121 DATA ISEED / 1988, 1989, 1990, 1991 / 122* 123 EPS = DLAMCH( 'Epsilon' ) 124 K = N 125 M2 = M+N 126 IF( M.GT.0 ) THEN 127 NP1 = N+1 128 ELSE 129 NP1 = 1 130 END IF 131 LWORK = M2*M2*NB 132* 133* Dynamically allocate all arrays 134* 135 ALLOCATE(A(M2,N),AF(M2,N),Q(M2,M2),R(M2,M2),RWORK(M2), 136 $ WORK(LWORK),T(NB,N),C(M2,N),CF(M2,N), 137 $ D(N,M2),DF(N,M2) ) 138* 139* Put random stuff into A 140* 141 LDT=NB 142 CALL ZLASET( 'Full', M2, N, CZERO, CZERO, A, M2 ) 143 CALL ZLASET( 'Full', NB, N, CZERO, CZERO, T, NB ) 144 DO J=1,N 145 CALL ZLARNV( 2, ISEED, J, A( 1, J ) ) 146 END DO 147 IF( M.GT.0 ) THEN 148 DO J=1,N 149 CALL ZLARNV( 2, ISEED, M-L, A( MIN(N+M,N+1), J ) ) 150 END DO 151 END IF 152 IF( L.GT.0 ) THEN 153 DO J=1,N 154 CALL ZLARNV( 2, ISEED, MIN(J,L), A( MIN(N+M,N+M-L+1), J ) ) 155 END DO 156 END IF 157* 158* Copy the matrix A to the array AF. 159* 160 CALL ZLACPY( 'Full', M2, N, A, M2, AF, M2 ) 161* 162* Factor the matrix A in the array AF. 163* 164 CALL ZTPQRT( M,N,L,NB,AF,M2,AF(NP1,1),M2,T,LDT,WORK,INFO) 165* 166* Generate the (M+N)-by-(M+N) matrix Q by applying H to I 167* 168 CALL ZLASET( 'Full', M2, M2, CZERO, ONE, Q, M2 ) 169 CALL ZGEMQRT( 'R', 'N', M2, M2, K, NB, AF, M2, T, LDT, Q, M2, 170 $ WORK, INFO ) 171* 172* Copy R 173* 174 CALL ZLASET( 'Full', M2, N, CZERO, CZERO, R, M2 ) 175 CALL ZLACPY( 'Upper', M2, N, AF, M2, R, M2 ) 176* 177* Compute |R - Q'*A| / |A| and store in RESULT(1) 178* 179 CALL ZGEMM( 'C', 'N', M2, N, M2, -ONE, Q, M2, A, M2, ONE, R, M2 ) 180 ANORM = ZLANGE( '1', M2, N, A, M2, RWORK ) 181 RESID = ZLANGE( '1', M2, N, R, M2, RWORK ) 182 IF( ANORM.GT.ZERO ) THEN 183 RESULT( 1 ) = RESID / (EPS*ANORM*MAX(1,M2)) 184 ELSE 185 RESULT( 1 ) = ZERO 186 END IF 187* 188* Compute |I - Q'*Q| and store in RESULT(2) 189* 190 CALL ZLASET( 'Full', M2, M2, CZERO, ONE, R, M2 ) 191 CALL ZHERK( 'U', 'C', M2, M2, DREAL(-ONE), Q, M2, DREAL(ONE), 192 $ R, M2 ) 193 RESID = ZLANSY( '1', 'Upper', M2, R, M2, RWORK ) 194 RESULT( 2 ) = RESID / (EPS*MAX(1,M2)) 195* 196* Generate random m-by-n matrix C and a copy CF 197* 198 DO J=1,N 199 CALL ZLARNV( 2, ISEED, M2, C( 1, J ) ) 200 END DO 201 CNORM = ZLANGE( '1', M2, N, C, M2, RWORK) 202 CALL ZLACPY( 'Full', M2, N, C, M2, CF, M2 ) 203* 204* Apply Q to C as Q*C 205* 206 CALL ZTPMQRT( 'L','N', M,N,K,L,NB,AF(NP1,1),M2,T,LDT,CF,M2, 207 $ CF(NP1,1),M2,WORK,INFO) 208* 209* Compute |Q*C - Q*C| / |C| 210* 211 CALL ZGEMM( 'N', 'N', M2, N, M2, -ONE, Q, M2, C, M2, ONE, CF, M2 ) 212 RESID = ZLANGE( '1', M2, N, CF, M2, RWORK ) 213 IF( CNORM.GT.ZERO ) THEN 214 RESULT( 3 ) = RESID / (EPS*MAX(1,M2)*CNORM) 215 ELSE 216 RESULT( 3 ) = ZERO 217 END IF 218* 219* Copy C into CF again 220* 221 CALL ZLACPY( 'Full', M2, N, C, M2, CF, M2 ) 222* 223* Apply Q to C as QT*C 224* 225 CALL ZTPMQRT( 'L','C',M,N,K,L,NB,AF(NP1,1),M2,T,LDT,CF,M2, 226 $ CF(NP1,1),M2,WORK,INFO) 227* 228* Compute |QT*C - QT*C| / |C| 229* 230 CALL ZGEMM('C','N',M2,N,M2,-ONE,Q,M2,C,M2,ONE,CF,M2) 231 RESID = ZLANGE( '1', M2, N, CF, M2, RWORK ) 232 IF( CNORM.GT.ZERO ) THEN 233 RESULT( 4 ) = RESID / (EPS*MAX(1,M2)*CNORM) 234 ELSE 235 RESULT( 4 ) = ZERO 236 END IF 237* 238* Generate random n-by-m matrix D and a copy DF 239* 240 DO J=1,M2 241 CALL ZLARNV( 2, ISEED, N, D( 1, J ) ) 242 END DO 243 DNORM = ZLANGE( '1', N, M2, D, N, RWORK) 244 CALL ZLACPY( 'Full', N, M2, D, N, DF, N ) 245* 246* Apply Q to D as D*Q 247* 248 CALL ZTPMQRT('R','N',N,M,N,L,NB,AF(NP1,1),M2,T,LDT,DF,N, 249 $ DF(1,NP1),N,WORK,INFO) 250* 251* Compute |D*Q - D*Q| / |D| 252* 253 CALL ZGEMM('N','N',N,M2,M2,-ONE,D,N,Q,M2,ONE,DF,N) 254 RESID = ZLANGE('1',N, M2,DF,N,RWORK ) 255 IF( CNORM.GT.ZERO ) THEN 256 RESULT( 5 ) = RESID / (EPS*MAX(1,M2)*DNORM) 257 ELSE 258 RESULT( 5 ) = ZERO 259 END IF 260* 261* Copy D into DF again 262* 263 CALL ZLACPY('Full',N,M2,D,N,DF,N ) 264* 265* Apply Q to D as D*QT 266* 267 CALL ZTPMQRT('R','C',N,M,N,L,NB,AF(NP1,1),M2,T,LDT,DF,N, 268 $ DF(1,NP1),N,WORK,INFO) 269 270* 271* Compute |D*QT - D*QT| / |D| 272* 273 CALL ZGEMM( 'N', 'C', N, M2, M2, -ONE, D, N, Q, M2, ONE, DF, N ) 274 RESID = ZLANGE( '1', N, M2, DF, N, RWORK ) 275 IF( CNORM.GT.ZERO ) THEN 276 RESULT( 6 ) = RESID / (EPS*MAX(1,M2)*DNORM) 277 ELSE 278 RESULT( 6 ) = ZERO 279 END IF 280* 281* Deallocate all arrays 282* 283 DEALLOCATE ( A, AF, Q, R, RWORK, WORK, T, C, D, CF, DF) 284 RETURN 285 END 286 287