1*> \brief \b CTSQR01 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 CTSQR01(TSSW, M,N, MB, NB, RESULT) 12* 13* .. Scalar Arguments .. 14* INTEGER M, N, MB 15* .. Return values .. 16* REAL RESULT(6) 17* 18* 19*> \par Purpose: 20* ============= 21*> 22*> \verbatim 23*> 24*> DTSQR01 tests DGEQR , DGELQ, DGEMLQ and DGEMQR. 25*> \endverbatim 26* 27* Arguments: 28* ========== 29* 30*> \param[in] TSSW 31*> \verbatim 32*> TSSW is CHARACTER 33*> 'TS' for testing tall skinny QR 34*> and anything else for testing short wide LQ 35*> \endverbatim 36*> \param[in] M 37*> \verbatim 38*> M is INTEGER 39*> Number of rows in test matrix. 40*> \endverbatim 41*> 42*> \param[in] N 43*> \verbatim 44*> N is INTEGER 45*> Number of columns in test matrix. 46*> \endverbatim 47*> \param[in] MB 48*> \verbatim 49*> MB is INTEGER 50*> Number of row in row block in test matrix. 51*> \endverbatim 52*> 53*> \param[in] NB 54*> \verbatim 55*> NB is INTEGER 56*> Number of columns in column block test matrix. 57*> \endverbatim 58*> 59*> \param[out] RESULT 60*> \verbatim 61*> RESULT is REAL array, dimension (6) 62*> Results of each of the six tests below. 63*> 64*> RESULT(1) = | A - Q R | or | A - L Q | 65*> RESULT(2) = | I - Q^H Q | or | I - Q Q^H | 66*> RESULT(3) = | Q C - Q C | 67*> RESULT(4) = | Q^H C - Q^H C | 68*> RESULT(5) = | C Q - C Q | 69*> RESULT(6) = | C Q^H - C Q^H | 70*> \endverbatim 71* 72* Authors: 73* ======== 74* 75*> \author Univ. of Tennessee 76*> \author Univ. of California Berkeley 77*> \author Univ. of Colorado Denver 78*> \author NAG Ltd. 79* 80* ===================================================================== 81 SUBROUTINE CTSQR01(TSSW, M, N, MB, NB, RESULT) 82 IMPLICIT NONE 83* 84* -- LAPACK test routine -- 85* -- LAPACK is a software package provided by Univ. of Tennessee, -- 86* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 87* 88* .. Scalar Arguments .. 89 CHARACTER TSSW 90 INTEGER M, N, MB, NB 91* .. Return values .. 92 REAL RESULT(6) 93* 94* ===================================================================== 95* 96* .. 97* .. Local allocatable arrays 98 COMPLEX, ALLOCATABLE :: AF(:,:), Q(:,:), 99 $ R(:,:), RWORK(:), WORK( : ), T(:), 100 $ CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:), LQ(:,:) 101* 102* .. Parameters .. 103 REAL ZERO 104 COMPLEX ONE, CZERO 105 PARAMETER( ZERO = 0.0, ONE = (1.0,0.0), CZERO=(0.0,0.0) ) 106* .. 107* .. Local Scalars .. 108 LOGICAL TESTZEROS, TS 109 INTEGER INFO, J, K, L, LWORK, TSIZE, MNB 110 REAL ANORM, EPS, RESID, CNORM, DNORM 111* .. 112* .. Local Arrays .. 113 INTEGER ISEED( 4 ) 114 COMPLEX TQUERY( 5 ), WORKQUERY( 1 ) 115* .. 116* .. External Functions .. 117 REAL SLAMCH, CLANGE, CLANSY 118 LOGICAL LSAME 119 INTEGER ILAENV 120 EXTERNAL SLAMCH, CLANGE, CLANSY, LSAME, ILAENV 121* .. 122* .. Intrinsic Functions .. 123 INTRINSIC MAX, MIN 124* .. Scalars in Common .. 125 CHARACTER*32 srnamt 126* .. 127* .. Common blocks .. 128 COMMON / srnamc / srnamt 129* .. 130* .. Data statements .. 131 DATA ISEED / 1988, 1989, 1990, 1991 / 132* 133* TEST TALL SKINNY OR SHORT WIDE 134* 135 TS = LSAME(TSSW, 'TS') 136* 137* TEST MATRICES WITH HALF OF MATRIX BEING ZEROS 138* 139 TESTZEROS = .FALSE. 140* 141 EPS = SLAMCH( 'Epsilon' ) 142 K = MIN(M,N) 143 L = MAX(M,N,1) 144 MNB = MAX ( MB, NB) 145 LWORK = MAX(3,L)*MNB 146* 147* Dynamically allocate local arrays 148* 149 ALLOCATE ( A(M,N), AF(M,N), Q(L,L), R(M,L), RWORK(L), 150 $ C(M,N), CF(M,N), 151 $ D(N,M), DF(N,M), LQ(L,N) ) 152* 153* Put random numbers into A and copy to AF 154* 155 DO J=1,N 156 CALL CLARNV( 2, ISEED, M, A( 1, J ) ) 157 END DO 158 IF (TESTZEROS) THEN 159 IF (M.GE.4) THEN 160 DO J=1,N 161 CALL CLARNV( 2, ISEED, M/2, A( M/4, J ) ) 162 END DO 163 END IF 164 END IF 165 CALL CLACPY( 'Full', M, N, A, M, AF, M ) 166* 167 IF (TS) THEN 168* 169* Factor the matrix A in the array AF. 170* 171 CALL CGEQR( M, N, AF, M, TQUERY, -1, WORKQUERY, -1, INFO ) 172 TSIZE = INT( TQUERY( 1 ) ) 173 LWORK = INT( WORKQUERY( 1 ) ) 174 CALL CGEMQR( 'L', 'N', M, M, K, AF, M, TQUERY, TSIZE, CF, M, 175 $ WORKQUERY, -1, INFO) 176 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 177 CALL CGEMQR( 'L', 'N', M, N, K, AF, M, TQUERY, TSIZE, CF, M, 178 $ WORKQUERY, -1, INFO) 179 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 180 CALL CGEMQR( 'L', 'C', M, N, K, AF, M, TQUERY, TSIZE, CF, M, 181 $ WORKQUERY, -1, INFO) 182 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 183 CALL CGEMQR( 'R', 'N', N, M, K, AF, M, TQUERY, TSIZE, DF, N, 184 $ WORKQUERY, -1, INFO) 185 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 186 CALL CGEMQR( 'R', 'C', N, M, K, AF, M, TQUERY, TSIZE, DF, N, 187 $ WORKQUERY, -1, INFO) 188 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 189 ALLOCATE ( T( TSIZE ) ) 190 ALLOCATE ( WORK( LWORK ) ) 191 srnamt = 'CGEQR' 192 CALL CGEQR( M, N, AF, M, T, TSIZE, WORK, LWORK, INFO ) 193* 194* Generate the m-by-m matrix Q 195* 196 CALL CLASET( 'Full', M, M, CZERO, ONE, Q, M ) 197 srnamt = 'CGEMQR' 198 CALL CGEMQR( 'L', 'N', M, M, K, AF, M, T, TSIZE, Q, M, 199 $ WORK, LWORK, INFO ) 200* 201* Copy R 202* 203 CALL CLASET( 'Full', M, N, CZERO, CZERO, R, M ) 204 CALL CLACPY( 'Upper', M, N, AF, M, R, M ) 205* 206* Compute |R - Q'*A| / |A| and store in RESULT(1) 207* 208 CALL CGEMM( 'C', 'N', M, N, M, -ONE, Q, M, A, M, ONE, R, M ) 209 ANORM = CLANGE( '1', M, N, A, M, RWORK ) 210 RESID = CLANGE( '1', M, N, R, M, RWORK ) 211 IF( ANORM.GT.ZERO ) THEN 212 RESULT( 1 ) = RESID / (EPS*MAX(1,M)*ANORM) 213 ELSE 214 RESULT( 1 ) = ZERO 215 END IF 216* 217* Compute |I - Q'*Q| and store in RESULT(2) 218* 219 CALL CLASET( 'Full', M, M, CZERO, ONE, R, M ) 220 CALL CHERK( 'U', 'C', M, M, REAL(-ONE), Q, M, REAL(ONE), R, M ) 221 RESID = CLANSY( '1', 'Upper', M, R, M, RWORK ) 222 RESULT( 2 ) = RESID / (EPS*MAX(1,M)) 223* 224* Generate random m-by-n matrix C and a copy CF 225* 226 DO J=1,N 227 CALL CLARNV( 2, ISEED, M, C( 1, J ) ) 228 END DO 229 CNORM = CLANGE( '1', M, N, C, M, RWORK) 230 CALL CLACPY( 'Full', M, N, C, M, CF, M ) 231* 232* Apply Q to C as Q*C 233* 234 srnamt = 'CGEMQR' 235 CALL CGEMQR( 'L', 'N', M, N, K, AF, M, T, TSIZE, CF, M, 236 $ WORK, LWORK, INFO) 237* 238* Compute |Q*C - Q*C| / |C| 239* 240 CALL CGEMM( 'N', 'N', M, N, M, -ONE, Q, M, C, M, ONE, CF, M ) 241 RESID = CLANGE( '1', M, N, CF, M, RWORK ) 242 IF( CNORM.GT.ZERO ) THEN 243 RESULT( 3 ) = RESID / (EPS*MAX(1,M)*CNORM) 244 ELSE 245 RESULT( 3 ) = ZERO 246 END IF 247* 248* Copy C into CF again 249* 250 CALL CLACPY( 'Full', M, N, C, M, CF, M ) 251* 252* Apply Q to C as QT*C 253* 254 srnamt = 'CGEMQR' 255 CALL CGEMQR( 'L', 'C', M, N, K, AF, M, T, TSIZE, CF, M, 256 $ WORK, LWORK, INFO) 257* 258* Compute |QT*C - QT*C| / |C| 259* 260 CALL CGEMM( 'C', 'N', M, N, M, -ONE, Q, M, C, M, ONE, CF, M ) 261 RESID = CLANGE( '1', M, N, CF, M, RWORK ) 262 IF( CNORM.GT.ZERO ) THEN 263 RESULT( 4 ) = RESID / (EPS*MAX(1,M)*CNORM) 264 ELSE 265 RESULT( 4 ) = ZERO 266 END IF 267* 268* Generate random n-by-m matrix D and a copy DF 269* 270 DO J=1,M 271 CALL CLARNV( 2, ISEED, N, D( 1, J ) ) 272 END DO 273 DNORM = CLANGE( '1', N, M, D, N, RWORK) 274 CALL CLACPY( 'Full', N, M, D, N, DF, N ) 275* 276* Apply Q to D as D*Q 277* 278 srnamt = 'CGEMQR' 279 CALL CGEMQR( 'R', 'N', N, M, K, AF, M, T, TSIZE, DF, N, 280 $ WORK, LWORK, INFO) 281* 282* Compute |D*Q - D*Q| / |D| 283* 284 CALL CGEMM( 'N', 'N', N, M, M, -ONE, D, N, Q, M, ONE, DF, N ) 285 RESID = CLANGE( '1', N, M, DF, N, RWORK ) 286 IF( DNORM.GT.ZERO ) THEN 287 RESULT( 5 ) = RESID / (EPS*MAX(1,M)*DNORM) 288 ELSE 289 RESULT( 5 ) = ZERO 290 END IF 291* 292* Copy D into DF again 293* 294 CALL CLACPY( 'Full', N, M, D, N, DF, N ) 295* 296* Apply Q to D as D*QT 297* 298 CALL CGEMQR( 'R', 'C', N, M, K, AF, M, T, TSIZE, DF, N, 299 $ WORK, LWORK, INFO) 300* 301* Compute |D*QT - D*QT| / |D| 302* 303 CALL CGEMM( 'N', 'C', N, M, M, -ONE, D, N, Q, M, ONE, DF, N ) 304 RESID = CLANGE( '1', N, M, DF, N, RWORK ) 305 IF( CNORM.GT.ZERO ) THEN 306 RESULT( 6 ) = RESID / (EPS*MAX(1,M)*DNORM) 307 ELSE 308 RESULT( 6 ) = ZERO 309 END IF 310* 311* Short and wide 312* 313 ELSE 314 CALL CGELQ( M, N, AF, M, TQUERY, -1, WORKQUERY, -1, INFO ) 315 TSIZE = INT( TQUERY( 1 ) ) 316 LWORK = INT( WORKQUERY( 1 ) ) 317 CALL CGEMLQ( 'R', 'N', N, N, K, AF, M, TQUERY, TSIZE, Q, N, 318 $ WORKQUERY, -1, INFO ) 319 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 320 CALL CGEMLQ( 'L', 'N', N, M, K, AF, M, TQUERY, TSIZE, DF, N, 321 $ WORKQUERY, -1, INFO) 322 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 323 CALL CGEMLQ( 'L', 'C', N, M, K, AF, M, TQUERY, TSIZE, DF, N, 324 $ WORKQUERY, -1, INFO) 325 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 326 CALL CGEMLQ( 'R', 'N', M, N, K, AF, M, TQUERY, TSIZE, CF, M, 327 $ WORKQUERY, -1, INFO) 328 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 329 CALL CGEMLQ( 'R', 'C', M, N, K, AF, M, TQUERY, TSIZE, CF, M, 330 $ WORKQUERY, -1, INFO) 331 LWORK = MAX( LWORK, INT( WORKQUERY( 1 ) ) ) 332 ALLOCATE ( T( TSIZE ) ) 333 ALLOCATE ( WORK( LWORK ) ) 334 srnamt = 'CGELQ' 335 CALL CGELQ( M, N, AF, M, T, TSIZE, WORK, LWORK, INFO ) 336* 337* 338* Generate the n-by-n matrix Q 339* 340 CALL CLASET( 'Full', N, N, CZERO, ONE, Q, N ) 341 srnamt = 'CGEMLQ' 342 CALL CGEMLQ( 'R', 'N', N, N, K, AF, M, T, TSIZE, Q, N, 343 $ WORK, LWORK, INFO ) 344* 345* Copy R 346* 347 CALL CLASET( 'Full', M, N, CZERO, CZERO, LQ, L ) 348 CALL CLACPY( 'Lower', M, N, AF, M, LQ, L ) 349* 350* Compute |L - A*Q'| / |A| and store in RESULT(1) 351* 352 CALL CGEMM( 'N', 'C', M, N, N, -ONE, A, M, Q, N, ONE, LQ, L ) 353 ANORM = CLANGE( '1', M, N, A, M, RWORK ) 354 RESID = CLANGE( '1', M, N, LQ, L, RWORK ) 355 IF( ANORM.GT.ZERO ) THEN 356 RESULT( 1 ) = RESID / (EPS*MAX(1,N)*ANORM) 357 ELSE 358 RESULT( 1 ) = ZERO 359 END IF 360* 361* Compute |I - Q'*Q| and store in RESULT(2) 362* 363 CALL CLASET( 'Full', N, N, CZERO, ONE, LQ, L ) 364 CALL CHERK( 'U', 'C', N, N, REAL(-ONE), Q, N, REAL(ONE), LQ, L) 365 RESID = CLANSY( '1', 'Upper', N, LQ, L, RWORK ) 366 RESULT( 2 ) = RESID / (EPS*MAX(1,N)) 367* 368* Generate random m-by-n matrix C and a copy CF 369* 370 DO J=1,M 371 CALL CLARNV( 2, ISEED, N, D( 1, J ) ) 372 END DO 373 DNORM = CLANGE( '1', N, M, D, N, RWORK) 374 CALL CLACPY( 'Full', N, M, D, N, DF, N ) 375* 376* Apply Q to C as Q*C 377* 378 CALL CGEMLQ( 'L', 'N', N, M, K, AF, M, T, TSIZE, DF, N, 379 $ WORK, LWORK, INFO) 380* 381* Compute |Q*D - Q*D| / |D| 382* 383 CALL CGEMM( 'N', 'N', N, M, N, -ONE, Q, N, D, N, ONE, DF, N ) 384 RESID = CLANGE( '1', N, M, DF, N, RWORK ) 385 IF( DNORM.GT.ZERO ) THEN 386 RESULT( 3 ) = RESID / (EPS*MAX(1,N)*DNORM) 387 ELSE 388 RESULT( 3 ) = ZERO 389 END IF 390* 391* Copy D into DF again 392* 393 CALL CLACPY( 'Full', N, M, D, N, DF, N ) 394* 395* Apply Q to D as QT*D 396* 397 CALL CGEMLQ( 'L', 'C', N, M, K, AF, M, T, TSIZE, DF, N, 398 $ WORK, LWORK, INFO) 399* 400* Compute |QT*D - QT*D| / |D| 401* 402 CALL CGEMM( 'C', 'N', N, M, N, -ONE, Q, N, D, N, ONE, DF, N ) 403 RESID = CLANGE( '1', N, M, DF, N, RWORK ) 404 IF( DNORM.GT.ZERO ) THEN 405 RESULT( 4 ) = RESID / (EPS*MAX(1,N)*DNORM) 406 ELSE 407 RESULT( 4 ) = ZERO 408 END IF 409* 410* Generate random n-by-m matrix D and a copy DF 411* 412 DO J=1,N 413 CALL CLARNV( 2, ISEED, M, C( 1, J ) ) 414 END DO 415 CNORM = CLANGE( '1', M, N, C, M, RWORK) 416 CALL CLACPY( 'Full', M, N, C, M, CF, M ) 417* 418* Apply Q to C as C*Q 419* 420 CALL CGEMLQ( 'R', 'N', M, N, K, AF, M, T, TSIZE, CF, M, 421 $ WORK, LWORK, INFO) 422* 423* Compute |C*Q - C*Q| / |C| 424* 425 CALL CGEMM( 'N', 'N', M, N, N, -ONE, C, M, Q, N, ONE, CF, M ) 426 RESID = CLANGE( '1', N, M, DF, N, RWORK ) 427 IF( CNORM.GT.ZERO ) THEN 428 RESULT( 5 ) = RESID / (EPS*MAX(1,N)*CNORM) 429 ELSE 430 RESULT( 5 ) = ZERO 431 END IF 432* 433* Copy C into CF again 434* 435 CALL CLACPY( 'Full', M, N, C, M, CF, M ) 436* 437* Apply Q to D as D*QT 438* 439 CALL CGEMLQ( 'R', 'C', M, N, K, AF, M, T, TSIZE, CF, M, 440 $ WORK, LWORK, INFO) 441* 442* Compute |C*QT - C*QT| / |C| 443* 444 CALL CGEMM( 'N', 'C', M, N, N, -ONE, C, M, Q, N, ONE, CF, M ) 445 RESID = CLANGE( '1', M, N, CF, M, RWORK ) 446 IF( CNORM.GT.ZERO ) THEN 447 RESULT( 6 ) = RESID / (EPS*MAX(1,N)*CNORM) 448 ELSE 449 RESULT( 6 ) = ZERO 450 END IF 451* 452 END IF 453* 454* Deallocate all arrays 455* 456 DEALLOCATE ( A, AF, Q, R, RWORK, WORK, T, C, D, CF, DF) 457* 458 RETURN 459 END 460