巨大な数を計算するJavaScriptプログラムの開発ノート。
処理内容 | Ver.0.3 | Ver.0.4 |
---|---|---|
2の10000乗 | 1300ms | 450ms |
1000の階乗 | 1000ms | 350ms |
400桁の10進数を2進数に変換 | 1000ms | 250ms |
ルート2小数点以下1024桁 | 2000ms | 600ms |
以下、巨大数演算Ver.0.4の開発記録。作ったプログラムのバージョン一覧。
通常乗算を行うときは文字列データを切り出して配列で計算し、計算後に配列を文字列に戻して出力する。 しかし連続で乗算を行う際は文字列に変換するのは無駄な処理だ。 そこで複数のデータを一気に入力しそれらを全て掛け合わせるlong_mul関数を作成。
累乗(pow)・剰余(pow_surp)の演算速度が1.5倍になった。両者共改良したのは後半だけで、全計算を配列でやらせるとさらに50%高速化ができると考えられる。
また、連続で除算をする場合は逆に文字列を維持したまま計算することで高速化できる。
新しいlong_mul関数にて、計算結果の末尾がゼロの時結果が狂うミスを発見。修正した。
あとdemoページに埋め込み易いように仕様を変更した。
速度は速くなっていないが、除算の15桁以上のアルゴリズムを改良。今まではわざわざ余りを被除数に還元して次のループで再度切り出していたが、必要な桁数だけをピンポイントに切り出すようにした。ソースコードを美しく書いた程度の改良だが、ソースの分かりやすさは後の改良に役立つ。
階乗を50%高速化。5000!の計算時間が5秒台になった。
通常乗算,2乗乗算,7桁以下の乗算も高速化。2^40000にかかるコストを6.5秒から5.3秒に短縮。円周率も1000桁の場合で11秒から10秒になった。
総合的なアルゴリズムではなく、式変形や代入の最適化だけでまだまだ高速化ができそうだ。
階乗デモを30%ほど高速化。デバッグ処理で使っていたif構文で処理速度が鈍っていたので、1000ごとにループ計算を区切ることによって最終的なif構文の回数を減らした。16326桁の5000!なら今まで14秒かかったがこれだと10秒で計算できる。
ちなみにVer.0.3時代は5000!に43秒もかかっていたりする。
toStringに100進数に変換するデモ(仮)を追加。巨大整数をの長さを半分に圧縮できる。
toString100 (169digits) i0!0>GD]dRJ&fQkCq>'Yr!1u$OQ'=!6J\G339L$]iS`}!3'rRV! 3QWo!2M/|OqU=@25=9!1!0B!19#7J}|#O!07|qIB!6J!5!3)po Y8o{!6].gJS0$!6(I]<P!6|;24;]!6#Vi?sY#o$TWn!7$)J,!4 vR1;*};Y8K4~#=+@W\ Reverse (300digits) 7013922736 3358654739 0367467232 7827045479 9382014446 0426983957 3616162241 0158704861 9095047947 5195465276 9442128944 7850262915 1826229392 3193220020 3990890044 9220897838 3198399795 0677765421 7688985811 6839481301 9805385825 4598892415 1724589800 5170288054 0076014952 7599010639 0996834714 2407902454 2140179100 2608295257
除算で同時に求める商の桁数を14桁に戻した。もし商に±2以上の誤差が出れば以下のようなwhile構文で補正。かなり無理矢理だが平方根も円周率も巨大剰余も演算ミスは起きていない。
while( broken_end.flag ) { broken_end = HugeInt.add( broken_end,Number2 ); divided_broken_Val--; } while( HugeInt.max_min( broken_end,Number2 ) ) { broken_end = HugeInt.sub( broken_end,Number2 ); divided_broken_Val++; }
除算単体で40%ほど高速化。平方根,円周率も速くなった。
7桁以下の乗算も高速化。2数のうち少なくとも片方が7桁以下のとき元々階乗用に作ったプログラムで計算する。 これだと単に2倍するときなどは14桁の配列になり演算効率がざっと倍になる。
コードを整理。alpha12で円周率の計算を達成したのでbetaにした。
単に雑なソースコードを整理しただけなので演算速度は同じ。
除算関数の295桁以上用のものを294桁のときと同じ方法で計算する方式にした。 オーバーフローした除数と同じ桁数だけ被除数もオーバーフローさせることでJavaScriptの有効範囲内(306桁)に収める。
手元では6秒近くかかった2*104000の平方根が2秒になった。 次のバージョンでガウス・ルジャンドルを用いた円周率計算の実験を行う予定。 そもそもこの企画は円周率を計算するために始めたものなので、円周率の計算に成功すれば晴れてベータ版となる。
マチンの公式を用いた円周率計算関数とガウス・ルジャンドルを用いた円周率計算関数を追加。 マチンの公式では1分弱で円周率1000桁を計算できる。これに対しガウス・ルジャンドルは14秒で1000桁を計算する、2000桁なら1分以内に計算できる。 なお、マチンの公式はπが完全に収束するまでループを繰り返すが、ガウス・ルジャンドルはループ前に必要なループ回数を概算するので最後数桁の値が実際のπの値と多少ズレる。
以下、ガウス・ルジャンドルによる円周率500桁の計算結果。6月下旬にやったときはこの500桁に150秒もかかった。
π 500digits π = 3. 1415926535 8979323846 2643383279 5028841971 6939937510 5820974944 5923078164 0628620899 8628034825 3421170679 8214808651 3282306647 0938446095 5058223172 5359408128 4811174502 8410270193 8521105559 6446229489 5493038196 4428810975 6659334461 2847564823 3786783165 2712019091 4564856692 3460348610 4543266482 1339360726 0249141273 7245870066 0631558817 4881520920 9628292540 9171536436 7892590360 0113305305 4882046652 1384146951 9415116094 3305727036 5759591953 0921861173 8193261179 3105118548 0744623799 6274956735 1885752724 8912279381 8301194870 (501digits)
_HugeInt_Version:Ver.0.4 alpha12 _HugeInt_Current:Wed Nov 30 16:21:35 UTC+0900 2005 _HugeInt_UserAgent:Mozilla/4.0 (compatible; MSIE 6.0; Windows NT 5.0; .NET CLR 1.1.4322) _HugeInt_Mode:pi DEBUG: HugeInt.pi_Gauss TEST DEBUG: [Gauss-Legendre] DEBUG:input: 500(3 digits) DEBUG: all frequency: 8 DEBUG: settling: 1/8(359ms) DEBUG: settling: 2/8(359ms) DEBUG: settling: 3/8(375ms) DEBUG: settling: 4/8(360ms) DEBUG: settling: 5/8(343ms) DEBUG: settling: 6/8(360ms) DEBUG: settling: 7/8(344ms) DEBUG: settling: 8/8(328ms) DEBUG: calculating... DEBUG: calculated (140ms) DEBUG:output: 31415926535897932384...(501 digits) DEBUG:Cost: 3375ms DEBUG: end
元々このプロジェクトは円周率を計算する目的で始まり、初めて100桁台の計算に成功したのが2005年6月。 しかし速度に不満がありその後半年間開発を続け今に至る。 ガウス・ルジャンドルが完成した今その目的は達成されたのだが、更なる高速化を目指して開発を続行する。
1万桁(精度9998桁)の計算に1900秒ほどかかった。
除算関数全体に致命的なバグが見つかった。 14桁以下では特定の条件下でゼロが付け足されないバグ。294桁以下では11桁同時除算後にフローした分の除算でゼロが付け足されないバグ。 これらが今回作った平方根プログラムの完成を遮っていた。
というわけで平方根。残念ながらあまり高速ではない。 しかし、古いVer.0.3では253以下の平方根の小数点以下を求めるものであったのに対し、今回のものは巨大整数の平方根を求めるものだ。 この場合ルート2の小数点以下1000桁を求めるには2*102000を入力する。巨大整数の平方根はガウス・ルジャンドルで円周率を求めるときに使用する。
平方根の高速化については除数が295桁以上の除算を高速化するのが手っ取り早いが、ループ回数を減らしても多少の高速化が望める。
binary.jsをtoString.jsに変更。 JavaScriptのtoStringメソッドの巨大整数版としてHugeInt.toString関数を追加した。 2〜32進数を比較的高速に展開できる。 しかし2進展開のみをやりたいときにはわずかに遅い。 そのため古いbinary関数も残している。
デバッグなど外部出力系もかなり弄った。 HugeInt._debug関数に渡す引数の順番はHugeInt._debug(種類,値)でだったので、単に文字列をデバッグ欄に表示させるときは種類の設定が面倒。 この順序を逆にした。
さらにMozilla用のプログレスバーを追加。 HugeInt._progress関数でIEならステータスバー,Mozillaならプログレスバーに経過表示を行える。今まではHugeInt._debugでもステータスバーに表示ができたが、無駄なデバッグ表示にいちいちステータスバーが上書きされるのは嫌いなので別々の関数にした。
除算アルゴリズムを改良。1回の概算で上位1桁しか信用してなかった部分を同時に14桁信用した結果である。 一度に計算できる商の桁数を14桁にしたことで15桁桁以上293桁以下において演算速度が飛躍的に上がった。
Ver.0.3-a9 | Ver.0.4-a6 | Ver.0.4-a7 |
---|---|---|
688ms | 625ms | 203ms |
また、pow_surp.jsにおいて巨大剰余演算の速度も一気に上がった。 以下は60桁の素数(979359608342168056772037984398918304826077541588952608951979)をフェルマーテストにかけるのにかかった時間。
処理/Ver. | Ver.0.3-a9 | Ver.0.4-a6 | Ver.0.4-a7 |
---|---|---|---|
二進展開 | 31ms | 0ms | 16ms |
剰余リスト生成 | 3187ms | 3031ms | 766ms |
リストの乗算 | 5891ms | 3109ms | 3062ms |
乗算結果を除算 | 1641ms | 1344ms | 203ms |
合計時間 | 10750ms | 7484ms | 4062ms |
この演算方式を294桁以上で使えるようにすれば、事実上全除算の実行速度が数倍に跳ね上がる。 しかし現段階では商がおかしくなる可能性も否定できない。
また、binary.js―二進展開も一度に行う除算を246から253に切り替えたことにより80%程度高速化。 ただし1000桁級の巨大整数のときだけしかこの高速化結果は出ない。小さな数の二進展開は他の処理で時間を食う。 手元では1000桁の整数を0.8秒で2進数に変換できる。alpha6ではこの処理に1.4秒、Ver.0.3にいたっては6秒もかかった。
古いVer.0.3で実装した素数検索とどれだけ差があるかの実験。 エラトステネスは判定する数が100桁になろうが一瞬で終わる。フェルマーテストは桁数が上がるごとに計算時間も伸びるが、4,50桁の素数でも短時間に検索できた。 RSAの鍵の生成には全く問題ない速度だ。
このデモのプログラムはつぎはぎだらけで汚いが気にしないように。
商を14桁同時計算すると特殊な場合において計算を間違える。処理する桁数を12桁(12,13桁でも特殊な場合で失敗)にしてこの問題を解決したalpha8。 alpha7の素数判定は無限ループを含む場合があるので通常はalpha8の素数判定を使ってください。
pow_surp.jsを追加。フェルマーテストで実行する巨大累乗の剰余を実装した。 しかしこれはver.0.3よりも速いとはいえない。
処理内容 | Ver.0.3 | Ver.0.4 |
---|---|---|
二進展開 | 15ms | 0ms |
剰余リスト生成 | 350ms | 70ms |
リストの乗算 | 15ms | 60ms |
乗算結果を除算 | 90ms | 0ms |
処理内容 | Ver.0.3 | Ver.0.4 |
---|---|---|
二進展開 | 15ms | 0ms |
剰余リスト生成 | 830ms | 890ms |
リストの乗算 | 270ms | 220ms |
乗算結果を除算 | 220ms | 220ms |
フェルマーテストで実際に使う15桁以上の乗数の場合で速度がでていない。 原因は15桁以上の除算がVer.0.3とあまり変わらないことだ。 また、Ver.0.3よりも複雑な実装であるVer.0.4は大規模な乗算では強力な演算能力を発揮するが、小規模な計算を複数回行う際はその複雑性が演算速度低下に直結する。
ようは除算を速くしましょうということだ。
平方根も作るつもりだったが、JavaScriptの精度を越える数値の平方根を計算する方法を考え中。 円周率を計算するには1000桁の数値を平方根に展開しないといけないからだ。
除算関数に新しいタイプを追加。今までは14桁以下か15桁以上の2種類であったが、15桁以上のアルゴリズムがかなりの低速。 そこで新に15桁以上306桁以下で動作する関数を追加。割る数の桁数が小さい順にlight,medium,heavyの3種類の関数を揃えた。 今回作ったmediumは圧倒的な速さではないものの、少なくともheavyよりは高速に動作する。
他にも2点の修正として、剰余を単純変数からオブジェクトに変更しResult.Surplusに格納した。また、商がゼロのとき空白が返るバグを修正。
もう一つ。任意の10進数を2進展開するbinary.jsを追加。繰り返す除算の回数では2^53で割るのが最短だが、割る数が15桁以上の場合は除算が極端に遅い。 そこであえて2^46で割ることによって最も速い除算関数を使って計算することにした。 結果Ver.0.3では400桁の10進数の展開に100ミリ秒かかっていたのが30ミリ秒。1000桁の整数の展開はVer.0.3が6秒に対し1.5秒という劇的な高速化を達成。
あと外部JSを数が多くなってきたのでINCLUDEされているファイル一覧をコンパクトにし、Demosでデモ一覧に戻れるようにした。またリリース日を指定する変数を追加。
不完全ながら除算関数を搭載。これで四則演算が揃った。
除算関数は14桁以下ではJavaScriptの精度範囲内でなので新しいアルゴリズムを採用し旧来のシステムより数十倍も速い。 しかし、15桁以上では"速い"方法が見つからなかったので旧来のアルゴリズムを採用している。 複雑な上に旧式なので15桁以上はかなり遅い。改善方法を考え中。
除算だけで5日以上かかったが、減算を行う際の符合処理が除算の関数にまで影響していたのが問題だった。
旧来のスクリプトは単純変数を用いていたため、単に代入する場合は変数の中身が移動していた。
しかし、今回はオブジェクトで処理しているので単に代入するだけではオブジェクトの中身ではなくオブジェクトのアドレスが代入される。
別名変数でもアドレスが同じだと片方が書き換えられればもう片方も書き換えられてしまう。
このため、HugeInt形式を判定する関数の戻り値をif(input.isHugeInt)return input;
からif(input.isHugeInt)return input.Copy();
に変更し、アドレスのみが代入されないようにした。
おまけ:40000!
乗算関数を追加。2000桁の乗算をVer.0.3alpha9と比べてみた。
計算速度がなんと2倍。やはり乗算結果を配列にしたのが大きい。期待通りの結果が出た。
また、2乗計算の場合10月27日の日記で記述していたように計算量を半分にできる。 今回これを実装したHugeInt.square関数を追加。pow.jsを組んで累乗計算も試した。
結果、Ver.0.3で33秒かかった2の5万乗を10秒で計算できるようになった。乗算の高速化が一番デカいが、square関数でも数秒短縮している。 既に2の1万乗1.3秒の記録を0.6秒にまで短縮した。
階乗(fact.js)を追加。今までの階乗計算は、まず精度限界までは自力乗算しその先は乗算関数を呼び出して計算するものだった。 つまりHugeInt.mulを連続で呼び出す方式。
今回はそうではなく、累乗のときに2乗専用のsquare関数を作ったのと同じように階乗専用に最適化された乗算プログラムを搭載した。 また、前回はn!ごとにまとめて計算して最後に全部掛けていたが、最終的な計算量はそんなに変わらないので今回はベタに積を積み重ねていく方式にした。 乗算方式は例によって配列型。専用の乗算を搭載したので、基礎プログラムを除けば階乗自体はスタンドアーロンで動く。
計算速度であるが、Ver.0.3とは桁違いに速い。1000!ならVer.0.3では1秒かかったが、今は大体半分の0.6秒で計算できる。 1000!は一瞬、2000!は2秒で楽勝。今まではかなりきつかった5000!も15秒程度で計算できる。 ちなみにデモでは10000!までのボタンを用意しているが、10000!は1分ちょっとで計算できる。 今までの最高記録が10000!だったので、今回300秒程度で20000!を計算し結果を掲載。末尾のゼロだけで4999桁ある。
おまけ:2500000。今までの記録は2222222なので一気に更新。目指せ100万乗。
Ver.0.3 betaでいこうと考えていたが、設計があまりに違うためVer.0.4とした。 また、アルゴリズムも純粋文字列型から配列を混ぜた複合型で行いたい。 文字列型アルゴリズムの強みは巨大整数を好きに切り分けられることであり、配列にはこのような柔軟な切り分けができない。 これに対し配列型アルゴリズムでは、変数の数の多さを生かし乗算などを文字列型より圧倒的に速く計算できる。 この2つの特性を使ってさらなる高速化を目指す。
脳内構想としては新しい除算アルゴリズム及び加算関数を全く呼び出さない乗算アルゴリズムがある。 以下、Ver.0.3 beta1 TESTからの修正点。
加算減算のデモが有効です。ランダムな1000桁の整数をインプットでき、1回では計算時間が0msなので100回ループ計算するデモも用意しました。