32ビット環境で64ビット整数を扱う (乗法編)

これまでの記事では、32ビット環境における64ビットの加法および減法について解説しました。いよいよ、ヤマ場である乗法、すなわち「掛け算」のやり方について解説しましょう。

なお、今回も64ビット整数値を表現するために、以下の QWORD 構造体を使用します。

typedef unsigned long DWORD;

typedef struct {
    DWORD dwLow;  //下位32ビット
    DWORD dwHigh; //上位32ビット
} QWORD;

基本的な考え方

コンピュータを使わずに掛け算を行う場合、殆どの人は筆算を用いるのではないでしょうか。今回の記事で紹介するプログラムにおいても、掛け算の実行手順はこの筆算とまったく同じです。

ただし、10進ではなく2進記法を用いる点が通常の筆算とは異なります。例えば、142×75を計算する場合は、オペランドをそれぞれ 2進数に直して、次のような「筆算」を行うことになります。

まず、第一オペランド (上段の 10001110 = 142₁₀) と第二オペランド (下段の 1001011 = 75₁₀) の一の位の数の積を計算します。積といっても、一の位の数字は 1 なので、第一オペランドの値をそのまま下に書き写すだけになります。

次に、第二オペランドの十の位 (二進法だから、正確には「二の位」) の数の積を計算します。通常の筆算と同じく、この積を書き出すときには、左に一桁ずらす (左に1ビットシフトする) 必要があります。

さらにその次、第二オペランドの百の位 (正確には「四の位」) の数字は 0 なので、第一オペランドとの積も当然ながら 0 になります。 0 は何ビットシフトしても 0 ですが、図には形式的に左に2ビットシフトしたものを載せました。

この動作を最上位の位まで繰り返すと、書き出される積の列は以下のようになります。あとはこれをすべて足し合わせれば求めるべき値、すなわち全体の積が得られます。 (加法についてはこちらの記事を参照してください。)

上の図では、スペースを節約するため、値が 0 の列を省略してありますが、実際の計算過程でも 0 の加算過程はスキップして差し支えありません。

32ビット乗算

さて、前節の説明をもとに、実際にこの計算を実行するプログラムを書いてみましょう。まずは、32ビット値の乗算をする関数を作るところから始めます。「32ビットの乗算なんて、普通に * 演算子を使えばできるじゃん。」と思われるかも知れませんが、二つの32ビット値の積の論理的な最大値は、

(2³² - 1)²
= 2⁶⁴ - 2³³ + 1
= 18,446,744,065,119,617,025

であり、32ビットの範囲には収まりません。そこで、この32ビットの範囲を超過した値を繰り上がり (carry) として取得することのできる乗算関数 mul32 を以下のように定義します。 (mul32 の中で加算を行うため使用されている add32 は、加法編で定義したものと同一。)

#define LOWORD(dw) (0xFFFFUL & dw)
#define HIWORD(dw) (0xFFFFUL & (dw >> 16))

typedef unsigned short WORD;

//加算 (32ビット)
DWORD add32(DWORD dw1, DWORD dw2, DWORD* lpdwCarry =NULL){

    DWORD dwTmp1 =LOWORD(dw1) + LOWORD(dw2);
    DWORD dwTmp2 =HIWORD(dw1) + HIWORD(dw2) + HIWORD(dwTmp1);

    if (lpdwCarry){
        *lpdwCarry =HIWORD(dwTmp2);
    }

    return LOWORD(dwTmp1) | (dwTmp2 << 16);
}

//乗算 (32ビット)
DWORD mul32(DWORD dw1, DWORD dw2, DWORD* lpdwCarry){

    //オペランドをそれぞれ上位, 下位16ビットに分解
    WORD wLow1  =LOWORD(dw1);
    WORD wLow2  =LOWORD(dw2);
    WORD wHigh1 =HIWORD(dw1);
    WORD wHigh2 =HIWORD(dw2);

    //16ビット乗算
    DWORD dwTmp1 =wLow1*wLow2;
    DWORD dwTmp2 =wLow1*wHigh2;
    DWORD dwTmp3 =wHigh1*wLow2;
    DWORD dwTmp4 =wHigh1*wHigh2;

    //戻り値 (積の下位32ビット) を計算
    DWORD dwRes, dwCarry1, dwCarry2;
    dwRes =add32(dwTmp1, LOWORD(dwTmp2) << 16, &dwCarry1);
    dwRes =add32(dwRes,  LOWORD(dwTmp3) << 16, &dwCarry2);

    //繰り上がり (積の上位32ビット) を計算・格納
    if (lpdwCarry){

        DWORD dwCarry;
        dwCarry =add32(dwTmp4,  HIWORD(dwTmp2));
        dwCarry =add32(dwCarry, HIWORD(dwTmp3));
        dwCarry =add32(dwCarry, dwCarry1);
        dwCarry =add32(dwCarry, dwCarry2);

        *lpdwCarry =dwCarry;
    }

    return dwRes;
}

積のデータ長は64ビット。そのうち下位32ビットは戻り値として返され、上位32ビットはポインタ lpdwCarry で指定されたアドレスに格納されます。上記のプログラムは、前節で紹介したように二つのオペランド dw1, dw2 に対して1ビットずつ計算していくのではなく、それぞれを下位・上位16ビットずつ (wLow1, wHigh1 と wLow2, wHigh2) に分けまとめて処理を単位で乗算を行っていますが、基本的な考え方に変わりはありません。(下図参照)

|0|0 codelogy|53|3|64ビット乗算|

さて、これでようやく準備が整いました。前節で定義した mul32 を使って、64ビット乗算を行う関数 mul64 を定義します。 ...といっても、この関数の定義は拍子抜けするほど簡単です。

//乗算
QWORD mul64(DWORD dwLow1, DWORD dwHigh1, DWORD dwLow2, DWORD dwHigh2){

    DWORD dwCarry;

    QWORD qw;
    qw.dwLow  =mul32(dwLow1, dwLow2, &dwCarry);
    qw.dwHigh =mul32(dwLow1, dwHigh2) + mul32(dwHigh1, dwLow2) + dwCarry;

    return qw;
}

二つの32ビット値の積を表現するには64ビットのデータサイズが必要であったのと同様、二つの64ビット値の積を表現するには128ビットのデータサイズが必要となります。しかし、このエントリで扱うのはあくまでの64ビット整数の範囲での演算であり、上位の64ビットの値は扱いようがないので、オーバーフローしたものとみなして捨てられてしまいます。コード中で、qw.dwHigh に代入する値の計算が add32 ではなく + 演算子を用いて行われているのは、この理由で繰り上がりを考慮する必要がないため。また、オペランドの上位16ビット (dwHigh1, dwHigh2) の積はすべてオーバーフロー領域上にあるため、その計算は省略することができます。

次回予告

ここまでで、64ビット整数の加法, 減法, 乗法についての解説が完了しました。この順番でいくと次は除法というのが順当に思われますが、その前に必要になると思われる機能があります。それは、「64ビットの値を10進で出力」する機能。 QWORD の形式で得られた数値は、2進数, 8進, 16進で出力するのはそれほど難しくはありません (8進はちょっとヒネリが必要) が、10進でとなるとそれなりに手間がかかります。

そこで次回は、この10進出力を行う方法についての解説をしたいと思います。お楽しみに。

成田（七の段が苦手）