Hundredsoft

多桁の整数値(s)の平方根 √s があったとして
ニュートン法（２次漸化式）を使って(s)の平方根の整数部分を求める方法について考えます。
与える引数(s)に基数のべき乗を掛けることで任意の精度に拡張できます。 __前回同様__ 高次漸化式にも応用できますが、ここでは２次漸化式を使っています。

√s を展開すると除算が必要になるので, (1/√s)を求め、これに(s)を乗じることとします。
(1/√s) の２次漸化式は、

X(i+1) = X(i) * (3 - s * X(i)^2) / 2

と表せます。
これは、

Bigint iSqrt(Bigint s){
  int n = s.length(); // bit数

  Bigint x = 1;
  x <<= (n / 2);  // Newton法の初期値として, 1 << (n/2) を与える。
  Bigint m(0);
  Bigint t3 = 3;
  t3 <<= n;

  while(m != x){
    m = x;
    x *= t3 - ((s * x * x) >> n);
    x >>= (n + 1);
  }
  x *= s;         // x = s / √s = √s
  x >>= n;

  // 誤差修正
  m = x + 1;
  if (s >= m * m){
     x = m;
  }
  return x;
}

と書けます。
考え方としては、計算各項を　(* 2ⁿ) することで整数を一時的に固定小数点化するものです。
演算量を減らす変形を重ねているのでわかりにくいかもしれませんが、非常に短いコードです。
初期値については、平方根は (n/2)桁程度になるので、その逆数もまた 2ⁿ倍することを考えれば、
(n/2)桁程度であることから決定していますが、初期値が甘いのでループ回数はやや多いです。
ニュートン法を使った除算については、__こちら__をご覧ください。

「Newton法初期値の改良バージョン」はこちら
ニュートン法を使った２進多桁の整数平方根(その２)

(注)ここでのBigintクラスは、多倍長整数を格納する仮想のクラスであり、実在するものではありません。

---

Tags: プログラムメモ