23

私は、大きな数を除算するアルゴリズムについて考えていました。剰余でbigintCをbigintDで除算します。ここで、基数bでのCの表現がわかっており、Dの形式はb^k-1です。例で示すのがおそらく最も簡単です。C=21979182173をD=999で割ってみましょう。

  • 数字を3桁のセットとして書き込みます：21 979 182 173
  • 左から順に、連続するセットの合計（999を法とする）を取ります：21 001 183 356
  • 「999を超えた」セットの前にあるセットに1を追加します：22 001 183 356

実際、21979182173/999=22001183および残りは356です。

複雑さを計算しました。間違いがなければ、アルゴリズムはO（n）で機能するはずです。ここで、nは基数b表現のCの桁数です。また、C ++で非常に大雑把で最適化されていないバージョンのアルゴリズム（b = 10のみ）を実行し、GMPの一般的な整数除算アルゴリズムに対してテストしましたが、実際にはGMPよりもうまくいくようです。私が見たところ、このような実装されたものはどこにも見つからなかったので、一般的な部門に対してテストすることに頼らざるを得ませんでした。

非常によく似た問題について説明している記事をいくつか見つけましたが、特に2とは異なるベースで、実際の実装に焦点を当てているものはありません。これは、数値が内部に格納される方法によるものと思われますが、前述のアルゴリズムは、たとえば、それを考慮しても、b=10です。私も他の人に連絡しようとしましたが、やはり役に立ちませんでした。

したがって、私の質問は次のようになります。前述のアルゴリズムが説明されている記事や本などがあり、実装について説明している可能性がありますか。そうでない場合、たとえばC / C ++でそのようなアルゴリズムを実装してテストすることは理にかなっていますか、それともこのアルゴリズムは本質的に悪いのでしょうか？

また、私はプログラマーではありません。プログラミングはかなり大丈夫ですが、確かにコンピューターの「内部」についてはあまり知識がありません。したがって、私の無知を許してください-この投稿には1つ以上の非常に愚かなことがある可能性が高いです。もう一度ごめんなさい。

どうもありがとう！

コメント/回答で提起されたポイントのさらなる明確化：

みなさん、ありがとうございました。すばらしい回答やアドバイスをすべて同じものでコメントしたくなかったので、多くの方が触れた1つのポイントについてお話ししたいと思います。

私は、ベース2 ^ nでの作業が、一般的に言って、明らかに最も効率的な方法であることを十分に認識しています。ほとんどすべてのbigintライブラリは2^32などを使用します。しかし、もし（そして、私が強調するのは、この特定のアルゴリズムにのみ役立つでしょう！）、基数bの数字の配列としてbigintsを実装するとどうなるでしょうか？もちろん、ここではbが「合理的」である必要があります。最も自然なケースであるb = 10は、十分に合理的であるように思われます。数値が内部に保存される方法を考慮すると、メモリと時間の両方を考慮すると、多かれ少なかれ非効率的であることはわかっていますが、私の（基本的でおそらく何らかの欠陥のある）テストが正しければ、GMPの一般的な区分よりも速く結果を生成することができました。これは、そのようなアルゴリズムを実装することに意味があります。

Ninefingersは、その場合、高価なモジュロ演算を使用する必要があることに気づきました。古い+新しい+1の桁数を見るだけで、古い+新しいがたとえば999と交差したかどうかを確認できます。4桁の場合は、これで完了です。さらに、old<999およびnew<= 999であるため、old + new + 1が4桁の場合（これ以上持つことはできません）、（old + new）％999は（の左端の桁を削除することになります） old + new + 1）、これは安くできると思います。

もちろん、私はこのアルゴリズムの明らかな制限に異議を唱えていませんし、改善できないと主張していません-それは特定のクラスの数でしか分割できず、ベースbでの配当の表現を事前に知っている必要があります。ただし、たとえばb = 10の場合、後者は自然に見えます。

さて、上で概説したように、bignumを実装したとしましょう。ベースbでC=（a_1a_2 ... a_n）、D = b^k-1と言います。アルゴリズム（おそらくはるかに最適化されている可能性があります）は次のようになります。タイプミスが少ないといいのですが。

  • k> nの場合、明らかに完了です
  • Cの先頭にゼロ（つまりa_0 = 0）を追加します（たとえば、9999を99で除算しようとする場合に備えて）
  • l = n％k （「通常の」整数のmod-高すぎないようにする必要があります）
  • old =（a_0 ... a_l）（最初の数字のセット、おそらくk桁未満）
  • for（i = l + 1; i <n; i = i + k）（floor（n / k）程度の反復があります）
    • new =（a_i ... a_（i + k-1））
    • new = new + old （これはbigintの追加であるため、O（k））
    • aux = new + 1 （繰り返しますが、bigintの追加-O（k）-私は満足していません）
    • auxがk桁を超える場合
      • 補助の最初の桁を削除します
      • old = old + 1 （もう一度bigint加算）
      • 古いものを最初にゼロで埋めて、必要な桁数になるようにします
      • （a_（ik）... a_（i-1））= old （i = l + 1の場合、（a _ 0 ... a _ l）= old）
      • new = aux
    • newを最初にゼロで埋めて、必要な桁数になるようにします
    • （a_i ... a_（i + k-1）= new
  • quot =（a_0 ... a_（n-k + 1））
  • rem = new

そこで、これについて私と話し合ってくれてありがとう-私が言ったように、これは、致命的な欠陥が見当たらない場合に、実装、テスト、および議論を試みる興味深い「特殊なケース」のアルゴリズムのようです。これまで広く議論されていないものであれば、さらに良いでしょう。ご意見をお聞かせください。長い投稿でごめんなさい。

また、もう少し個人的なコメント：

@Ninefingers：私は実際にGMPがどのように機能するか、それが何をするか、そして一般的なbigint除算アルゴリズムについてある程度の（非常に基本的な！）知識を持っているので、あなたの議論の多くを理解することができました。また、GMPが高度に最適化されており、さまざまなプラットフォームに合わせてカスタマイズされていることも承知しています。そのため、一般的に「打ち負かそう」とはしていません。これは、先のとがった棒で戦車を攻撃するのと同じくらい実り多いようです。ただし、これはこのアルゴリズムの考え方ではありません。非常に特殊な場合に機能します（GMPではカバーされていないようです）。無関係なことに、一般的な分割はO（n）で行われていますか？私が見た中で最も多いのはM（n）です。（そして、私が正しく理解していれば、実際には（Schönhage–Strassenなど）O（n）に到達しない可能性があります）。それでもO（n）に到達しないフューラーのアルゴリズムは、私が正しければ、ほぼ純粋です。理論的。）

@Avi Berger：考え方は似ていますが、これは実際には「九去法」とまったく同じではないようです。しかし、私が間違っていなければ、前述のアルゴリズムは常に機能するはずです。

4

3 に答える 3

12

このアルゴリズムは、「キャスト アウト ナイン」として知られる 10 進法アルゴリズムのバリエーションです。あなたの例は、ベース1000と「キャストアウト」999（ベースより1少ない）を使用しています。これは、手計算の簡単なチェックを行う方法として小学校で教えられていました。高校の数学の先生がいて、もう教えられていないことを知ってぞっとし、私たちにそれを記入してくれました。

基数 1000 で 999 をキャストアウトしても、一般的な除算アルゴリズムとしては機能しません。実際の値ではなく、実際の商と余りに 999 を法とする合同の値を生成します。あなたのアルゴリズムは少し異なり、動作するかどうかは確認していませんが、基数 1000 を効果的に使用し、除数が基数より 1 小さいことに基づいています。47 で割るために試してみたい場合は、最初に基数 48 の数値システムに変換する必要があります。

詳細については、Google の「キャスト アウト ナイン」を参照してください。

編集：私はもともとあなたの投稿を少し早く読みすぎましたが、これが機能するアルゴリズムであることを知っています. @Ninefingers と @Karl Bielefeldt がコメントで私よりも明確に述べているように、パフォーマンスの見積もりに含まれていないのは、目前の特定の除数に適したベースへの変換です。

于 2011-02-23T22:54:07.403 に答える
5

私のコメントに基づいて、これに追加する必要があると感じています。これは答えではなく、背景についての説明です。

bignum ライブラリは、リムと呼ばれるものを使用します。gmp ソースで mp_limb_t を検索します。これは通常、固定サイズの整数フィールドです。

足し算のようなことをするとき、(非効率ではありますが) アプローチする 1 つの方法は、次のようにすることです。

doublelimb r = limb_a + limb_b + carryfrompreviousiteration

この 2 倍のサイズのリムは、合計がリムのサイズよりも大きい場合に、limb_a + limb_b のオーバーフローをキャッチします。そのため、四肢のサイズとして uint32_t を使用している場合、合計が 2^32 より大きい場合、オーバーフローをキャッチできます。

なぜ私たちはこれが必要なのですか？さて、あなたが通常行うことは、すべてのリムをループすることです-整数を分割してそれぞれを通過することでこれを自分で行いました-しかし、算術を行うのと同じように、LSL を最初に実行します (したがって、最小のリムを最初に実行します)。手で。

これは非効率に見えるかもしれませんが、これは単に C のやり方です。大きな銃を本当に打破するために、x86 にはadc指示として - キャリーを追加します。これが行うことは、フィールドに対する算術演算であり、算術演算がレジスタのサイズをオーバーフローした場合にキャリー ビットを設定します。次にaddまたはを実行するadcと、プロセッサはキャリー ビットも考慮に入れます。減算では、借用フラグと呼ばれます。

これはシフト操作にも当てはまります。そのため、プロセッサのこの機能は、bignum を高速にするために非常に重要です。実際のところ、チップにはこのようなことを行うための電子回路があります。ソフトウェアで行うと、常に遅くなります。

あまり詳細には触れませんが、操作はこの加算、シフト、減算などの機能から構築されます。これらは非常に重要です。ああ、それを正しく行っている場合は、リムごとにプロセッサのレジスタの全幅を使用します。

2 番目のポイント - ベース間の変換。元の基数でその下の桁からのオーバーフローを説明することはできず、その数値は下の桁からのオーバーフローを説明できないため、数値の途中で値を取得してその基数を変更することはできません。 .. 等々。つまり、基数を変更するたびに、bignum 全体を元の基数から新しい基数に再度変換する必要があります。したがって、ビッグナム (すべての手足) を少なくとも 3 回歩く必要があります。または、代わりに、他のすべての操作でオーバーフローを高価に検出します...プロセッサが私たちのためにそれを行う前に、オーバーフローした場合に解決するためにモジュロ操作を行う必要があることを思い出してください。

また、この場合はおそらく高速ですが、bignum ライブラリ gmp がメモリ管理などのかなりの作業を行うことを覚えておいてください。を使用している場合mpz_は、まず、ここで説明した以上の抽象化を使用しています。最後に、gmp は、これまでに聞いたことのあるほぼすべてのプラットフォームに加えて、展開されたループを使用して手動で最適化されたアセンブリを使用します。Mathematica に同梱されているのには十分な理由がある、Maple など。

さて、参考までに、読み物をいくつか。

  • Modern Computer Arithmeticは、任意精度ライブラリの Knuth に似た作業です。
  • Donald Knuth 著、Seminumerical Algorithms (The Art of Computer Programming Volume II)。
  • bsdntのアルゴリズムの実装に関するWilliam Hart のブログでは、さまざまな除算アルゴリズムについて説明しています。bignum ライブラリに興味がある場合、これは優れたリソースです。この種のものを追い始めるまで、私は自分自身を良いプログラマーだと思っていました...

要約すると、除算のアセンブリ命令は最悪なので、剰余算術で除算を定義するときに行うように、逆数を計算して代わりに乗算するのが一般的です。存在するさまざまな手法 (MCA を参照) は、ほとんどが O(n) です。

編集: OK、すべてのテクニックが O(n) であるとは限りません。div1 と呼ばれるテクニックのほとんどは O(n) です。大きくすると、O(n^2) の複雑さになります。これを避けるのは困難です。

さて、bigints を数字の配列として実装できますか? はい、もちろんできます。ただし、追加の直前のアイデアを検討してください

/* you wouldn't do this just before add, it's just to 
   show you the declaration.
 */
uint32_t* x = malloc(num_limbs*sizeof(uint32_t));
uint32_t* y = malloc(num_limbs*sizeof(uint32_t));
uint32_t* a = malloc(num_limbs*sizeof(uint32_t));
uint32_t m;

for ( i = 0; i < num_limbs; i++ )
{
    m = 0;
    uint64_t t = x[i] + y[i] + m;
    /* now we need to work out if that overflowed at all */
    if ( (t/somebase) >= 1 ) /* expensive division */
    {
        m = t % somebase; /* get the overflow */
    }
}

/* frees somewhere */

これは、スキームを介して追加するために見ているもののラフスケッチです。したがって、塩基間の変換を実行する必要があります。したがって、ベースの表現への変換が必要になり、完了したら元に戻す必要があります。これは、このフォームが他の場所で非常に遅いためです。ここでは O(n) と O(n^2) の違いについて話しているわけではありませんが、リムごとの高価な除算命令または除算するたびに高価な変換について話しているのです。これを参照してください。

次に、一般的なケース部門の部門をどのように拡張しますか? つまり、上記のコードからxとyの 2 つの数値を除算したい場合です。高価な Bignum ベースの施設に頼らないと、答えはありません。クヌートを参照してください。自分のサイズよりも大きいモジュロをとることはうまくいきません。

説明させてください。21979182173 mod 1099 を試してください。ここでは、簡単にするために、最大サイズのフィールドが 3 桁であると仮定します。これは不自然な例ですが、gcc 拡張機能を使用して 128 ビットを使用する場合、私が知っている最大のフィールド サイズです。とにかく、ポイントはあなたです：

21 979 182 173

あなたの数を手足に分割します。次に、モジュロと合計を取ります。

21 1000 1182 1355

うまくいきません。これは Avi が正しい場所です。これは 9 をキャストする形式またはその適応であるためです。ただし、最初はフィールドがオーバーフローしているため、ここでは機能しません。モジュロを使用して、各フィールドが確実に範囲内に収まるようにします。その四肢/フィールドサイズ。

それで、解決策は何ですか？数値を適切なサイズの一連の bignum に分割しますか? bignum 関数を使用して、必要なすべての計算を開始しますか? これは、フィールドを直接操作する既存の方法よりもはるかに遅くなります。

おそらく、bignum ではなく四肢で除算するためのこのケースを提案しているにすぎません。このアルゴリズムがヘンセル除算よりも高速かどうかはわかりません。興味深い比較になります。この問題は、bignum ライブラリ全体で共通の表現に伴います。既存の bignum ライブラリで選択された表現は、私が拡張した理由によるものです。最初に行われたアセンブリ レベルでは意味があります。

補足として; uint32_t手足を表すために使用する必要はありません。アセンブリに最適化されたバージョンを利用できるように、理想的にはシステムのレジスタのサイズ (uint64_t など) のサイズを使用します。そのため、64 ビット システムadc rax, rbxでは、結果が 2^64 ビットを超える場合にのみ、オーバーフロー (CF) を設定します。

tl;drバージョン: 問題はアルゴリズムやアイデアではありません。アルゴリズムに必要な表現は、add/sub/mul などで最も効率的な方法ではないため、基数間の変換の問題です。knuth を言い換えると、これは数学的エレガンスと計算効率の違いを示しています。

于 2011-02-24T00:05:45.330 に答える
0

同じ除数で頻繁に除算する必要がある場合は、それ(またはその累乗) を基数として使用すると、基数 2 の 2 進整数のビット シフトと同じくらい安く除算できます。

必要に応じて、ベース 999 を使用できます。10 の累乗基数を使用することについて特別なことはありませんが、10 進整数への変換が非常に安価になります。(整数全体で完全な除算を行う代わりに、一度に 1 つのリムを処理できます。これは、2 進整数を 10 進数に変換することと、4 ビットごとに 16 進数に変換することの違いのようなものです。2 進数 -> 16 進数で開始できます。ただし、2 のべき乗以外の基数への変換は、除算を使用して LSB ファーストにする必要があります。)

たとえば、パフォーマンス要件のあるコード ゴルフの質問でフィボナッチ (10 9 )の最初の 1000 桁を計算するには、私の 105 バイトの x86 マシン コードの回答で、この Python の回答と同じアルゴリズムを使用しました: 通常のa+=b; b+=aフィボナッチ反復ですが、毎回 10 (の累乗) で割るとa大きくなりすぎます。

フィボナッチは、キャ​​リーが伝播するよりも速く成長するため、10 進数の下位桁を破棄しても、長期的には上位桁が変わらないことがあります。（必要な精度を超えていくつか余分に保持します）。

最終的な 2 進数 -> 10 進数の変換はそれに依存するため、破棄した 2 の累乗の数を追跡しない限り、2 の累乗で除算することはできません。

したがって、このアルゴリズムでは、拡張精度の加算と 10 (または任意の 10 の累乗) による除算を行う必要があります。

基数 10 の9リムを 32 ビット整数要素に格納しました。10 9で除算するのは簡単です。下肢をスキップするためのポインタのインクリメントだけです。実際に を実行する代わりにmemmove、次の追加反復で使用されるポインターをオフセットするだけです。

10^9 以外の 10 のべき乗による除算はやや安上がりだと思いますが、各肢で実際の除算を行い、残りを次の肢に伝播する必要があります。

この方法では拡張精度の加算はバイナリ リムよりもいくらかコストがかかります。これは、compare: sum[i] = a[i] + b[i]; carry = sum < a;(符号なし比較) を使用してキャリーアウトを手動で生成する必要があるためです。また、条件付き移動命令を使用して、その比較に基づいて手動で 10^9 にラップします。しかし、そのキャリーアウトをadc(x86 add-with-carry 命令) への入力として使用することができました。

最大で 1 回ラップしたことがわかっているため、加算時にラップを処理するために完全なモジュロは必要ありません。

これは、各 32 ビット リムの 2 ビット強を無駄にします: 代わりに 10^9 です2^32 = 4.29... * 10^9。基数 10 の数字を 1 バイトに 1 つずつ格納すると、スペース効率が大幅に低下し、パフォーマンスが大幅に低下します。これは、最新の 64 ビット CPU での 8 ビット バイナリ加算のコストが 64 ビット バイナリ加算のコストと同じであるためです。

私はコードサイズを目指していました: 純粋なパフォーマンスのために、基数 10^19 の「数字」を保持する 64 ビットリムを使用していました。( 、したがって、これは 64 ごとに 1 ビット未満を無駄にします。) これにより、各ハードウェア命令2^64 = 1.84... * 10^19で 2 倍の作業を行うことができます。addうーん、実際にはこれが問題になる可能性があります。2 つのリムの合計が 64 ビット整数をラップする可能性があるため、チェックするだけで> 10^19はもはや十分ではありません。5*10^18baseまたは base で作業する10^18か、バイナリキャリーと手動キャリーをチェックするより複雑なキャリーアウト検出を行うことができます。

1 バイト内で 1 つのニブルから次のニブルへのキャリーをブロックするためのハードウェア サポートがないため、4 ビット ニブルごとに 1 桁でパックされた BCD を格納すると、パフォーマンスがさらに低下します。

全体として、私のバージョンは、同じハードウェア上で Python 拡張精度バージョンよりも約 10 倍高速に実行されました (ただし、分割の頻度を減らすことで、速度を大幅に最適化する余地がありました)。(70 秒または 80 秒対 12 分)

それでも、そのアルゴリズムのこの特定の実装 (足し算と割り算だけが必要で、いくつかの足し算のたびに割り算が行われる) では、基数 10^9 リムの選択は非常に良かったと思います。N 番目のフィボナッチ数に対して、10 億回の拡張精度の加算を行う必要のない、はるかに効率的なアルゴリズムがあります。

于 2017-11-16T21:41:11.907 に答える