|
PUBLIC NUMERIC BIAS,KETA,SIGN,EPS
LET BIAS=0
LET KETA=2500 !'桁数=KETA*4
LET SIGN=-BIAS-1
LET EPS=0
LET KETA=KETA+EPS
DIM X(-BIAS-1 TO KETA),S(-BIAS-1 TO KETA)
!'π/4=2805*ATN(1/5257)-398*ATN(1/9466)+1950*ATN(1/12943)+1850*ATN(1/34208)+2021*ATN(1/44179)+2097*ATN(1/85353)+1484*ATN(1/114669)+1389*ATN(1/330182)+808*ATN(1/485298)
PRINT(KETA-EPS)*4;"桁の計算開始 ";TIME$
LET T=TIME
PRINT "2805*ATN(1/5257) ";TIME$
CALL SATN(X,5257)
CALL SMUL(X,2805*4)
CALL LCOPY(S,X)
PRINT "398*ATN(1/9466) ";TIME$
CALL SATN(X,9466)
CALL SMUL(X,398*4)
CALL LSUB2(S,X)
PRINT "1950*ATN(1/12943) ";TIME$
CALL SATN(X,12943)
CALL SMUL(X,1950*4)
CALL LADD2(S,X)
PRINT "1850*ATN(1/34208) ";TIME$
CALL SATN(X,34208)
CALL SMUL(X,1850*4)
CALL LADD2(S,X)
PRINT "2021*ATN(1/44179) ";TIME$
CALL SATN(X,44179)
CALL SMUL(X,2021*4)
CALL LADD2(S,X)
PRINT "2097*ATN(1/85353) ";TIME$
CALL SATN(X,85353)
CALL SMUL(X,2097*4)
CALL LADD2(S,X)
PRINT "1484*ATN(1/114669) ";TIME$
CALL SATN(X,114669)
CALL SMUL(X,1484*4)
CALL LADD2(S,X)
PRINT "1389*ATN(1/330182) ";TIME$
CALL SATN(X,330182)
CALL SMUL(X,1389*4)
CALL LADD2(S,X)
PRINT "808*ATN(1/485298) ";TIME$
CALL SATN(X,485298)
CALL SMUL(X,808*4)
CALL LADD2(S,X)
PRINT "計算終了 ";TIME$;TIME-T
CALL DISPLAY(S,"")
END
EXTERNAL SUB SATN(S(),XX)!'ATN(1/X)
!'ATN(1/X)=X/(X^2+1)*(1+(2)/(3)*1/(X^2+1)+(2*4)/(3*5)*(1/(X^2+1)^2+(2*4*6)/(3*5*7)*(1/(X^2+1))^3+...
DIM SS(-BIAS-1 TO KETA),Y(-BIAS-1 TO KETA)
LET Y(0)=1
LET Y(SIGN)=1
CALL SDIV(Y,XX*XX+1)
CALL SMUL(Y,XX)
CALL LCOPY(S,Y)
DO
LET N=N+1
CALL LCOPY(SS,S)
CALL SMUL(Y,2*N)
CALL SDIV(Y,2*N+1)
CALL SDIV(Y,XX*XX+1)
CALL LADD2(S,Y)
LOOP UNTIL EQUAL(SS,S)<>0
END SUB
EXTERNAL SUB DISPLAY(X(),N$)
IF N$="" THEN
OPEN #1:TEXTWINDOW1
ELSE
OPEN #1:NAME N$
END IF
ERASE #1
FOR K=-BIAS TO 0
IF X(K)<>0 THEN EXIT FOR
NEXT K
IF X(SIGN)=-1 THEN PRINT #1:"- ";
IF K>=0 THEN
LET K=0
PRINT #1:STR$(X(0));"."
ELSE
PRINT #1:STR$(X(K));
FOR I=K+1 TO 0
LET A$=A$&RIGHT$("000"&STR$(X(I)),4)
IF LEN(A$)=100 THEN
PRINT #1:A$
LET A$=""
END IF
NEXT I
IF LEN(A$)>0 THEN
PRINT #1:A$;"."
LET A$=""
END IF
END IF
LET S=0
FOR I=1 TO KETA-EPS
LET A$=A$&RIGHT$("000"&STR$(X(I)),4)
IF LEN(A$)=100 THEN
LET S=S+100
FOR J=1 TO 10
PRINT #1:LEFT$(A$,10);" ";
IF J=5 THEN PRINT #1:" ";
LET A$=RIGHT$(A$,LEN(A$)-10)
NEXT J
PRINT #1:":";S
LET A$=""
IF MOD(S,1000)=0 THEN PRINT #1
END IF
NEXT I
IF LEN(A$)>0 THEN
LET S=S+LEN(A$)
LET A$=A$&REPEAT$(" ",10)
FOR J=1 TO 9
PRINT #1:RTRIM$(LEFT$(A$,10));" ";
IF J=5 THEN PRINT #1:" ";
LET A$=RIGHT$(A$,LEN(A$)-10)
IF RTRIM$(A$)="" THEN EXIT FOR
NEXT J
IF A$<>"" THEN PRINT #1:A$
END IF
CLOSE #1
END SUB
EXTERNAL SUB LCLR(X())
MAT X=ZER
LET X(SIGN)=1
END SUB
EXTERNAL SUB LCOPY(X(),Y())
MAT X=Y
END SUB
EXTERNAL SUB LADD(A(),B(),C())
LET SIGNA=A(SIGN)
LET SIGNB=B(SIGN)
IF SIGNA=1 AND SIGNB=-1 THEN
LET B(SIGN)=1
CALL LSUB(A,B,C)
LET B(SIGN)=-1
EXIT SUB
ELSEIF SIGNA=-1 AND SIGNB=1 THEN
LET A(SIGN)=1
CALL LSUB(B,A,C)
LET A(SIGN)=-1
EXIT SUB
END IF
MAT C=A+B
FOR I=KETA TO-BIAS+1 STEP-1
IF C(I)>=10000 THEN
LET C(I)=MOD(C(I),10000)
LET C(I-1)=C(I-1)+1
END IF
NEXT I
IF C(-BIAS)>=10000 THEN PRINT "OVER FLOW in LADD"
IF SIGNA=-1 AND SIGNB=-1 THEN LET C(SIGN)=-1 ELSE LET C(SIGN)=1
END SUB
EXTERNAL SUB LADD2(A(),B())
DIM C(-BIAS-1 TO KETA)
CALL LADD(A,B,C)
CALL LCOPY(A,C)
END SUB
EXTERNAL SUB LSUB(A(),B(),C())
LET SIGNA=A(SIGN)
LET SIGNB=B(SIGN)
LET A(SIGN)=1
LET B(SIGN)=1
IF SIGNA*SIGNB=-1 THEN
CALL LADD(A,B,C)
LET C(SIGN)=SIGNA
LET A(SIGN)=SIGNA
LET B(SIGN)=SIGNB
EXIT SUB
END IF
LET GR=GREAT(A,B)
IF SIGNA=1 AND SIGNB=1 THEN
IF GR<>0 THEN
MAT C=A-B
LET C(SIGN)=1
ELSE
MAT C=B-A
LET C(SIGN)=-1
END IF
ELSE
IF GR<>0 THEN
MAT C=B-A
LET C(SIGN)=1
ELSE
MAT C=A-B
LET C(SIGN)=-1
END IF
END IF
FOR I=KETA TO-BIAS+1 STEP-1
IF C(I)<0 THEN
LET C(I)=C(I)+10000
LET C(I-1)=C(I-1)-1
END IF
NEXT I
LET A(SIGN)=SIGNA
LET B(SIGN)=SIGNB
END SUB
EXTERNAL SUB LSUB2(A(),B())
DIM C(-BIAS-1 TO KETA)
CALL LSUB(A,B,C)
CALL LCOPY(A,C)
END SUB
EXTERNAL SUB SMUL(A(),XA)
LET SIGNA=A(SIGN)
LET SG=SGN(XA)
LET XA=ABS(XA)
MAT A=XA*A
FOR I=KETA TO-BIAS+1 STEP-1
IF A(I)>=10000 THEN
LET R=INT(A(I)/10000)
LET A(I)=MOD(A(I),10000)
LET A(I-1)=A(I-1)+R
END IF
NEXT I
IF A(-BIAS)>=10000 THEN PRINT "OVER FLOW in SMUL"
IF SIGNA*SG=-1 THEN LET A(SIGN)=-1 ELSE LET A(SIGN)=1
END SUB
EXTERNAL SUB SDIV(A(),XA)
LET SIGNA=A(SIGN)
LET SG=SGN(XA)
LET XA=ABS(XA)
FOR I=-BIAS TO KETA-1
LET R=A(I)-INT(A(I)/XA)*XA
LET A(I)=INT(A(I)/XA)
LET A(I+1)=A(I+1)+R*10000
NEXT I
LET A(KETA)=INT(A(KETA)/XA)
IF SIGNA*SG=-1 THEN LET A(SIGN)=-1
END SUB
EXTERNAL FUNCTION GREAT(A(),B())
LET SIGNA=A(SIGN)
LET SIGNB=B(SIGN)
IF SIGNA=-1 AND SIGNB=1 THEN
LET GREAT=0
EXIT FUNCTION
END IF
IF SIGNA=1 AND SIGNB=-1 THEN
LET GREAT=-1
EXIT FUNCTION
END IF
FOR I=-BIAS TO KETA
IF SIGNA=-1 AND SIGNB=-1 THEN
IF A(I)<B(I)THEN
LET GREAT=-1
EXIT FUNCTION
END IF
IF A(I)>B(I)THEN
LET GREAT=0
EXIT FUNCTION
END IF
ELSE
IF A(I)>B(I)THEN
LET GREAT=-1
EXIT FUNCTION
END IF
IF A(I)<B(I)THEN
LET GREAT=0
EXIT FUNCTION
END IF
END IF
NEXT I
LET GREAT=0
END FUNCTION
EXTERNAL FUNCTION EQUAL(A(),B())
FOR I=-BIAS-1 TO KETA-EPS
IF A(I)<>B(I)THEN
LET EQUAL=0
EXIT FUNCTION
END IF
NEXT I
LET EQUAL=-1
END FUNCTION
|
|