アルゴリズム論特論(塩田)第10回 (2) 多倍長整数の平方根
多倍長整数の平方根
多倍長整数の平方根は倍精度実数の sqrt では求められませんので、
「開平法」と呼ばれる方法を用いたり、
次のニュートン法を用いて計算します。
Algorithm 5(ニュートン法)
- 入力:自然数 $n$
- 出力:$n$ の平方根の切り捨て
def mysqrt(n):
if n < = 1:
return n
s = n
while s * s > n:
s = (s + n//s) // 2
return s
数値解析で実数のニュートン法を習いましたね。方程式 $f(x)=0$ の近似解は数列
$\displaystyle{x_{k+1}=x_k - \frac{x_k}{f'(x_k)}}$
を作って求めました。$\sqrt{n}$ は $f(x)=x^2-n=0$ の解で、そのとき上の数列は
$\displaystyle{x_{k+1}=x_k + \frac{x_k^2-n}{2x_k}=\frac{1}{2}\left(x_k+\frac{n}{x_k}\right)}$
となります。除算を整数商に代えたものが
Alg.5 です。
Ex.6
n = 9956894514697920852 ( 64 bit )
Step 0 : s = 9956894514697920852
Step 1 : s = 4978447257348960426
Step 2 : s = 2489223628674480214
Step 3 : s = 1244611814337240108
Step 4 : s = 622305907168620057
Step 5 : s = 311152953584310036
Step 6 : s = 155576476792155033
Step 7 : s = 77788238396077548
Step 8 : s = 38894119198038837
Step 9 : s = 19447059599019546
Step 10 : s = 9723529799510028
Step 11 : s = 4861764899755525
Step 12 : s = 2430882449878786
Step 13 : s = 1215441224941440
Step 14 : s = 607720612474815
Step 15 : s = 303860306245599
Step 16 : s = 151930153139183
Step 17 : s = 75965076602359
Step 18 : s = 37982538366715
Step 19 : s = 18991269314429
Step 20 : s = 9495634919358
Step 21 : s = 4747817983966
Step 22 : s = 2373910040558
Step 23 : s = 1186957117429
Step 24 : s = 593482753008
Step 25 : s = 296749765032
Step 26 : s = 148391659099
Step 27 : s = 74229378923
Step 28 : s = 37181757882
Step 29 : s = 18724773827
Step 30 : s = 9628261807
Step 31 : s = 5331196948
Step 32 : s = 3599431402
Step 33 : s = 3182836172
Step 34 : s = 3155572503
Step 35 : s = 3155454726
Step 36 : s = 3155454723
answer = 3155454723
Prop.7 Alg.5 の計算量は $O(\log^3 n)$.
証明 終盤までは $s$ がほとんど半分半分になっていきますので (
Ex.6 参照 )、ループ回数は $O(\log n)$ で、
ループの中の計算は乗算と除算1回ずつゆえ。(証明終)