「なんちゃって」BCD演算クラスを作って、nCr(組み合わせ数)計算した話

発端

本質は、「ワークシート関数「COMBIN」の結果をmodに入力すると、ありえない数値が返ってくる。」だったのです。

が、なぜかそこから、「必ず整数になるはずなので、誤差ができるのはおかしい」

という話になって、VBAで計算し始める方がちらほら。

もっと居たような気がしたのですけど、私とはけたさん以外、どなたがいたのでしたっけ。

はけたさんが書かれてるのは「パスカルの三角形」と呼ばれる二項展開の手法であるようです。。

恥ずかしながら、知らなかったのです。中高大と数学好きで通してたつもりだったのですが、サッカリンレベルで甘かったです。

でも、桁数足りないよね

で、私も作ってみたのですが、かなりぜい肉だらけのロジック。

でも、分母は毎回きちんと約分完了してたし、Excelが計算してくれる桁数には追い付けたし、まあ、いいかなあ、あとは桁数ですね。BCDライブラリ欲しいですね、というところまで進みました。

漸化式

この発言で様相が一変したと言っても過言ではないと思うのです。

実はこの計算方法、はけたさんがツイートしたものと本質的には同じ考え方なのですが、大きなメリットとして、「再起呼び出しが不要」という点があります。

VBAは再起呼び出ししやすい言語とは言えないので、この差は大きいと思うのです。

これまでのやりとりで、「BCDライブラリ必要だよね」というところまで議論は進んでいたのです。

が、8bit機でマシン語やってた方は理解できると思うのですが、掛け算とか割り算って、計算コストがそこそこ大きいのです。

CPUの演算命令に用意されてなかったら自作しないといけなかったし、CPU組み込みも昔は遅かったのです。

それをVBAでやるとなると、結構大変なのです。

まあ、やることは筆算と変わらなくて、桁上がりを検知して上の桁に1を加算したりするだけなのですが、単純な加算の方が恐ろしく単純に解決できます。

で、黒猫と香辛料さんの発言のおかげで、「加算さえ用意すれば演算できる」ことが分かってしまったわけです。

そしてBCDへ

加算だけで良いなら、作らないと損ですね、ってわけで、晩酌の肴に組んだのが以下のクラス。

Option Explicit

'クラス変数定義
Dim byteDigit_() As Byte     'BCDデータ
Dim lngMax_ As Long          '桁数

'--------------------------------------------------------------------------------
'疑似コンストラクタ
'引数:strArg As String
'・引数strArgの文字数からlngArg算出
'・クラス変数lngMaxに値セット
'・クラス変数byteDigit()を初期化
'--------------------------------------------------------------------------------
Public Sub Initialize(strArg As String)
Dim lngArg As Long
Dim lngIdx As Long

    lngArg = Len(strArg)

    '・クラス変数lngMaxに値セット
    Me.lngMax = lngArg
    
    '・クラス変数byteDigit()を初期化
    ReDim byteDigit_(lngArg - 1) As Byte
    
    'BCD変換しつつ代入
    '文字列の1桁目→配列子0の配列へ、and so on.
    For lngIdx = 0 To Me.lngMax - 1
        byteDigit_(lngIdx) = Asc(Mid(strArg, Me.lngMax - lngIdx, 1)) - Asc("0")
    Next
End Sub

'--------------------------------------------------------------------------------
'lngMax プロパティ
'--------------------------------------------------------------------------------
Public Property Get lngMax() As Long
    lngMax = lngMax_
End Property
 
Public Property Let lngMax(lngArg As Long)
    lngMax_ = lngArg
End Property
 
'--------------------------------------------------------------------------------
'byteDigit プロパティ
'--------------------------------------------------------------------------------
Public Property Get byteDigit() As Byte()
    byteDigit = byteDigit_
End Property
 
Public Property Let byteDigit(byteArg() As Byte)
    byteDigit_ = byteArg
End Property
 

'--------------------------------------------------------------------------------
'toString
'BCD表現のByte列を文字列へ変換
'--------------------------------------------------------------------------------
Public Function toString(byteArg() As Byte) As String
Dim strTemp As String
Dim lngIdx As Long

    strTemp = ""

    'BCD変換しつつ代入
    '文字列の1桁目→配列子0の配列へ、and so on.
    For lngIdx = 0 To UBound(byteArg)
        strTemp = Chr(byteArg(lngIdx) + Asc("0")) & strTemp
    Next
    toString = strTemp
End Function

'--------------------------------------------------------------------------------
'toBCD
'文字列をBCD表現のByte列へ変換
'--------------------------------------------------------------------------------
Public Function toBCD(strArg As String) As Byte()
Dim lngArg As Long
Dim lngIdx As Long
Dim lngMax As Long
Dim byteDigit() As Byte
    lngArg = Len(strArg)

    lngMax = lngArg
    
    'byteDigit()を初期化
    ReDim byteDigit(lngArg - 1) As Byte
    
    'BCD変換しつつ代入
    '文字列の1桁目→配列子0の配列へ、and so on.
    For lngIdx = 0 To lngMax - 1
        byteDigit(lngIdx) = Asc(Mid(strArg, lngMax - lngIdx - 1, 1)) - Asc("0")
    Next
End Function

'--------------------------------------------------------------------------------
'AddBCD
'メンバ変数のBCD表現のByte列に、引数に指定したBCD表現のByte列を加算する
'--------------------------------------------------------------------------------
Public Sub addBCD(byteArg() As Byte)
Dim lngIdx As Long
Dim lngMax As Long
Dim lngTemp As Long
Dim lngArgMax As Long

Dim lngArgA As Long
Dim lngArgB As Long

lngArgMax = UBound(byteArg)

    'お互いの桁数で小さい方の最大桁数-1まで計算
    For lngIdx = 0 To Max(lngArgMax, Me.lngMax - 1)
        '加算
        lngArgA = 0
        lngArgB = 0
        If lngIdx > lngArgMax Then
            '引数に指定したBCD表現のByte列を加算し終わっているので、なにもしない
        Else
            lngArgA = byteArg(lngIdx)
        End If
        
        If lngIdx > Me.lngMax - 1 Then
            'BCD表現のByte列が足りなくなったので桁を追加
            ReDim Preserve byteDigit_(Me.lngMax) As Byte
            Me.lngMax = Me.lngMax + 1
        Else
            lngArgB = byteDigit_(lngIdx)
        End If
        
        lngTemp = lngArgA + lngArgB
        
        If lngTemp > 9 Then
            '桁上がり発生
            '上位の桁が今の桁数より大きい場合、Byte列を増やす
            If (lngIdx + 1) >= Me.lngMax Then
                ReDim Preserve byteDigit_(Me.lngMax) As Byte
                Me.lngMax = Me.lngMax + 1
            End If
            
            '現在桁へ-10した値を設定
            byteDigit_(lngIdx) = lngTemp - 10
            
            '上の桁へ1を加算
            byteDigit_(lngIdx + 1) = byteDigit_(lngIdx + 1) + 1
        Else
            '桁上がり発生せず
            '現在桁へ値を設定
            byteDigit_(lngIdx) = lngTemp
        End If
    Next
End Sub

Private Function Min(lngA As Long, lngB As Long) As Long
    Min = IIf(lngA >= lngB, lngB, lngA)
End Function

Private Function Max(lngA As Long, lngB As Long) As Long
    Max = IIf(lngA >= lngB, lngA, lngB)
End Function

そして、実際に足し算を行うときに、入力を「文字列」にするためのラッパー関数を用意しました。

Function addBCD(strArgA As String, strArgB As String) As String
Dim bcdA As New bcd
Dim bcdB As New bcd

    Call bcdA.Initialize(strArgA)
    Call bcdB.Initialize(strArgB)
    
    '加算
    Call bcdA.addBCD(bcdB.byteDigit)
    
    '戻り値を返す
    addBCD = bcdA.toString(bcdA.byteDigit)
End Function

自画自賛ですが、このラッピングの仕方はすごく理想的で、「データ構造を無視できる」のが強みです。

内部処理でBCD演算のためにByte型配列を扱っていて、かつ、桁数が足りなくなると自動で桁数増やしている、なんてことを全く意識せずに足し算だけ任せられるのは便利なのです。

Function nCrBCD(lngN As Long, lngR As Long) As String
Dim strArg() As String   '演算用の文字列
Dim lngIdx As Long
Dim lngIdx2 As Long
Dim dateStart

dateStart = Now()

ReDim strArg(lngR) As String

    '演算用文字列初期化
    For lngIdx = 0 To lngR
        strArg(lngIdx) = "1"
    Next

    For lngIdx2 = 1 To lngN - lngR
        For lngIdx = 1 To lngR
            strArg(lngIdx) = addBCD(strArg(lngIdx), strArg(lngIdx - 1))
        Next
    Next

Debug.Print CDate(Now() - dateStart)

    nCrBCD = strArg(lngR)
End Function

そして最後に、黒猫と香辛料さんのコードを、BCD演算で動くようにして完了。

配列確保も動的にしたので、理論上は、「コンピューターのスペックが許す限り、何桁でも計算できるようになりました。」

で、どこまで計算できるの

確認した範囲でこんな感じでした。

1000C50:1秒くらい

10000C50:26秒くらい

100000C50:5分59秒くらい

ちなみに、100000C50は186桁で、VBAとCONMINで値はそれぞれ以下の通りです。

VBA:324791116448528873575523182071046314176283476648027244927270999288005964035859548411037117880276576299533569666297949933317955323578121238572954027153746965881539586782982104673414248000

CONBIN:324791116448529000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

CONBIN関数が上15桁で切れてるのがわかりますね。BCD大好きってなりました。

こんな演算を見ると、「COBOL無くならないな」って感想が出てきたりします。

CやC++にはBCDライブラリありますけどデフォルトではないのでデフォルトで使えるメリットはかなり大きいのです。

次へ

BCDライブラリ、作るしかないんじゃない?と、思っています。

整数演算だけでなく小数点演算も必要なので工夫が必要ですが。

ちなみに、そういえば、BCD使った演算ルーチンの考察って、25年前にやってたな、と、思い出したので、週末当たり昔のノートを見てみようと思います。

あと、今回掲載したクラスは出来が良くないです。

速度を考えてもっとギークなつくりにできるはずなのです。

例えば、BCDでの1桁加算演算、マシン語だと単純にレジスター足し算してキャリフラグで分岐、なのですけど、VBAだと、同様に、select Case文で実現可能なのです。

速度検証してないのですけど、理論上、10進数1桁の演算は、結果が100通りしかないのですよね。10*10の演算なので。

その辺りも含めて、クラスは作り直したいと思うのです。わくわくしてきますね!

久しぶりのブログはこんな感じで締めたいと思います。

なお、まだお酒が抜けてない酔っ払いです。ぷはー。