VBAHaskellの紹介 その21（多倍長整数 bigInt）

前回の記事は有理数計算のサンブルだったので、今回は多倍長整数の計算。

有理数の計算で、オーバーフローに対するガードがないことを不満点としてあげたが、既約分数のべき乗を計算したりすると分母や分子がLongの範囲を超えてしまうのは仕方がない。4バイトとか8バイトといったサイズの制約を超えるためには、配列A(i)を用意して `Σ A(i)・N^i` という形で整数を表現すればいいだけなので簡単だ。Nを 2^15 以下にしておけば掛け算する時の一時的なオーバーフローの心配もない。

今回もデータ構造をクラス化したりせず、単純な配列で値を保持した。**関数一覧**にある機能はだいたい明らかだと思う。


```vbnet
b1 = str2bigInt("0xABCDEF0123456789")    ' 文字列から構築
?bigInt2str(b1)
12379813738877118345
b2 = bigInt_pow(b1, 22)                  ' べき乗
?bigInt2str(b2)
1095747356866942874670236654359969216276597683447788174955606986733509423670022225232131983743614090882854465815739052044839779267369899822789127883550951615180272664979938646672434148953351414478944827697602226184465600916188807253315306160536630641445850300617883726624620685918980495021734095162050648727094223906841315178756746911585634217866565736702066664882279583761024961817642917680220796189606573581695556640625

' 階乗  1000!
fac1000 = foldl1(p_bigInt_mult, mapF(p_long2bigInt, iota(1,1000)))
fac1000s = bigInt2str(fac1000)
?fac1000s
402387260077093773543702433923003985719374864210714632543799910429938512398629020592044208486969404800479988610197196058631666872994808558901323829669944590997424504087073759918823627727188732519779505950995276120874975462497043601418278094646496291056393887437886487337119181045825783647849977012476632889835955735432513185323958463075557409114262417474349347553428646576611667797396668820291207379143853719588249808126867838374559731746136085379534524221586593201928090878297308431392844403281231558611036976801357304216168747609675871348312025478589320767169132448426236131412508780208000261683151027341827977704784635868170164365024153691398281264810213092761244896359928705114964975419909342221566832572080821333186116811553615836546984046708975602900950537616475847728421889679646244945160765353408198901385442487984959953319101723355556602139450399736280750137837615307127761926849034352625200015888535147331611702103968175921510907788019393178114194545257223865541461062892187960223838971476088506276862967146674697562911234082439208160153780889893964518263243671616762179168909779911903754031274622289988005195444414282012187361745992642956581746628302955570299024324153181617210465832036786906117260158783520751516284225540265170483304226143974286933061690897968482590125458327168226458066526769958652682272807075781391858178889652208164348344825993266043367660176999612831860788386150279465955131156552036093988180612138558600301435694527224206344631797460594682573103790084024432438465657245014402821885252470935190620929023136493273497565513958720559654228749774011413346962715422845862377387538230483865688976461927383814900140767310446640259899490222221765904339901886018566526485061799702356193897017860040811889729918311021171229845901641921068884387121855646124960798722908519296819372388642614839657382291123125024186649353143970137428531926649875337218940694281434118520158014123344828015051399694290153483077644569099073152433278288269864602789864321139083506217095002597389863554277196742822248757586765752344220207573630569498825087968928162753848863396909959826280956121450994871701244516461260379029309120889086942028510640182154399457156805941872748998094254742173582401063677404595741785160829230135358081840096996372524230560855903700624271243416909004153690105933983835777939410970027753472000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
' 個々の数字をひとつずつ切り離して配列化
fac1000m = mapF(p_mid(fac1000s), zip(iota(1, Len(fac1000s)), repeat(1, Len(fac1000s))))
digit = Array("0", "1", "2", "3", "4", "5", "6", "7", "8", "9")
' 0,1,2,3,4,5,6,7,8,9それぞれの数字の出現する回数を調べる
printM mapF(p_count_if(lambdaExpr(p_equal,1,ph_1), fac1000m), digit)
  472  239  248  216  229  213  231  217  257  246
' "0"が断然多くてそれ以外は出現回数にそれほど差はなかった
```
（bigInt2str 関数による表示はとてつもなく遅かったが改善させた。）

**関数一覧**

|   関数      |    内容   |        |
|:----------------|:-------- |:-------- |
|long2bigInt|LongからbigIntを生成  ||
|bigInt2double|bigIntからDoubleを生成||
|double2bigInt|DoubleからbigIntを生成||
|bigInt_log   |bigIntの対数||
|bigInt_sgn|bigIntの符号||
|bigInt_base|bigIntの基数(2 ^15)||
|bigInt_abs|bigIntの絶対値||
|bigInt_plus|bigIntの加算||
|bigInt_minus     |bigIntの減算||
|bigInt_mult|bigIntの乗算||
|bigInt_divide_mod  |bigIntの除算 (商とMod)||
|bigInt_divide|bigIntの除算||
|bigInt_mod|bigIntのMod||
|bigInt_pow|bigIntのベキ乗||
|bigInt2str|bigIntからStringへの変換（10進表示）||
|str2bigInt|StringからbigIntへの変換||
|bigInt_equal|bigIntの比較 (a = b)||
|bigInt_not_equal|bigIntの比較 (a <> b)||
|bigInt_less|bigIntの比較 (a < b)||
|bigInt_less_equal|bigIntの比較 (a <= b)||
|bigInt_greater|bigIntの比較 (a > b)||
|bigInt_greater_equal|bigIntの比較 (a >= b)||
|bigInt_gcd|最大公約数||

Double値からbigIntを生成する関数 `double2bigInt` は正確な値を求めるものではない。この逆の操作である `bigInt2double` はオーバーフローに対するガードはしていないので、Double の最大値である ±e+308 を超えている場合の挙動はわからない。有理数の分母と分子をこれで表現する対応はまだやっていない。

多倍長整数 x を配列 m で表現する内部構造は以下の通り。

|配列要素|m(0)|m(1)|m(2)|・・・|m(N+1)|
|:------|:---|:---|:---|:---|:---|
|内容|±基数B|A(0)|A(1)|・・・|　A(N)|
このとき、
x = ±A(0) + A(1)・B^1 + A(2)・B^2 + ・・・ + A(N)・B^N
である。基数Bは 2 ^ 15 = 32768 固定だが、10進表示の関数 `bigInt2str` を実装する際に一時的に B = 10000 のものを作るときがある。
「±基数B」 とあるのは、基数の値と数値そのものの符号を兼ねていて、たとえば 32770 と -32770 の表現はこうなる。
 32770 → Array(32768 , 2, 1)　・・・　　2 + 1 * 32768 ^ 1
-32770 → Array(-32768 , 2, 1)  ・・・　-(2 + 1 * 32768 ^ 1)


今回もHaskellとは全然関係ない話題だし、bigIntの実装にもほとんど使っていない。
VBAソースコードは [misc_ratio.bas](https://github.com/mYmd/VBA/blob/master/misc_ratio.bas) にある。

* * *

[VBAHaskellの紹介 その20（有理数の計算）](http://qiita.com/mmYYmmdd/items/cc7515e18338ef0fc765)
[VBAHaskellの紹介 その1（最初はmapF）](http://qiita.com/mmYYmmdd/items/c731edf943acc0a0ebe9)
ソース https://github.com/mYmd/VBA