C*********************************************************************C C** Floating-Point Multiple-Precision Arithmetic Test Program **C C** for Massively Parallel Processor HITACHI SR2201 **C C** **C C** All rights reserved, Copyright (C) 1998, Daisuke TAKAHASHI **C C** **C C*********************************************************************C c/*- c * Copyright (c) 1998 Information Promotion Agency of Japan and D.Takahashi c * All rights reserved. c * c * Redistribution and non-commercial use in source and binary forms, c * with or without modification, are permitted provided that the following c * conditions are met: c * c * 1. Redistributions of source code must retain the above copyright c * notice, this list of conditions and the following disclaimer. c * c * 2. Redistributions in binary form must reproduce the above copyright c * notice, this list of conditions and the following disclaimer in the c * documentation and/or other materials provided with the distribution. c * c * 3. All advertising materials mentioning features or use of this c * software must display the following acknowledgement: c * c * This product includes software developed by Information Promotion c * Agency of Japan and D.Takahashi c * c * 4. The name "Information Promotion Agency" or "D.Takahashi" should not be c * used to endorse or promote products derived from this software c * without specific prior written permission. c * c * 5. Any use of the source code or the binary in a commercial product, c * whether may it be the origial representation or in some modified form, c * is not permitted without specific prior written permission. c * c * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND c * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE c * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE c * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE c * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL c * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS c * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) c * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT c * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY c * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF c * SUCH DAMAGE. c */ IMPLICIT REAL*8 (A-H,O-Z) INCLUDE 'mpif.h' PARAMETER (NDA=1048576,NDB=NDA/4) DIMENSION IA(NDB+2),IB(NDB+2),IS(NDB+2),IT(NDB+2) DIMENSION IV(NDB+2),IW(NDB+2),IX(NDB+2),IY(NDB+2),IZ(NDB+2) COMMON /IPROC/ME/NPROC/NPU REAL*8 N C CALL MPI_INIT(IERR) CALL MPI_COMM_RANK(MPI_COMM_WORLD,ME,IERR) CALL MPI_COMM_SIZE(MPI_COMM_WORLD,NPU,IERR) C IF (ME .EQ. 0) THEN READ(5,*) M END IF CALL MPI_BCAST(M,1,MPI_INTEGER,0,MPI_COMM_WORLD,IERR) N=2.0D0**(M-3) C CALL MPINIT(ME,NPU) C CALL MPI_BARRIER(MPI_COMM_WORLD,IERR) IF (ME .EQ. 0) THEN TIME1=MPI_WTIME() END IF C CALL GAUSS(IA,IB,IS,IT,IV,IW,IX,IY,IZ,N) IF (ME .EQ. 0) THEN WRITE(6,*) ' GAUSS-LEGENDRE ALGORITHM' END IF CALL OUTPUT(IB(3),N) C CALL BORWEIN(IA,IB,IS,IV,IW,IX,IY,IZ,N) IF (ME .EQ. 0) THEN WRITE(6,*) ' BORWEIN 4-TH CONVERGENT ALGORITHM' END IF CALL OUTPUT(IB(3),N) C IF (ME .EQ. 0) THEN TIME2=MPI_WTIME() TIME0=TIME2-TIME1 WRITE(6,*) ' TIME=',TIME0 END IF C CALL MPI_FINALIZE(IERR) STOP END SUBROUTINE GAUSS(IA,IB,IS,IT,IV,IW,IX,IY,IZ,N) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IA(*),IB(*),IS(*),IT(*),IV(*),IW(*),IX(*),IY(*),IZ(*) COMMON /IPROC/ME/NPROC/NPU REAL*8 N C M=IDINT(DLOG10(N*8.0D0)/DLOG10(2.0D0)) C CALL MPCLEAR(IA,N,IERR) CALL MPCLEAR(IB,N,IERR) CALL MPCLEAR(IS,N,IERR) CALL MPCLEAR(IT,N,IERR) CALL MPCLEAR(IV,N,IERR) CALL MPCLEAR(IW,N,IERR) CALL MPCLEAR(IX,N,IERR) CALL MPCLEAR(IY,N,IERR) CALL MPCLEAR(IZ,N,IERR) C IA(1)=1 IA(2)=1 IF (ME .EQ. 0) THEN IA(3)=10000000 END IF IW(1)=1 IW(2)=0 IF (ME .EQ. 0) THEN IW(3)=50000000 END IF CALL MPSQRT(IW,IB,IV,IZ,N,IERR) IT(1)=1 IT(2)=0 IF (ME .EQ. 0) THEN IT(3)=25000000 END IF IS(1)=1 IS(2)=0 IF (ME .EQ. 0) THEN IS(3)=50000000 END IF IX(1)=1 IX(2)=1 IF (ME .EQ. 0) THEN IX(3)=10000000 END IF C DO 10 I=1,M CALL MPMOVE(IA,IY,N,IERR) CALL MPADD(IA,IB,IA,N,IERR) CALL MPMUL(IA,IS,IA,N,IERR) CALL MPMUL(IB,IY,IW,N,IERR) CALL MPSQRT(IW,IB,IV,IZ,N,IERR) CALL MPSUB(IY,IA,IW,N,IERR) CALL MPMUL(IW,IW,IW,N,IERR) CALL MPMUL(IX,IW,IW,N,IERR) CALL MPSUB(IT,IW,IT,N,IERR) CALL MPADD(IX,IX,IX,N,IERR) 10 CONTINUE CALL MPADD(IA,IB,IA,N,IERR) CALL MPMUL(IA,IA,IA,N,IERR) IS(1)=1 IS(2)=0 IF (ME .EQ. 0) THEN IS(3)=25000000 END IF CALL MPMUL(IA,IS,IA,N,IERR) CALL MPDIV(IA,IT,IB,IZ,N,IERR) RETURN END SUBROUTINE BORWEIN(IA,IB,IS,IV,IW,IX,IY,IZ,N) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IA(*),IB(*),IS(*),IV(*),IW(*),IX(*),IY(*),IZ(*) COMMON /IPROC/ME/NPROC/NPU REAL*8 N C M=IDINT(DLOG10(N*8.0D0)/DLOG10(2.0D0)) C CALL MPCLEAR(IA,N,IERR) CALL MPCLEAR(IB,N,IERR) CALL MPCLEAR(IS,N,IERR) CALL MPCLEAR(IV,N,IERR) CALL MPCLEAR(IW,N,IERR) CALL MPCLEAR(IX,N,IERR) CALL MPCLEAR(IY,N,IERR) CALL MPCLEAR(IZ,N,IERR) C IA(1)=1 IA(2)=1 IF (ME .EQ. 0) THEN IA(3)=60000000 END IF IW(1)=1 IW(2)=1 IF (ME .EQ. 0) THEN IW(3)=20000000 END IF CALL MPSQRT(IW,IB,IV,IZ,N,IERR) CALL MPADD(IB,IB,IB,N,IERR) CALL MPADD(IB,IB,IB,N,IERR) CALL MPSUB(IA,IB,IA,N,IERR) CALL MPADD(IB,IB,IY,N,IERR) CALL MPADD(IB,IY,IY,N,IERR) IW(1)=1 IW(2)=2 IF (ME .EQ. 0) THEN IW(3)=17000000 END IF CALL MPSUB(IW,IY,IY,N,IERR) IX(1)=1 IX(2)=1 IF (ME .EQ. 0) THEN IX(3)=20000000 END IF IS(1)=1 IS(2)=1 IF (ME .EQ. 0) THEN IS(3)=10000000 END IF C DO 10 I=1,M/2 CALL MPSUB(IS,IY,IY,N,IERR) CALL MPSQRT(IY,IB,IV,IZ,N,IERR) CALL MPSQRT(IB,IW,IV,IZ,N,IERR) CALL MPSUB(IS,IW,IB,N,IERR) CALL MPADD(IS,IW,IW,N,IERR) CALL MPDIV(IB,IW,IY,IZ,N,IERR) CALL MPMUL(IY,IY,IB,N,IERR) CALL MPADD(IY,IY,IW,N,IERR) CALL MPMUL(IB,IB,IY,N,IERR) CALL MPADD(IW,IB,IW,N,IERR) CALL MPADD(IS,IW,IW,N,IERR) CALL MPMUL(IW,IW,IW,N,IERR) CALL MPMUL(IA,IW,IA,N,IERR) CALL MPADD(IB,IB,IB,N,IERR) CALL MPADD(IS,IB,IB,N,IERR) CALL MPADD(IB,IY,IB,N,IERR) CALL MPSUB(IW,IB,IW,N,IERR) CALL MPMUL(IX,IW,IW,N,IERR) CALL MPSUB(IA,IW,IA,N,IERR) CALL MPADD(IX,IX,IX,N,IERR) CALL MPADD(IX,IX,IX,N,IERR) 10 CONTINUE CALL MPDIV(IS,IA,IB,IZ,N,IERR) RETURN END