|
投稿者：しばっち
投稿日：2016年10月23日(日)17時45分27秒
|
|
|
|
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
|
|
|
|
投稿者：たろさ
投稿日：2016年10月25日(火)20時07分53秒
|
|
|
|
投稿者：しばっち
投稿日：2016年10月25日(火)20時57分17秒
|
|
|
|
たろささんへのお返事です。
> 何時もお世話になります。最初に、
> 十進 OPTION ARITHMETIC DECIMAL (Full BASIC互換)
> 計算して、高速だ!!と感動しました。
>
> 次に
> 有理数 OPTION ARITHMETIC RATIONAL
> 計算して、高速だ!!と感動しました。
多倍長計算ルーチンは通常「2進モード」で実行します。
「ヘルプ(H)」の下あたりにある「2」のボタンを押してから実行してください。
※起動時は「10」のボタンが押されています。10進モードから変更してください
> 最終桁の微調整を教えて頂けないでしょうか?
LET KETA=2500 <-----ここの数字で桁数を調整してください。ここでは2500*4で1万桁です
LET SIGN=-BIAS-1
LET EPS=0 <---------ここを増やしてください。1上げる毎に4桁分増えます。
LET KETA=KETA+EPS
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
|
|
|
|
投稿者：たろさ
投稿日：2016年10月25日(火)21時38分33秒
|
|
|
|
> No.4156[元記事へ]
しばっちさんへのお返事です。
> たろささんへのお返事です。
>
> LET KETA=2500 <-----ここの数字で桁数を調整してください。ここでは2500*4で1万桁です
> LET SIGN=-BIAS-1
> LET EPS=0 <---------ここを増やしてください。1上げる毎に4桁分増えます。
> LET KETA=KETA+EPS
>
> 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
>
とてもよくわかりました。精度確認完了
実は、最初に二進モードでエラーが出ました。
オプションを外して、ボタンで実行すると、エラーが出ません。
■ パソコン環境
ウィンドウズ : Microsoft Windows 8.1
サービスパック : なし 自作機
システムの種類 : 32 ビット
プロセッサー : AMD Athlon(tm) 64 Processor 3800+
周波数 : 2.40 GHz
メモリー : 3.00 GB
Re: データの並べ替えについて
#4151
わたしは、「-,24 ,25, , 26,,,,,=27M,28\, 29,,.,」このような数列は読めないと
諦めていたので、感動しました。BASIC Acc 二進モードで円周率を書き出して
0～9までの、個数を数えるprogramも出来そうです。
まだ、BASIC Accで書き出しに失敗してます。
http://blogs.yahoo.co.jp/donald_stinger
|
|
|
|
投稿者：Yoshi
投稿日：2016年12月 6日(火)11時46分8秒
|
|
|
|
> No.4138[元記事へ]
しばっちさんへのお返事です。
失礼します。
このプログラムを動かしてみたのですが、
結果の数値が正確な円周率と違うようです。
3.
1415926535 8979323846 2643383279 5028841971 6945427326
ネットで見た他の方々の表示
3.
1415926535 8979323846 2643383279 5028841971 6939937510
> 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
>
|
|
|
|
投稿者：たろさ
投稿日：2016年12月 6日(火)20時29分16秒
|
|
|
戻る
