Hundredsoft

除算 x = p / q があったとして
ニュートン法（２次漸化式）を使って(q)の逆数を求め、それに(p)を乗算する方法は多桁演算の定石です。
この場合、多桁の浮動小数点型に変換して演算することが多いでしょう。
ここでは、ニュートン法を整数演算のみで実装する方法について考えます。
高次漸化式にも応用できますが、ここでは２次漸化式を使っています。

1/q の２次漸化式は、

X(i+1) = X(i) * (2 - q * X(i))

と表せます。
ここで被除数のBit数をpb, 除数のBit数をqbとして、n = pb + qb とすれば

Bigint iDiv(Bigint p, Bigint q){
  int pb = p.length(); // 被除数のbit数
  int qb = q.length(); // 除数のbit数
  int n = pb + qb;

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

  while(m != x){
    m = x;
    x *= c2 - q * x;
    x >>= n;
  }
  x *= p;        // x = p / q
  x >>= n;

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

と書けます。
考え方としては、計算各項を　(* 2ⁿ) することで整数を一時的に固定小数点化するものです。
演算量を減らす変形を重ねているのでわかりにくいかもしれませんが、非常に短いコードです。
但し、初期値が甘いのでループ回数はやや多いです。
被除数と除数の桁差が少ない場合は、回復法を使った２進多桁の除算が有利です。


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

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

---

Tags: プログラムメモ