平方根 √A のニュートン法による漸化式
Xn+1 ：= (Xn + A / Xn) / 2
これを希望する桁になるまで繰り返し計算する。
繰り返す回数は、１００万桁≒２の２０乗 ２０回。
１億桁なら、≒２の２６．５乗→２７回。
A / Xnは、多倍長÷多倍長の形なので計算が複雑になる。
1 / √A のニュートン法による漸化式
Xn+1： = Xn * (3 – A * Xn2) / 2 <乗算回数 5回 係数含>
但し、Ａ＝２のときには Xn * (１．５ – Xn2) <乗算回数 ２回>
これだと回数は多いが 多倍長÷多倍長はなくなる。
※２の割算があるが、０．５倍と考えれば、すべて乗算になる。
√A＝ A ＊ （1 / √A） であるから、
1 / √Aがわかれば、√Aが 乗算でもとめられる。
√２ １００万桁.txt (約100kＢ）：但し、抜粋。
計算時間 １１時間４０分５秒 (10進ＢＡＳＩＣによる１００桁区切り計算にて)
ＣＰＵ：Celelon ９００ＭＨz
1.
4142135623 7309504880 1688724209 6980785696 7187537694 8073176679 7379907324 7846210703 8850387534 3276415727
3501384623 0912297024 9248360558 5073721264 4121497099 9358314132 2266592750 5592755799 9505011527 8206057147
0109559971 6059702745 3459686201 4728517418 6408891986 0955232923 0484308714 3214508397 6260362799 5251407989
6872533965 4633180882 9640620615 2583523950 5474575028 7759961729 8355752203 3753185701 1354374603 4084988471
6038689997 0699004815 0305440277 9031645424 7823068492 9369186215 8057846311 1596668713 0130156185 6898723723
5288509264 8612494977 1542183342 0428568606 0146824720 7714358548 7415565706 9677653720 2264854470 1585880162
0758474922 6572260020 8558446652 1458398893 9443709265 9180031138 8246468157 0826301005 9485870400 3186480342
1948972782 9064104507 2636881313 7398552561 1732204024 5091227700 2269411275 7362728049 5738108967 5040183698
6836845072 5799364729 0607629969 4138047565 4823728997 1803268024 7442062926 9124859052 1810044598 4215059112
0249441341 7285314781 0580360337 1077309182 8693147101 7111168391 6581726889 4197587165 8215212822 9518488472
1000 桁目
! √２の計算 百桁区切り 十進ＢＡＳＩＣ
! １／√２ ： Xn+1 = Xn * (3 - A * Xn^2) / 2
!
OPTION ARITHMETIC DECIMAL_HIGH
DO
DO
input PROMPT "桁数は？ (1000単位で)":KTA
loop while KTA<2000
loop until KTA/1000=int(KTA/1000)
print "計算開始 ";date$;" ";time$
print "桁数=";KTA;"桁"
! 100桁区切り
LET BLOCK=int(KTA*1.1/100) ! マージン10％とする
print "BLOCK=";BLOCK
option base 0
dim XN0(BLOCK),XN1(BLOCK),A(BLOCK),XN02(BLOCK)
dim AXN02(BLOCK),AXN02m3(BLOCK),AXN02m3X(BLOCK)
dim ZPT5(BLOCK)
LET KS=0
gosub 20000 ! まず初期化
10 !
DO
LET KS=KS+1
loop while KTA > 2^(KS-1)
print "KS=";KS+1
for K=0 to KS+1
print "計算";K;"回目";TIME$
gosub 100
next K
!最終計算 √２＝２×（１／√２）
FOR I=0 to BLOCK
LET XN0(I)=2*XN0(I)
next I
20 ! 桁処理 --桁上げ処理
Do
LET OBF=0
for I=BLOCK TO 1 step -1
LET XN0(I-1)=XN0(I-1)+int(XN0(I)/10^100)
LET XN0(I)=XN0(I)-10^100*int(XN0(I)/10^100)
if XN0(I)>10^100 then LET OBF=1
NEXT I
loop while OBF=1
!印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷
print "終了！"
print "計算終了 ";date$;" ";time$
print str$(XN0(0));"."
for I=1 TO BLOCK
LET p$=repeat$("0",100-len(str$(XN0(I))))&str$(XN0(I))
print MID$(p$,1,10);" ";MID$(p$,11,10);" ";MID$(p$,21,10);" ";MID$(p$,31,10);" ";MID$(p$,41,10);" ";
print MID$(p$,51,10);" ";MID$(p$,61,10);" ";MID$(p$,71,10);" ";MID$(p$,81,10);" ";MID$(p$,91,10);" "
if int(I/10)=I/10 then print I*100;"桁目"
next I
!印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷印刷
stop
100 !
! Xn+1 = Xn * (3 - A * Xn^2) / 2
! ①Xn^2の計算
for I=0 TO BLOCK
for J=0 TO BLOCK-I
LET XN02(I+J)=XN0(I)*XN0(J)+XN02(I+J)
next J
next I
! ②A * Xn^2の計算
for I=0 TO BLOCK
for J=0 TO BLOCK-I
LET AXN02(I+J)=A(I)*XN02(J)+AXN02(I+J)
next J
next I
! ③３－A * Xn^2の計算
LET AXN02m3(0)=3
for I=0 TO BLOCK
LET AXN02m3(I)=AXN02m3(I)-AXN02(I)
next I
! ④Xn *（３－A * Xn^2）の計算
for I=0 TO BLOCK
for J=0 TO BLOCK-I
LET AXN02m3X(I+J)=XN0(I)*AXN02m3(J)+AXN02m3X(I+J)
next J
next I
! ⑤Xn *（３－A * Xn^2）／２の計算
for I=0 TO BLOCK
for J=0 TO BLOCK-I
LET XN1(I+J)=AXN02m3X(I)*ZPT5(J)+XN1(I+J)
next J
next I
1000 ! 桁処理１ --負記号の削除
DO
LET MNSF=0
if XN1(1)<0 then
LET XN1(0)=XN1(0)+int(XN1(1)/10^100)-1
LET XN1(1)=XN1(1)-10^100*int(XN1(1)/10^100)+10^100
end if
for I=BLOCK TO 1 step -1
if XN1(I)<0 then
LET MNSF=1
LET XN1(I-1)=XN1(I-1)+int(XN1(I)/10^100)-1
LET XN1(I)=XN1(I)-10^100*int(XN1(I)/10^100)+10^100
end if
NEXT I
LOOP while MNSF=1
2000 ! 桁処理２ --桁上げ処理
DO
LET OBF=0
for I=BLOCK TO 1 step -1
LET XN1(I-1)=XN1(I-1)+int(XN1(I)/10^100)
LET XN1(I)=XN1(I)-10^100*int(XN1(I)/10^100)
if XN1(I)>10^100 then
LET OBF=1
end if
NEXT I
loop while OBF=1
3000 ! XN1 => XN0 へ入替
for I=0 to BLOCK
LET XN0(I)=XN1(I)
NEXT I
!２回目以降の初期化
for I=0 To BLOCK
LET A(I)=0
! XN0(I)=0 ２回目の初期化では『0』にしない。
LET XN1(I)=0
LET XN02(I)=0
LET AXN02(I)=0
LET AXN02m3(I)=0
LET AXN02m3X(I)=0
LET ZPT5(I)=0
next I
LET A(0)=2
LET ZPT5(1)=5*10^99 ! 0.5
return
20000 !初期設定
for I=0 To BLOCK
LET A(I)=0
LET XN0(I)=0
LET XN1(I)=0
LET XN02(I)=0
LET AXN02(I)=0
LET AXN02m3(I)=0
LET AXN02m3X(I)=0
LET ZPT5(I)=0
next I
LET XN0(0)=1
LET A(0)=2
LET ZPT5(1)=5*10^99 ! 0.5
return
END