//=================================================================C // DRM公式によるπ計算プログラム ( ATAN 版 , C ) // 2003/02/15 Yasunori Ushiro ( 後 保範 ) //=================================================================C // Pi Value Computation Program by Atan with DRM Method C // Pi = M1*Atan(1/X1) + M2*Atan(1/X2) + M3*Atan(1/X3) C // Atan(1/X) = 1/X - 1/(3*X**3) + 1/(5*X**5) - 1/(7*X**7) +.. C // Using FMT for High-Precision Multiplication C //-----------------------------------------------------------------C // Written by Yasunori Ushiro, 2003/02/15 C // (Hitachi Ltd. and Waseda University) C //=================================================================C #include #include #include // Global Value #define ND 270001 #define L10V 5 double A[ND], B[ND], E[ND], X[ND], Y[ND], PI[ND] ; double W[4*ND], WK[6*ND] ; int NPI[ND+3] ; FILE *FT1 ; double BASE=1.0e5 ; int NORID=-1 ; double VS ; // Number of Current Term in Atan int INO[5] = {2,3,3,4,4}; // Number of Atan for Pi // Constant of Atan formulation ; (M,X) int IDATA[5][4][2] = { 16, 5, -4, 239, 0,0, 0,0, // ID=1 48,18, 32, 57, -20,239, 0,0, // ID=2 32,10, -4, 239, -16,515, 0,0, // ID=3 176,57, 28,239, -48,682, 96,12943, // ID=4 48,49, 128,57, -20,239, 48,110443 } ; // ID=5 // Function Definition in the Program double DRMATN(int IX, int NN, int MLP) ; double DRMAT2(int IX, int NN) ; double DRMAT1(int N, int M) ; // Function Definition in the PIOUT int XPIIN(int MID, int NX, int *N, int *NC, int *NN, int *L10, int *LBASE, int *ISB) ; void XPIOUT(double *PI, int NC, int L10, int LBASE, int ISB, int T1, double ERR) ; // Function Definition in the MULT double MTCFM(double *A, double *B, double *C, int N, double *WK) ; double FMTCFM(double *A, double *B, double *C,int N, double *WK) ; double FDIV(double *A, double *B, double *C, int N, double *FW, double *WK) ; double FCSQRT(int IA, double *C, int N, int ID, double *FW, double *WK) ; double FSQRT(double *A, double *C, int N, int ID, double *FW, double *WK) ; void CFMT4(double *A, int N, int ID, double *WK) ; void CFMT(double *A, int N, int ID, double *WK) ; void CFMTTB(int N2, int ID, double *TBL) ; void CFMTSB(double *A, int M, int L, double *TB, double *WK) ; // Function definition in the ADD void ADD(double *A, double *B, double *C, int N, int ID) ; void FADD(double *A, double *B, double *C, int N, int ID) ; void FCADD(int IA, double *B, double *C, int N, int ID) ; void FCMT(double *A, int IB, double *C, int N) ; void FCDIV(double *A, int IB, double *C, int N) ; void NORM(double *A, int N) ; void FNORM(double *A, int N) ; void FLOAT(double *A, int N) ; void COPY(double *A, double *C, int N) ; void SETV(int IV, double *C, int N) ; void FVSET(int IV, double *C, int N) ; double fmax(double A1, double A2) ; double fmax3(double A1, double A2, double A3) ; int MIN(int N1, int N2) ; int LOG2F(int N) ; int power(int N1, int N2) ; //=================================================================C // Main Program C //-----------------------------------------------------------------C void main () { int k, ID, N, NN, NC, L10, LBASE, ISB, MLP ; int IX, IB, T1, T2, T3 ; double XX, ALG, ERR, ER1, TT ; char FNAME[] = "AtanDRM?.txt" ; char FID[] = "12345" ; // Parameter Read and Set ID = XPIIN (2, ND, &N, &NC, &NN, &L10, &LBASE, &ISB) ; T1 = clock() ; // Start Time ERR = 0.0 ; SETV(0, PI, N+1) ; // PI=float(0) // Main Loop for (k=0; k