1 PROGRAM PCQRDRIVER 2* 3* -- ScaLAPACK testing driver (version 1.7) -- 4* University of Tennessee, Knoxville, Oak Ridge National Laboratory, 5* and University of California, Berkeley. 6* May 28, 2001 7* 8* Purpose 9* ======= 10* 11* PCQRDRIVER is the main test program for the COMPLEX 12* SCALAPACK QR factorization routines. This test driver performs a QR 13* QL, LQ, RQ, QP (QR factorization with column pivoting) or TZ 14* (complete unitary factorization) factorization and checks the 15* results. 16* 17* The program must be driven by a short data file. An annotated 18* example of a data file can be obtained by deleting the first 3 19* characters from the following 16 lines: 20* 'ScaLAPACK QR factorizations input file' 21* 'PVM machine' 22* 'QR.out' output file name (if any) 23* 6 device out 24* 6 number of factorizations 25* 'QR' 'QL' 'LQ' 'RQ' 'QP' 'TZ' factorization: QR, QL, LQ, RQ, QP, TZ 26* 4 number of problems sizes 27* 55 17 31 201 values of M 28* 5 71 31 201 values of N 29* 3 number of MB's and NB's 30* 4 3 5 values of MB 31* 4 7 3 values of NB 32* 7 number of process grids (ordered P & Q) 33* 1 2 1 4 2 3 8 values of P 34* 7 2 4 1 3 2 1 values of Q 35* 1.0 threshold 36* 37* Internal Parameters 38* =================== 39* 40* TOTMEM INTEGER, default = 2000000 41* TOTMEM is a machine-specific parameter indicating the 42* maximum amount of available memory in bytes. 43* The user should customize TOTMEM to his platform. Remember 44* to leave room in memory for the operating system, the BLACS 45* buffer, etc. For example, on a system with 8 MB of memory 46* per process (e.g., one processor on an Intel iPSC/860), the 47* parameters we use are TOTMEM=6200000 (leaving 1.8 MB for OS, 48* code, BLACS buffer, etc). However, for PVM, we usually set 49* TOTMEM = 2000000. Some experimenting with the maximum value 50* of TOTMEM may be required. 51* 52* INTGSZ INTEGER, default = 4 bytes. 53* REALSZ INTEGER, default = 4 bytes. 54* CPLXSZ INTEGER, default = 8 bytes. 55* INTGSZ, REALSZ and CPLXSZ indicate the length in bytes on 56* the given platform for an integer, a single precision real 57* and a single precision complex. 58* MEM COMPLEX array, dimension ( TOTMEM / CPLXSZ ) 59* 60* All arrays used by SCALAPACK routines are allocated from 61* this array and referenced by pointers. The integer IPA, 62* for example, is a pointer to the starting element of MEM for 63* the matrix A. 64* 65* ===================================================================== 66* 67* .. Parameters .. 68 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_, 69 $ LLD_, MB_, M_, NB_, N_, RSRC_ 70 PARAMETER ( BLOCK_CYCLIC_2D = 1, DLEN_ = 9, DTYPE_ = 1, 71 $ CTXT_ = 2, M_ = 3, N_ = 4, MB_ = 5, NB_ = 6, 72 $ RSRC_ = 7, CSRC_ = 8, LLD_ = 9 ) 73 INTEGER CPLXSZ, INTGSZ, MEMSIZ, NTESTS, REALSZ, TOTMEM 74 COMPLEX PADVAL 75 PARAMETER ( CPLXSZ = 8, INTGSZ = 4, REALSZ = 4, 76 $ TOTMEM = 2000000, MEMSIZ = TOTMEM / CPLXSZ, 77 $ NTESTS = 20, 78 $ PADVAL = ( -9923.0E+0, -9923.0E+0 ) ) 79* .. 80* .. Local Scalars .. 81 CHARACTER*2 FACT 82 CHARACTER*6 PASSED 83 CHARACTER*7 ROUT 84 CHARACTER*8 ROUTCHK 85 CHARACTER*80 OUTFILE 86 LOGICAL CHECK 87 INTEGER I, IAM, IASEED, ICTXT, IMIDPAD, INFO, IPA, 88 $ IPOSTPAD, IPPIV, IPREPAD, IPTAU, IPRW, IPW, J, 89 $ K, KFAIL, KPASS, KSKIP, KTESTS, L, LIPIV, 90 $ LRWORK, LTAU, LWORK, M, MAXMN, MB, MINMN, MNP, 91 $ MNQ, MP, MYCOL, MYROW, N, NB, NFACT, NGRIDS, 92 $ NMAT, NNB, NOUT, NPCOL, NPROCS, NPROW, NQ, 93 $ WORKFCT, WORKRFCT, WORKSIZ 94 REAL ANORM, FRESID, THRESH 95 DOUBLE PRECISION NOPS, TMFLOPS 96* .. 97* .. Arrays .. 98 CHARACTER*2 FACTOR( NTESTS ) 99 INTEGER DESCA( DLEN_ ), IERR( 1 ), MBVAL( NTESTS ), 100 $ MVAL( NTESTS ), NBVAL( NTESTS ), 101 $ NVAL( NTESTS ), PVAL( NTESTS ), QVAL( NTESTS ) 102 DOUBLE PRECISION CTIME( 1 ), WTIME( 1 ) 103 COMPLEX MEM( MEMSIZ ) 104* .. 105* .. External Subroutines .. 106 EXTERNAL BLACS_BARRIER, BLACS_EXIT, BLACS_GET, 107 $ BLACS_GRIDEXIT, BLACS_GRIDINFO, BLACS_GRIDINIT, 108 $ BLACS_PINFO, DESCINIT, IGSUM2D, PCCHEKPAD, 109 $ PCFILLPAD, PCGELQF, PCGELQRV, 110 $ PCGEQLF, PCGEQLRV, PCGEQPF, 111 $ PCQPPIV, PCGEQRF, PCGEQRRV, 112 $ PCGERQF, PCGERQRV, PCTZRZRV, 113 $ PCMATGEN, PCLAFCHK, PCQRINFO, 114 $ PCTZRZF, SLBOOT, SLCOMBINE, SLTIMER 115* .. 116* .. External Functions .. 117 LOGICAL LSAMEN 118 INTEGER ICEIL, NUMROC 119 REAL PCLANGE 120 EXTERNAL ICEIL, LSAMEN, NUMROC, PCLANGE 121* .. 122* .. Intrinsic Functions .. 123 INTRINSIC DBLE, MAX, MIN 124* .. 125* .. Data Statements .. 126 DATA KTESTS, KPASS, KFAIL, KSKIP /4*0/ 127* .. 128* .. Executable Statements .. 129* 130* Get starting information 131* 132 CALL BLACS_PINFO( IAM, NPROCS ) 133 IASEED = 100 134 CALL PCQRINFO( OUTFILE, NOUT, NFACT, FACTOR, NTESTS, NMAT, MVAL, 135 $ NTESTS, NVAL, NTESTS, NNB, MBVAL, NTESTS, NBVAL, 136 $ NTESTS, NGRIDS, PVAL, NTESTS, QVAL, NTESTS, 137 $ THRESH, MEM, IAM, NPROCS ) 138 CHECK = ( THRESH.GE.0.0E+0 ) 139* 140* Loop over the different factorization types 141* 142 DO 40 I = 1, NFACT 143* 144 FACT = FACTOR( I ) 145* 146* Print headings 147* 148 IF( IAM.EQ.0 ) THEN 149 WRITE( NOUT, FMT = * ) 150 IF( LSAMEN( 2, FACT, 'QR' ) ) THEN 151 ROUT = 'PCGEQRF' 152 ROUTCHK = 'PCGEQRRV' 153 WRITE( NOUT, FMT = 9986 ) 154 $ 'QR factorization tests.' 155 ELSE IF( LSAMEN( 2, FACT, 'QL' ) ) THEN 156 ROUT = 'PCGEQLF' 157 ROUTCHK = 'PCGEQLRV' 158 WRITE( NOUT, FMT = 9986 ) 159 $ 'QL factorization tests.' 160 ELSE IF( LSAMEN( 2, FACT, 'LQ' ) ) THEN 161 ROUT = 'PCGELQF' 162 ROUTCHK = 'PCGELQRV' 163 WRITE( NOUT, FMT = 9986 ) 164 $ 'LQ factorization tests.' 165 ELSE IF( LSAMEN( 2, FACT, 'RQ' ) ) THEN 166 ROUT = 'PCGERQF' 167 ROUTCHK = 'PCGERQRV' 168 WRITE( NOUT, FMT = 9986 ) 169 $ 'RQ factorization tests.' 170 ELSE IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 171 ROUT = 'PCGEQPF' 172 ROUTCHK = 'PCGEQRRV' 173 WRITE( NOUT, FMT = 9986 ) 174 $ 'QR factorization with column pivoting tests.' 175 ELSE IF( LSAMEN( 2, FACT, 'TZ' ) ) THEN 176 ROUT = 'PCTZRZF' 177 ROUTCHK = 'PCTZRZRV' 178 WRITE( NOUT, FMT = 9986 ) 179 $ 'Complete unitary factorization tests.' 180 END IF 181 WRITE( NOUT, FMT = * ) 182 WRITE( NOUT, FMT = 9995 ) 183 WRITE( NOUT, FMT = 9994 ) 184 WRITE( NOUT, FMT = * ) 185 END IF 186* 187* Loop over different process grids 188* 189 DO 30 J = 1, NGRIDS 190* 191 NPROW = PVAL( J ) 192 NPCOL = QVAL( J ) 193* 194* Make sure grid information is correct 195* 196 IERR( 1 ) = 0 197 IF( NPROW.LT.1 ) THEN 198 IF( IAM.EQ.0 ) 199 $ WRITE( NOUT, FMT = 9999 ) 'GRID', 'nprow', NPROW 200 IERR( 1 ) = 1 201 ELSE IF( NPCOL.LT.1 ) THEN 202 IF( IAM.EQ.0 ) 203 $ WRITE( NOUT, FMT = 9999 ) 'GRID', 'npcol', NPCOL 204 IERR( 1 ) = 1 205 ELSE IF( NPROW*NPCOL.GT.NPROCS ) THEN 206 IF( IAM.EQ.0 ) 207 $ WRITE( NOUT, FMT = 9998 ) NPROW*NPCOL, NPROCS 208 IERR( 1 ) = 1 209 END IF 210* 211 IF( IERR( 1 ).GT.0 ) THEN 212 IF( IAM.EQ.0 ) 213 $ WRITE( NOUT, FMT = 9997 ) 'grid' 214 KSKIP = KSKIP + 1 215 GO TO 30 216 END IF 217* 218* Define process grid 219* 220 CALL BLACS_GET( -1, 0, ICTXT ) 221 CALL BLACS_GRIDINIT( ICTXT, 'Row-major', NPROW, NPCOL ) 222 CALL BLACS_GRIDINFO( ICTXT, NPROW, NPCOL, MYROW, MYCOL ) 223* 224* Go to bottom of loop if this case doesn't use my process 225* 226 IF( MYROW.GE.NPROW .OR. MYCOL.GE.NPCOL ) 227 $ GO TO 30 228* 229 DO 20 K = 1, NMAT 230* 231 M = MVAL( K ) 232 N = NVAL( K ) 233* 234* Make sure matrix information is correct 235* 236 IERR(1) = 0 237 IF( M.LT.1 ) THEN 238 IF( IAM.EQ.0 ) 239 $ WRITE( NOUT, FMT = 9999 ) 'MATRIX', 'M', M 240 IERR( 1 ) = 1 241 ELSE IF( N.LT.1 ) THEN 242 IF( IAM.EQ.0 ) 243 $ WRITE( NOUT, FMT = 9999 ) 'MATRIX', 'N', N 244 IERR( 1 ) = 1 245 END IF 246* 247* Make sure no one had error 248* 249 CALL IGSUM2D( ICTXT, 'All', ' ', 1, 1, IERR, 1, -1, 0 ) 250* 251 IF( IERR( 1 ).GT.0 ) THEN 252 IF( IAM.EQ.0 ) 253 $ WRITE( NOUT, FMT = 9997 ) 'matrix' 254 KSKIP = KSKIP + 1 255 GO TO 20 256 END IF 257* 258* Loop over different blocking sizes 259* 260 DO 10 L = 1, NNB 261* 262 MB = MBVAL( L ) 263 NB = NBVAL( L ) 264* 265* Make sure mb is legal 266* 267 IERR( 1 ) = 0 268 IF( MB.LT.1 ) THEN 269 IERR( 1 ) = 1 270 IF( IAM.EQ.0 ) 271 $ WRITE( NOUT, FMT = 9999 ) 'MB', 'MB', MB 272 END IF 273* 274* Check all processes for an error 275* 276 CALL IGSUM2D( ICTXT, 'All', ' ', 1, 1, IERR, 1, -1, 277 $ 0 ) 278* 279 IF( IERR( 1 ).GT.0 ) THEN 280 IF( IAM.EQ.0 ) 281 $ WRITE( NOUT, FMT = 9997 ) 'MB' 282 KSKIP = KSKIP + 1 283 GO TO 10 284 END IF 285* 286* Make sure nb is legal 287* 288 IERR( 1 ) = 0 289 IF( NB.LT.1 ) THEN 290 IERR( 1 ) = 1 291 IF( IAM.EQ.0 ) 292 $ WRITE( NOUT, FMT = 9999 ) 'NB', 'NB', NB 293 END IF 294* 295* Check all processes for an error 296* 297 CALL IGSUM2D( ICTXT, 'All', ' ', 1, 1, IERR, 1, -1, 298 $ 0 ) 299* 300 IF( IERR( 1 ).GT.0 ) THEN 301 IF( IAM.EQ.0 ) 302 $ WRITE( NOUT, FMT = 9997 ) 'NB' 303 KSKIP = KSKIP + 1 304 GO TO 10 305 END IF 306* 307* Padding constants 308* 309 MP = NUMROC( M, MB, MYROW, 0, NPROW ) 310 NQ = NUMROC( N, NB, MYCOL, 0, NPCOL ) 311 MNP = NUMROC( MIN( M, N ), MB, MYROW, 0, NPROW ) 312 MNQ = NUMROC( MIN( M, N ), NB, MYCOL, 0, NPCOL ) 313 IF( CHECK ) THEN 314 IPREPAD = MAX( MB, MP ) 315 IMIDPAD = NB 316 IPOSTPAD = MAX( NB, NQ ) 317 ELSE 318 IPREPAD = 0 319 IMIDPAD = 0 320 IPOSTPAD = 0 321 END IF 322* 323* Initialize the array descriptor for the matrix A 324* 325 CALL DESCINIT( DESCA, M, N, MB, NB, 0, 0, ICTXT, 326 $ MAX( 1, MP ) + IMIDPAD, IERR( 1 ) ) 327* 328* Check all processes for an error 329* 330 CALL IGSUM2D( ICTXT, 'All', ' ', 1, 1, IERR, 1, -1, 331 $ 0 ) 332* 333 IF( IERR( 1 ).LT.0 ) THEN 334 IF( IAM.EQ.0 ) 335 $ WRITE( NOUT, FMT = 9997 ) 'descriptor' 336 KSKIP = KSKIP + 1 337 GO TO 10 338 END IF 339* 340* Assign pointers into MEM for ScaLAPACK arrays, A is 341* allocated starting at position MEM( IPREPAD+1 ) 342* 343 IPA = IPREPAD+1 344 IPTAU = IPA + DESCA( LLD_ ) * NQ + IPOSTPAD + IPREPAD 345* 346 IF( LSAMEN( 2, FACT, 'QR' ) ) THEN 347* 348 LTAU = MNQ 349 IPW = IPTAU + LTAU + IPOSTPAD + IPREPAD 350* 351* Figure the amount of workspace required by the QR 352* factorization 353* 354 LWORK = DESCA( NB_ ) * ( MP + NQ + DESCA( NB_ ) ) 355 WORKFCT = LWORK + IPOSTPAD 356 WORKSIZ = WORKFCT 357* 358 IF( CHECK ) THEN 359* 360* Figure the amount of workspace required by the 361* checking routines PCLAFCHK, PCGEQRRV and 362* PCLANGE 363* 364 WORKSIZ = LWORK + MP*DESCA( NB_ ) + IPOSTPAD 365* 366 END IF 367* 368 ELSE IF( LSAMEN( 2, FACT, 'QL' ) ) THEN 369* 370 LTAU = NQ 371 IPW = IPTAU + LTAU + IPOSTPAD + IPREPAD 372* 373* Figure the amount of workspace required by the QL 374* factorization 375* 376 LWORK = DESCA( NB_ ) * ( MP + NQ + DESCA( NB_ ) ) 377 WORKFCT = LWORK + IPOSTPAD 378 WORKSIZ = WORKFCT 379* 380 IF( CHECK ) THEN 381* 382* Figure the amount of workspace required by the 383* checking routines PCLAFCHK, PCGEQLRV and 384* PCLANGE 385* 386 WORKSIZ = LWORK + MP*DESCA( NB_ ) + IPOSTPAD 387* 388 END IF 389* 390 ELSE IF( LSAMEN( 2, FACT, 'LQ' ) ) THEN 391* 392 LTAU = MNP 393 IPW = IPTAU + LTAU + IPOSTPAD + IPREPAD 394* 395* Figure the amount of workspace required by the LQ 396* factorization 397* 398 LWORK = DESCA( MB_ ) * ( MP + NQ + DESCA( MB_ ) ) 399 WORKFCT = LWORK + IPOSTPAD 400 WORKSIZ = WORKFCT 401* 402 IF( CHECK ) THEN 403* 404* Figure the amount of workspace required by the 405* checking routines PCLAFCHK, PCGELQRV and 406* PCLANGE 407* 408 WORKSIZ = LWORK + 409 $ MAX( MP*DESCA( NB_ ), NQ*DESCA( MB_ ) 410 $ ) + IPOSTPAD 411* 412 END IF 413* 414 ELSE IF( LSAMEN( 2, FACT, 'RQ' ) ) THEN 415* 416 LTAU = MP 417 IPW = IPTAU + LTAU + IPOSTPAD + IPREPAD 418* 419* Figure the amount of workspace required by the QR 420* factorization 421* 422 LWORK = DESCA( MB_ ) * ( MP + NQ + DESCA( MB_ ) ) 423 WORKFCT = LWORK + IPOSTPAD 424 WORKSIZ = WORKFCT 425* 426 IF( CHECK ) THEN 427* 428* Figure the amount of workspace required by the 429* checking routines PCLAFCHK, PCGERQRV and 430* PCLANGE 431* 432 WORKSIZ = LWORK + 433 $ MAX( MP*DESCA( NB_ ), NQ*DESCA( MB_ ) 434 $ ) + IPOSTPAD 435* 436 END IF 437* 438 ELSE IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 439* 440 LTAU = MNQ 441 IPPIV = IPTAU + LTAU + IPOSTPAD + IPREPAD 442 LIPIV = ICEIL( INTGSZ*NQ, CPLXSZ ) 443 IPW = IPPIV + LIPIV + IPOSTPAD + IPREPAD 444* 445* Figure the amount of workspace required by the 446* factorization i.e from IPW on. 447* 448 LWORK = MAX( 3, MP + MAX( 1, NQ ) ) 449 WORKFCT = LWORK + IPOSTPAD 450 LRWORK = MAX( 1, 2 * NQ ) 451 WORKRFCT = ICEIL( LRWORK*REALSZ, CPLXSZ ) + 452 $ IPOSTPAD 453 IPRW = IPW + WORKFCT + IPREPAD 454 WORKSIZ = WORKFCT + IPREPAD + WORKRFCT 455* 456 IF( CHECK ) THEN 457* 458* Figure the amount of workspace required by the 459* checking routines PCLAFCHK, PCGEQRRV, 460* PCLANGE. 461* 462 WORKSIZ = MAX( WORKSIZ - IPOSTPAD, 463 $ DESCA( NB_ )*( 2*MP + NQ + DESCA( NB_ ) ) ) + 464 $ IPOSTPAD 465 END IF 466* 467 ELSE IF( LSAMEN( 2, FACT, 'TZ' ) ) THEN 468* 469 LTAU = MP 470 IPW = IPTAU + LTAU + IPOSTPAD + IPREPAD 471* 472* Figure the amount of workspace required by the TZ 473* factorization 474* 475 LWORK = DESCA( MB_ ) * ( MP + NQ + DESCA( MB_ ) ) 476 WORKFCT = LWORK + IPOSTPAD 477 WORKSIZ = WORKFCT 478* 479 IF( CHECK ) THEN 480* 481* Figure the amount of workspace required by the 482* checking routines PCLAFCHK, PCTZRZRV and 483* PCLANGE 484* 485 WORKSIZ = LWORK + 486 $ MAX( MP*DESCA( NB_ ), NQ*DESCA( MB_ ) 487 $ ) + IPOSTPAD 488* 489 END IF 490* 491 END IF 492* 493* Check for adequate memory for problem size 494* 495 IERR( 1 ) = 0 496 IF( IPW+WORKSIZ.GT.MEMSIZ ) THEN 497 IF( IAM.EQ.0 ) 498 $ WRITE( NOUT, FMT = 9996 ) 499 $ FACT // ' factorization', 500 $ ( IPW+WORKSIZ )*CPLXSZ 501 IERR( 1 ) = 1 502 END IF 503* 504* Check all processes for an error 505* 506 CALL IGSUM2D( ICTXT, 'All', ' ', 1, 1, IERR, 1, -1, 507 $ 0 ) 508* 509 IF( IERR( 1 ).GT.0 ) THEN 510 IF( IAM.EQ.0 ) 511 $ WRITE( NOUT, FMT = 9997 ) 'MEMORY' 512 KSKIP = KSKIP + 1 513 GO TO 10 514 END IF 515* 516* Generate the matrix A 517* 518 CALL PCMATGEN( ICTXT, 'N', 'N', DESCA( M_ ), 519 $ DESCA( N_ ), DESCA( MB_ ), 520 $ DESCA( NB_ ), MEM( IPA ), 521 $ DESCA( LLD_ ), DESCA( RSRC_ ), 522 $ DESCA( CSRC_ ), IASEED, 0, MP, 0, NQ, 523 $ MYROW, MYCOL, NPROW, NPCOL ) 524* 525* Need the Infinity of A for checking 526* 527 IF( CHECK ) THEN 528 CALL PCFILLPAD( ICTXT, MP, NQ, MEM( IPA-IPREPAD ), 529 $ DESCA( LLD_ ), IPREPAD, IPOSTPAD, 530 $ PADVAL ) 531 IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 532 CALL PCFILLPAD( ICTXT, LIPIV, 1, 533 $ MEM( IPPIV-IPREPAD ), LIPIV, 534 $ IPREPAD, IPOSTPAD, PADVAL ) 535 END IF 536 CALL PCFILLPAD( ICTXT, LTAU, 1, 537 $ MEM( IPTAU-IPREPAD ), LTAU, 538 $ IPREPAD, IPOSTPAD, PADVAL ) 539 CALL PCFILLPAD( ICTXT, WORKSIZ-IPOSTPAD, 1, 540 $ MEM( IPW-IPREPAD ), 541 $ WORKSIZ-IPOSTPAD, 542 $ IPREPAD, IPOSTPAD, PADVAL ) 543 ANORM = PCLANGE( 'I', M, N, MEM( IPA ), 1, 1, 544 $ DESCA, MEM( IPW ) ) 545 CALL PCCHEKPAD( ICTXT, 'PCLANGE', MP, NQ, 546 $ MEM( IPA-IPREPAD ), DESCA( LLD_ ), 547 $ IPREPAD, IPOSTPAD, PADVAL ) 548 CALL PCCHEKPAD( ICTXT, 'PCLANGE', 549 $ WORKSIZ-IPOSTPAD, 1, 550 $ MEM( IPW-IPREPAD ), 551 $ WORKSIZ-IPOSTPAD, IPREPAD, 552 $ IPOSTPAD, PADVAL ) 553 IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 554 CALL PCFILLPAD( ICTXT, WORKRFCT-IPOSTPAD, 1, 555 $ MEM( IPRW-IPREPAD ), 556 $ WORKRFCT-IPOSTPAD, 557 $ IPREPAD, IPOSTPAD, PADVAL ) 558 END IF 559 CALL PCFILLPAD( ICTXT, WORKFCT-IPOSTPAD, 1, 560 $ MEM( IPW-IPREPAD ), 561 $ WORKFCT-IPOSTPAD, 562 $ IPREPAD, IPOSTPAD, PADVAL ) 563 END IF 564* 565 CALL SLBOOT() 566 CALL BLACS_BARRIER( ICTXT, 'All' ) 567* 568* Perform QR factorizations 569* 570 IF( LSAMEN( 2, FACT, 'QR' ) ) THEN 571 CALL SLTIMER( 1 ) 572 CALL PCGEQRF( M, N, MEM( IPA ), 1, 1, DESCA, 573 $ MEM( IPTAU ), MEM( IPW ), LWORK, 574 $ INFO ) 575 CALL SLTIMER( 1 ) 576 ELSE IF( LSAMEN( 2, FACT, 'QL' ) ) THEN 577 CALL SLTIMER( 1 ) 578 CALL PCGEQLF( M, N, MEM( IPA ), 1, 1, DESCA, 579 $ MEM( IPTAU ), MEM( IPW ), LWORK, 580 $ INFO ) 581 CALL SLTIMER( 1 ) 582 ELSE IF( LSAMEN( 2, FACT, 'LQ' ) ) THEN 583 CALL SLTIMER( 1 ) 584 CALL PCGELQF( M, N, MEM( IPA ), 1, 1, DESCA, 585 $ MEM( IPTAU ), MEM( IPW ), LWORK, 586 $ INFO ) 587 CALL SLTIMER( 1 ) 588 ELSE IF( LSAMEN( 2, FACT, 'RQ' ) ) THEN 589 CALL SLTIMER( 1 ) 590 CALL PCGERQF( M, N, MEM( IPA ), 1, 1, DESCA, 591 $ MEM( IPTAU ), MEM( IPW ), LWORK, 592 $ INFO ) 593 CALL SLTIMER( 1 ) 594 ELSE IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 595 CALL SLTIMER( 1 ) 596 CALL PCGEQPF( M, N, MEM( IPA ), 1, 1, DESCA, 597 $ MEM( IPPIV ), MEM( IPTAU ), 598 $ MEM( IPW ), LWORK, MEM( IPRW ), 599 $ LRWORK, INFO ) 600 CALL SLTIMER( 1 ) 601 ELSE IF( LSAMEN( 2, FACT, 'TZ' ) ) THEN 602 CALL SLTIMER( 1 ) 603 IF( N.GE.M ) 604 $ CALL PCTZRZF( M, N, MEM( IPA ), 1, 1, DESCA, 605 $ MEM( IPTAU ), MEM( IPW ), LWORK, 606 $ INFO ) 607 CALL SLTIMER( 1 ) 608 END IF 609* 610 IF( CHECK ) THEN 611* 612* Check for memory overwrite in factorization 613* 614 CALL PCCHEKPAD( ICTXT, ROUT, MP, NQ, 615 $ MEM( IPA-IPREPAD ), DESCA( LLD_ ), 616 $ IPREPAD, IPOSTPAD, PADVAL ) 617 CALL PCCHEKPAD( ICTXT, ROUT, LTAU, 1, 618 $ MEM( IPTAU-IPREPAD ), LTAU, 619 $ IPREPAD, IPOSTPAD, PADVAL ) 620 IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 621 CALL PCCHEKPAD( ICTXT, ROUT, LIPIV, 1, 622 $ MEM( IPPIV-IPREPAD ), LIPIV, 623 $ IPREPAD, IPOSTPAD, PADVAL ) 624 CALL PCCHEKPAD( ICTXT, ROUT, WORKRFCT-IPOSTPAD, 625 $ 1, MEM( IPRW-IPREPAD ), 626 $ WORKRFCT-IPOSTPAD, 627 $ IPREPAD, IPOSTPAD, PADVAL ) 628 END IF 629 CALL PCCHEKPAD( ICTXT, ROUT, WORKFCT-IPOSTPAD, 1, 630 $ MEM( IPW-IPREPAD ), 631 $ WORKFCT-IPOSTPAD, IPREPAD, 632 $ IPOSTPAD, PADVAL ) 633 CALL PCFILLPAD( ICTXT, WORKSIZ-IPOSTPAD, 1, 634 $ MEM( IPW-IPREPAD ), 635 $ WORKSIZ-IPOSTPAD, 636 $ IPREPAD, IPOSTPAD, PADVAL ) 637* 638 IF( LSAMEN( 2, FACT, 'QR' ) ) THEN 639* 640* Compute residual = ||A-Q*R|| / (||A||*N*eps) 641* 642 CALL PCGEQRRV( M, N, MEM( IPA ), 1, 1, DESCA, 643 $ MEM( IPTAU ), MEM( IPW ) ) 644 CALL PCLAFCHK( 'No', 'No', M, N, MEM( IPA ), 1, 645 $ 1, DESCA, IASEED, ANORM, FRESID, 646 $ MEM( IPW ) ) 647 ELSE IF( LSAMEN( 2, FACT, 'QL' ) ) THEN 648* 649* Compute residual = ||A-Q*L|| / (||A||*N*eps) 650* 651 CALL PCGEQLRV( M, N, MEM( IPA ), 1, 1, DESCA, 652 $ MEM( IPTAU ), MEM( IPW ) ) 653 CALL PCLAFCHK( 'No', 'No', M, N, MEM( IPA ), 1, 654 $ 1, DESCA, IASEED, ANORM, FRESID, 655 $ MEM( IPW ) ) 656 ELSE IF( LSAMEN( 2, FACT, 'LQ' ) ) THEN 657* 658* Compute residual = ||A-L*Q|| / (||A||*N*eps) 659* 660 CALL PCGELQRV( M, N, MEM( IPA ), 1, 1, DESCA, 661 $ MEM( IPTAU ), MEM( IPW ) ) 662 CALL PCLAFCHK( 'No', 'No', M, N, MEM( IPA ), 1, 663 $ 1, DESCA, IASEED, ANORM, FRESID, 664 $ MEM( IPW ) ) 665 ELSE IF( LSAMEN( 2, FACT, 'RQ' ) ) THEN 666* 667* Compute residual = ||A-R*Q|| / (||A||*N*eps) 668* 669 CALL PCGERQRV( M, N, MEM( IPA ), 1, 1, DESCA, 670 $ MEM( IPTAU ), MEM( IPW ) ) 671 CALL PCLAFCHK( 'No', 'No', M, N, MEM( IPA ), 1, 672 $ 1, DESCA, IASEED, ANORM, FRESID, 673 $ MEM( IPW ) ) 674 ELSE IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 675* 676* Compute residual = ||AP-Q*R|| / (||A||*N*eps) 677* 678 CALL PCGEQRRV( M, N, MEM( IPA ), 1, 1, DESCA, 679 $ MEM( IPTAU ), MEM( IPW ) ) 680 ELSE IF( LSAMEN( 2, FACT, 'TZ' ) ) THEN 681* 682* Compute residual = ||A-T*Z|| / (||A||*N*eps) 683* 684 IF( N.GE.M ) THEN 685 CALL PCTZRZRV( M, N, MEM( IPA ), 1, 1, DESCA, 686 $ MEM( IPTAU ), MEM( IPW ) ) 687 END IF 688 CALL PCLAFCHK( 'No', 'No', M, N, MEM( IPA ), 1, 689 $ 1, DESCA, IASEED, ANORM, FRESID, 690 $ MEM( IPW ) ) 691 END IF 692* 693* Check for memory overwrite 694* 695 CALL PCCHEKPAD( ICTXT, ROUTCHK, MP, NQ, 696 $ MEM( IPA-IPREPAD ), DESCA( LLD_ ), 697 $ IPREPAD, IPOSTPAD, PADVAL ) 698 CALL PCCHEKPAD( ICTXT, ROUTCHK, LTAU, 1, 699 $ MEM( IPTAU-IPREPAD ), LTAU, 700 $ IPREPAD, IPOSTPAD, PADVAL ) 701 CALL PCCHEKPAD( ICTXT, ROUTCHK, WORKSIZ-IPOSTPAD, 702 $ 1, MEM( IPW-IPREPAD ), 703 $ WORKSIZ-IPOSTPAD, IPREPAD, 704 $ IPOSTPAD, PADVAL ) 705* 706 IF( LSAMEN( 2, FACT, 'QP' ) ) THEN 707* 708 CALL PCQPPIV( M, N, MEM( IPA ), 1, 1, DESCA, 709 $ MEM( IPPIV ) ) 710* 711* Check for memory overwrite 712* 713 CALL PCCHEKPAD( ICTXT, 'PCQPPIV', MP, NQ, 714 $ MEM( IPA-IPREPAD ), 715 $ DESCA( LLD_ ), 716 $ IPREPAD, IPOSTPAD, PADVAL ) 717 CALL PCCHEKPAD( ICTXT, 'PCQPPIV', LIPIV, 1, 718 $ MEM( IPPIV-IPREPAD ), LIPIV, 719 $ IPREPAD, IPOSTPAD, PADVAL ) 720* 721 CALL PCLAFCHK( 'No', 'No', M, N, MEM( IPA ), 1, 722 $ 1, DESCA, IASEED, ANORM, FRESID, 723 $ MEM( IPW ) ) 724* 725* Check for memory overwrite 726* 727 CALL PCCHEKPAD( ICTXT, 'PCLAFCHK', MP, NQ, 728 $ MEM( IPA-IPREPAD ), 729 $ DESCA( LLD_ ), 730 $ IPREPAD, IPOSTPAD, PADVAL ) 731 CALL PCCHEKPAD( ICTXT, 'PCLAFCHK', 732 $ WORKSIZ-IPOSTPAD, 1, 733 $ MEM( IPW-IPREPAD ), 734 $ WORKSIZ-IPOSTPAD, IPREPAD, 735 $ IPOSTPAD, PADVAL ) 736 END IF 737* 738* Test residual and detect NaN result 739* 740 IF( LSAMEN( 2, FACT, 'TZ' ) .AND. N.LT.M ) THEN 741 KSKIP = KSKIP + 1 742 PASSED = 'BYPASS' 743 ELSE 744 IF( FRESID.LE.THRESH .AND. 745 $ (FRESID-FRESID).EQ.0.0E+0 ) THEN 746 KPASS = KPASS + 1 747 PASSED = 'PASSED' 748 ELSE 749 KFAIL = KFAIL + 1 750 PASSED = 'FAILED' 751 END IF 752 END IF 753* 754 ELSE 755* 756* Don't perform the checking, only timing 757* 758 KPASS = KPASS + 1 759 FRESID = FRESID - FRESID 760 PASSED = 'BYPASS' 761* 762 END IF 763* 764* Gather maximum of all CPU and WALL clock timings 765* 766 CALL SLCOMBINE( ICTXT, 'All', '>', 'W', 1, 1, WTIME ) 767 CALL SLCOMBINE( ICTXT, 'All', '>', 'C', 1, 1, CTIME ) 768* 769* Print results 770* 771 IF( MYROW.EQ.0 .AND. MYCOL.EQ.0 ) THEN 772* 773 MINMN = MIN( M, N ) 774 MAXMN = MAX( M, N ) 775* 776 IF( LSAMEN( 2, FACT, 'TZ' ) ) THEN 777 IF( M.GE.N ) THEN 778 NOPS = 0.0D+0 779 ELSE 780* 781* 9 ( M^2 N - M^3 ) + 13 M N - M^2 for 782* complete unitary factorization (M <= N). 783* 784 NOPS = 9.0D+0 * ( 785 $ DBLE( N )*( DBLE( M )**2 ) - 786 $ DBLE( M )**3 ) + 787 $ 13.0D+0*DBLE( N )*DBLE( M ) - 788 $ DBLE( M )**2 789 END IF 790* 791 ELSE 792* 793* 8 M N^2 - 8/3 N^2 + 6 M N + 8 N^2 for QR type 794* factorization when M >= N. 795* 796 NOPS = 8.0D+0 * ( DBLE( MINMN )**2 ) * 797 $ ( DBLE( MAXMN )-DBLE( MINMN ) / 3.0D+0 ) + 798 $ ( 6.0D+0 * DBLE( MAXMN ) + 799 $ 8.0D+0 * DBLE( MINMN ) ) * 800 $ DBLE( MINMN ) 801 END IF 802* 803* Print WALL time 804* 805 IF( WTIME( 1 ).GT.0.0D+0 ) THEN 806 TMFLOPS = NOPS / ( WTIME( 1 ) * 1.0D+6 ) 807 ELSE 808 TMFLOPS = 0.0D+0 809 END IF 810 IF( WTIME( 1 ).GE.0.0D+0 ) 811 $ WRITE( NOUT, FMT = 9993 ) 'WALL', M, N, MB, NB, 812 $ NPROW, NPCOL, WTIME( 1 ), TMFLOPS, 813 $ PASSED, FRESID 814* 815* Print CPU time 816* 817 IF( CTIME( 1 ).GT.0.0D+0 ) THEN 818 TMFLOPS = NOPS / ( CTIME( 1 ) * 1.0D+6 ) 819 ELSE 820 TMFLOPS = 0.0D+0 821 END IF 822 IF( CTIME( 1 ).GE.0.0D+0 ) 823 $ WRITE( NOUT, FMT = 9993 ) 'CPU ', M, N, MB, NB, 824 $ NPROW, NPCOL, CTIME( 1 ), TMFLOPS, 825 $ PASSED, FRESID 826* 827 END IF 828* 829 10 CONTINUE 830* 831 20 CONTINUE 832* 833 CALL BLACS_GRIDEXIT( ICTXT ) 834* 835 30 CONTINUE 836* 837 40 CONTINUE 838* 839* Print out ending messages and close output file 840* 841 IF( IAM.EQ.0 ) THEN 842 KTESTS = KPASS + KFAIL + KSKIP 843 WRITE( NOUT, FMT = * ) 844 WRITE( NOUT, FMT = 9992 ) KTESTS 845 IF( CHECK ) THEN 846 WRITE( NOUT, FMT = 9991 ) KPASS 847 WRITE( NOUT, FMT = 9989 ) KFAIL 848 ELSE 849 WRITE( NOUT, FMT = 9990 ) KPASS 850 END IF 851 WRITE( NOUT, FMT = 9988 ) KSKIP 852 WRITE( NOUT, FMT = * ) 853 WRITE( NOUT, FMT = * ) 854 WRITE( NOUT, FMT = 9987 ) 855 IF( NOUT.NE.6 .AND. NOUT.NE.0 ) 856 $ CLOSE ( NOUT ) 857 END IF 858* 859 CALL BLACS_EXIT( 0 ) 860* 861 9999 FORMAT( 'ILLEGAL ', A6, ': ', A5, ' = ', I3, 862 $ '; It should be at least 1' ) 863 9998 FORMAT( 'ILLEGAL GRID: nprow*npcol = ', I4, '. It can be at most', 864 $ I4 ) 865 9997 FORMAT( 'Bad ', A6, ' parameters: going on to next test case.' ) 866 9996 FORMAT( 'Unable to perform ', A, ': need TOTMEM of at least', 867 $ I11 ) 868 9995 FORMAT( 'TIME M N MB NB P Q Fact Time ', 869 $ ' MFLOPS CHECK Residual' ) 870 9994 FORMAT( '---- ------ ------ --- --- ----- ----- --------- ', 871 $ '----------- ------ --------' ) 872 9993 FORMAT( A4, 1X, I6, 1X, I6, 1X, I3, 1X, I3, 1X, I5, 1X, I5, 1X, 873 $ F9.2, 1X, F11.2, 1X, A6, 2X, G8.1 ) 874 9992 FORMAT( 'Finished ', I6, ' tests, with the following results:' ) 875 9991 FORMAT( I5, ' tests completed and passed residual checks.' ) 876 9990 FORMAT( I5, ' tests completed without checking.' ) 877 9989 FORMAT( I5, ' tests completed and failed residual checks.' ) 878 9988 FORMAT( I5, ' tests skipped because of illegal input values.' ) 879 9987 FORMAT( 'END OF TESTS.' ) 880 9986 FORMAT( A ) 881* 882 STOP 883* 884* End of PCQRDRIVER 885* 886 END 887* 888 SUBROUTINE PCQPPIV( M, N, A, IA, JA, DESCA, IPIV ) 889* 890* -- ScaLAPACK routine (version 1.7) -- 891* University of Tennessee, Knoxville, Oak Ridge National Laboratory, 892* and University of California, Berkeley. 893* May 1, 1997 894* 895* .. Scalar Arguments .. 896 INTEGER IA, JA, M, N 897* .. 898* .. Array Arguments .. 899 INTEGER DESCA( * ), IPIV( * ) 900 COMPLEX A( * ) 901* .. 902* 903* Purpose 904* ======= 905* 906* PCQPPIV applies to sub( A ) = A(IA:IA+M-1,JA:JA+N-1) the pivots 907* returned by PCGEQPF in reverse order for checking purposes. 908* 909* Notes 910* ===== 911* 912* Each global data object is described by an associated description 913* vector. This vector stores the information required to establish 914* the mapping between an object element and its corresponding process 915* and memory location. 916* 917* Let A be a generic term for any 2D block cyclicly distributed array. 918* Such a global array has an associated description vector DESCA. 919* In the following comments, the character _ should be read as 920* "of the global array". 921* 922* NOTATION STORED IN EXPLANATION 923* --------------- -------------- -------------------------------------- 924* DTYPE_A(global) DESCA( DTYPE_ )The descriptor type. In this case, 925* DTYPE_A = 1. 926* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating 927* the BLACS process grid A is distribu- 928* ted over. The context itself is glo- 929* bal, but the handle (the integer 930* value) may vary. 931* M_A (global) DESCA( M_ ) The number of rows in the global 932* array A. 933* N_A (global) DESCA( N_ ) The number of columns in the global 934* array A. 935* MB_A (global) DESCA( MB_ ) The blocking factor used to distribute 936* the rows of the array. 937* NB_A (global) DESCA( NB_ ) The blocking factor used to distribute 938* the columns of the array. 939* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first 940* row of the array A is distributed. 941* CSRC_A (global) DESCA( CSRC_ ) The process column over which the 942* first column of the array A is 943* distributed. 944* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local 945* array. LLD_A >= MAX(1,LOCr(M_A)). 946* 947* Let K be the number of rows or columns of a distributed matrix, 948* and assume that its process grid has dimension p x q. 949* LOCr( K ) denotes the number of elements of K that a process 950* would receive if K were distributed over the p processes of its 951* process column. 952* Similarly, LOCc( K ) denotes the number of elements of K that a 953* process would receive if K were distributed over the q processes of 954* its process row. 955* The values of LOCr() and LOCc() may be determined via a call to the 956* ScaLAPACK tool function, NUMROC: 957* LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ), 958* LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ). 959* An upper bound for these quantities may be computed by: 960* LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A 961* LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A 962* 963* Arguments 964* ========= 965* 966* M (global input) INTEGER 967* The number of rows to be operated on, i.e. the number of rows 968* of the distributed submatrix sub( A ). M >= 0. 969* 970* N (global input) INTEGER 971* The number of columns to be operated on, i.e. the number of 972* columns of the distributed submatrix sub( A ). N >= 0. 973* 974* A (local input/local output) COMPLEX pointer into the 975* local memory to an array of dimension (LLD_A, LOCc(JA+N-1)). 976* On entry, the local pieces of the M-by-N distributed matrix 977* sub( A ) which is to be permuted. On exit, the local pieces 978* of the distributed permuted submatrix sub( A ) * Inv( P ). 979* 980* IA (global input) INTEGER 981* The row index in the global array A indicating the first 982* row of sub( A ). 983* 984* JA (global input) INTEGER 985* The column index in the global array A indicating the 986* first column of sub( A ). 987* 988* DESCA (global and local input) INTEGER array of dimension DLEN_. 989* The array descriptor for the distributed matrix A. 990* 991* IPIV (local input) INTEGER array, dimension LOCc(JA+N-1). 992* On exit, if IPIV(I) = K, the local i-th column of sub( A )*P 993* was the global K-th column of sub( A ). IPIV is tied to the 994* distributed matrix A. 995* 996* ===================================================================== 997* 998* .. Parameters .. 999 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_, 1000 $ LLD_, MB_, M_, NB_, N_, RSRC_ 1001 PARAMETER ( BLOCK_CYCLIC_2D = 1, DLEN_ = 9, DTYPE_ = 1, 1002 $ CTXT_ = 2, M_ = 3, N_ = 4, MB_ = 5, NB_ = 6, 1003 $ RSRC_ = 7, CSRC_ = 8, LLD_ = 9 ) 1004* .. 1005* .. Local Scalars .. 1006 INTEGER IACOL, ICOFFA, ICTXT, IITMP, IPVT, IPCOL, 1007 $ IPROW, ITMP, J, JJ, JJA, KK, MYCOL, MYROW, 1008 $ NPCOL, NPROW, NQ 1009* .. 1010* .. External Subroutines .. 1011 EXTERNAL BLACS_GRIDINFO, IGEBR2D, IGEBS2D, IGERV2D, 1012 $ IGESD2D, IGAMN2D, INFOG1L, PCSWAP 1013* .. 1014* .. External Functions .. 1015 INTEGER INDXL2G, NUMROC 1016 EXTERNAL INDXL2G, NUMROC 1017* .. 1018* .. Intrinsic Functions .. 1019 INTRINSIC MIN, MOD 1020* .. 1021* .. Executable Statements .. 1022* 1023* Get grid parameters 1024* 1025 ICTXT = DESCA( CTXT_ ) 1026 CALL BLACS_GRIDINFO( ICTXT, NPROW, NPCOL, MYROW, MYCOL ) 1027 CALL INFOG1L( JA, DESCA( NB_ ), NPCOL, MYCOL, DESCA( CSRC_ ), JJA, 1028 $ IACOL ) 1029 ICOFFA = MOD( JA-1, DESCA( NB_ ) ) 1030 NQ = NUMROC( N+ICOFFA, DESCA( NB_ ), MYCOL, IACOL, NPCOL ) 1031 IF( MYCOL.EQ.IACOL ) 1032 $ NQ = NQ - ICOFFA 1033* 1034 DO 20 J = JA, JA+N-2 1035* 1036 IPVT = JA+N-1 1037 ITMP = JA+N 1038* 1039* Find first the local minimum candidate for pivoting 1040* 1041 CALL INFOG1L( J, DESCA( NB_ ), NPCOL, MYCOL, DESCA( CSRC_ ), 1042 $ JJ, IACOL ) 1043 DO 10 KK = JJ, JJA+NQ-1 1044 IF( IPIV( KK ).LT.IPVT )THEN 1045 IITMP = KK 1046 IPVT = IPIV( KK ) 1047 END IF 1048 10 CONTINUE 1049* 1050* Find the global minimum pivot 1051* 1052 CALL IGAMN2D( ICTXT, 'Rowwise', ' ', 1, 1, IPVT, 1, IPROW, 1053 $ IPCOL, 1, -1, MYCOL ) 1054* 1055* Broadcast the corresponding index to the other process columns 1056* 1057 IF( MYCOL.EQ.IPCOL ) THEN 1058 ITMP = INDXL2G( IITMP, DESCA( NB_ ), MYCOL, DESCA( CSRC_ ), 1059 $ NPCOL ) 1060 CALL IGEBS2D( ICTXT, 'Rowwise', ' ', 1, 1, ITMP, 1 ) 1061 IF( IPCOL.NE.IACOL ) THEN 1062 CALL IGERV2D( ICTXT, 1, 1, IPIV( IITMP ), 1, MYROW, 1063 $ IACOL ) 1064 ELSE 1065 IF( MYCOL.EQ.IACOL ) 1066 $ IPIV( IITMP ) = IPIV( JJ ) 1067 END IF 1068 ELSE 1069 CALL IGEBR2D( ICTXT, 'Rowwise', ' ', 1, 1, ITMP, 1, MYROW, 1070 $ IPCOL ) 1071 IF( MYCOL.EQ.IACOL .AND. IPCOL.NE.IACOL ) 1072 $ CALL IGESD2D( ICTXT, 1, 1, IPIV( JJ ), 1, MYROW, IPCOL ) 1073 END IF 1074* 1075* Swap the columns of A 1076* 1077 CALL PCSWAP( M, A, IA, ITMP, DESCA, 1, A, IA, J, DESCA, 1 ) 1078* 1079 20 CONTINUE 1080* 1081* End of PCQPPIV 1082* 1083 END 1084