多倍長加減算
多倍長整数
多倍長整数での加減算はそれこそ筆算通りの計算でできる。 気持ち程度の選択肢として、計算途中に正規化するものと最後に正規化するものとがある。 また、符号を別に持つ形で処理を行う場合は両者の符号関係や絶対値の大きさを比較して処理する必要がある一方、 符号を持たない、非負整数を扱う場合には減算結果が負になったらどうするか決めておくか、 必ず負にならないように演算を適用させる必要がある。
\[ A = a_{n_A-1}R^{n_A-1} + a_{n_A-2}R^{n_A-2} + \cdots + a_1R + a_0 \] \[ B = b_{n_B-1}R^{n_B-1} + b_{n_B-2}R^{n_B-2} + \cdots + b_1R + b_0 \] \[ \Rightarrow C = A \pm B = c_{n_C-1}R^{n_C-1} + c_{n_C-2}R^{n_C-2} + \cdots + c_1R + c_0 \]上記のような計算モデルで計算を行う。
後で正規化を行う
\[ \begin{array}{rrrrrr} &31&41&59&26&53\\ +& & &58&97&93\\ \hline &31&41&117&123&146\\ \hline \text{(正規化)}&31&42&18&24&46 \end{array} \]まずは後で正規化を行う手法から。 上の計算イメージのようにとりあえず要素ごとの和を求め、後で正規化を行う。 この方式では複数の正規化ルールを取り換える場合に加算ルーチンが 1 つで済む他、 横並びのデータに依存がないので並列化、ベクトル化しやすい。 また本来やりたいことだけ記述するのでコードの可読性も高くなりやすい。
2 つの多倍長数 $a$ と $b$ の和を求めるプログラムイメージを書く。
BigUInt UnsignedAdd(const BigUInt& a, const BigUInt& b) {
BigUInt c(std::max(a.size(), b.size()));
// 桁毎に足す
for (size_t i = 0; i < c.size(); ++i) {
Digit ai = (i < a.size()) ? a[i] : 0; // 桁数を超えると 0 にする
Digit bi = (i < b.size()) ? b[i] : 0;
c[i] = ai + bi;
}
// 正規化
Digit carry = 0;
for (Digit& ci : c) {
ci = ci + carry;
carry = ci / R;
ci = ci % R;
}
if (carry)
c.push_back(carry); // 最後に繰り上がりがあれば桁数を増やす
return c;
}
このプログラムでは $R$ を少し超えた数でも Digit
型がオーバーフローしない前提で書いているので
例えば 64 ビット計算機で $R=2^{64}$ とするなど
オーバーフローすることが前提となる場合には使えない。
後で正規化を行うやり方はこのような形でサイズ制限がある。
計算中に正規化を行う
続いて、計算途中で正規化を行うバージョン。 こちらは上記のようなサイズ制限を回避しやすい他、 計算データを複数回読み込まなくて良いのでキャッシュ的に嬉しい。
\[ \begin{array}{rrrrrr} &31&41&59&26&53\\ +& & &58&97&93\\ \hline &&&&(1)&46 \end{array} \] \[ \begin{array}{rrrrrr} &31&41&59&26&53\\ +& & &58&97&93\\ \hline &&&(1)&24&46 \end{array} \] \[ \begin{array}{rrrrrr} &31&41&59&26&53\\ +& & &58&97&93\\ \hline &&(1)&18&24&46 \end{array} \] \[ \begin{array}{rrrrrr} &31&41&59&26&53\\ +& & &58&97&93\\ \hline &31&42&18&24&46 \end{array} \]この場合のプログラムは以下の通り。
BigUInt UnsignedAdd(const BigUInt& a, const BigUInt& b) {
BigUInt c(max(a.size(), b.size()));
Digit carry = 0;
for (size_t i = 0; i < c.size(); ++i) {
Digit ai = (i < a.size()) ? a[i] : 0;
Digit bi = (i < b.size()) ? b[i] : 0;
Digit ci = ai + bi + carry;
carry = ci / R;
c[i] = ci % R;
}
if (carrry)
c.push_back(carry);
return c;
}
こちらの手法であれば Digit が 64bit 符号なし整数で $R=2^{64}$
というような場合でもループ内の計算を例えば
Digit ci = ai + carry; carry = (ci < ai) ? 1 : 0; ci = ci + bi; if (ci < bi) carry = 1;
とすることで計算できるようになる。
多倍長浮動小数点数
浮動小数では仮数部の有効桁数と指数部との関係で足し合わせるインデックスにズレが出てくるが、 基本的には整数と違いが無い。 指数が異なる場合には片方の数値の末尾が端数となるが、好きな方法で丸めると良い。