BigInteger が使えると簡単なんですけどね・・。
独自で BigInteger の割り算を実装しようかと思ったのですが、
分数の約分で実装してみました。
int の配列を用意して 2 * 3 * 4 * 5 を
int[] a = {0,0,1,1,1,1}
と表します。
このa を 3 倍する時は a[3] = a[3] + 1;
5倍する時は a[5] = a[5] + 1;
と表します。
/** 掛け算 */
void mul(int ar[], int n) {
if (n > 0 && n < MAX_AR) {
ar[n] ++;
}
}
また 約分をするための関数を作ります。
64 は
a[64] = 1;
で表し
これを 約分すると 2 * 2 * 2 * 2 * 2 * 2 なので 2 を 6 回 掛ける事なので
a[2] = 6 と表します。
void yakubun(int x[]) {
int prime[] = { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47 };
for (int i = 0; i < 15 ; i++) {
int pi = prime[i];
for (int j = (MAX_AR-1)/i ; j > 1; j --) {
int ij = i * j;
if (x[ij] > 0) {
// 例:3の5倍 (15) に 2 が 入っている場合 x[3] と x[5] に 2 を足す;
x[i] += x[ij];
x[j] += x[ij];
x[ij] = 0;
}
}
}
}
次に 割り算を作ります。
分数の約分を使って、
必ず割り切れるという前提での実装です。
汎用的には使えません。
/** 2つの数字の割り算 */
long long div(int b[], int a[]) {
int c[MAX_AR];
for (int i = 0; i < MAX_AR; i++) {
c[i] = 0;
}
// b / a を実施する 分母と分子で同じ物を約分する
// 配列に入っている数字を約分する
yakubun(a);
yakubun(b);
for (int i = 0; i < MAX_AR; i++) {
if (b[i] >= a[i]) {
c[i] = b[i] - a[i];
}
}
// 最終的な計算結果
long long s = 1;
for (int i = 1; i < MAX_AR; i++) {
if (c[i] > 0) {
for (int j = 0; j < c[i]; j++) {
s = s * i;
// 掛け算の余りを使うとオーバーフローせず 正しい答えが出せる
if (s >= 1000000000) {
s = s % 1000000000;
}
}
}
}
return s;
}
最後に combi は 再帰を使わず 掛け算と割り算に展開します。
nCr = ( n から 数を下げながら r 個の掛け算) / (1 から 数を上げながら r 個の掛け算)
long long combi(long l, long r) {
int a[MAX_AR];
int b[MAX_AR];
// 変数初期化
for (int j = 0; j < MAX_AR; j++) {
a[j] = 0;
b[j] = 0;
}
for (int i = 0; i < r; i++) {
// 1 から r までの掛け算
mul(a, i + 1); // a = a * (i + 1);
// l から r + 1 までの掛け算
mul(b, l - i); // b = b * (64 - i);
}
return div(b, a);
}
main 関数と、定数宣言等
#define MAX_AR 65
int main()
{
char ans[32];
char last9[32];
sprintf_s(ans, "%lld", combi(64, 32));
strcpy_s(last9 , &ans[strlen(ans) - 9]);
printf("%s\n", last9);
return 0;
}
割り算を 手本なしに作る自信がなかったので 約分という方法を使ってみました。
かける数も 64 までの 小さな数という事が前提なので 汎用的な BigInteger の掛け算
を作っていません。
C# で BigInteger を使えるなら
// 64C32 を計算する
BigInteger a = 1, b = 1;
for(int i=0;i<32;i++)
{
a = a * ( i + 1);
b = b * (64 - i);
}
BigInteger c = b / a;
たったこれだけです。
482219923991114978843459072919892677776312893440000000 / 263130836933693530167218012160000000 = 1832624140942590534
となります。
上記プログラムで MAX_AR を大きくして それに合わせた素数を計算する事で
かなり大きな数でも 高速に計算できることを確認しました。
combi(60, 28)/(5*29)*(2*21*61)
ま、参考までにどうぞ。