著作権フリーのC++高速汎用多倍長演算を作るｽﾚ

まずは高速フーリエ変換を実装する
ネット上でここが詳しい
http://homepage1.nifty.com/~umetani/ims/Algebra.html

てめぇ一人でやれ

日本じゃ著作権フリーは無理だから、外人に頼むしかないな。

ここに書き込んだ人が著作権を主張しないで自由に使って良いっていえばいいだけではないの?

作るスレじゃなくて作ってくださいスレの間違いだろ。
単発依頼スレを自分だけ特別に立てて良いと思った理由を述べろ。

ちゃんと活動するよ　まずは離散フーリエ変換のノートを書き込むからまっててくれ

やっぱり面倒になった

作るやつはスレ立てる前に作り始める
できんやつはスレ立てる．
若干反例があるが，これは普遍の真理．

なぜかというと理論的な部分は簡単なのだが、高速化がする部分が複雑で既にあるやつ使った方が良いと思うからだ

著作権フリーなものを海外から日本へ持ち込んだ場合、
その著作権は誰のものになるんだ？

まさか、FFTをやってくれるフリーなライブラリが欲しいだけってことはないよね？
ちなみに、それはもう存在してて、ググれば見つかるよ。
商用利用もOK。GPLじゃない。

Ｒ　環に対して
Ｒ［Ｘ］/（Ｘ＾ｎ-1）　→　R[G]　（同型のフーリエ変換）

という写像を作って
多項式環の乗法が、　右側だと簡単で計算出来て、
それを同型で戻すという寸法だ

>>11
では、高速離散フーリエ変換のライブラリでフリーで高速にうごくやつを組み込んで作ろう
どこかにいいやつありますか?

ググレカス

>>14
それより厳密な整数演算に使える高速フーリエ変換ってあるんですか?
誤差が出るやつでは駄目です　整数演算のみでやるやつです

多倍長じゃなかったのかよ

多倍長って2バイトや4バイトや8バイト以上の整数、小数を扱うって事だよ

つくってみるわ
とりあえず
基幹部分はC
インターフェイスはPython
でいいんだよな

ttp://homepage1.nifty.com/~umetani/ims/Tools.html

FFT乗算の原理


数値を多項式と見立てて、多項式の乗算、即ち、離散系の畳み込みをフーリエ空間で行い、
それを実空間に戻して、数値の乗算をしたことにするのである。
今、4096個の一次元データに4096個データの一次元フィルタを掛けるとすれば、
実空間では、4096*4096 約1678万回の乗算が必要となるが、
DFTでは、畳み込みだけでは、0.4万回で済む。前後に3回のDFTが必要となるが、
この分は、FFTで行えば8万回程度なので、約8.4万回で済むことになる。約200倍のスピードとなる。
http://www.nextftp.com/swlabo/m0_pctech/hp_ultraprecision/up_815_1.htm

高速フーリエ変換にかかるコストを除けば、Ｎ*Ｎの乗法がＮに出来る最強アルゴリズムだよな
3/4　Ｎ*Ｎにできるアルゴリズムがあるが断然違うよな

lispでも使ってれば？

せめてメール蘭にsageって入れようぜ。

このスレはFFTを実装するスレだったのか。
多倍長整数の加減乗除を実装するスレなのかと思った。

スレタイに汎用って入れているくらいだしね。

>>9
じゃあ低速でいいからとっとと作れよ、簡単なんだろ。

高速化するのはその後でいいよ。
最低限、自分の発言に責任を持て。

FFTが多倍長の乗算に使われる事を知らない奴ばかりで、
1がDQNあつかいされてるな。
フーリエ変換って、多倍長乗算やCTとか、
周波数とか関係のない、意外なとこで利用されてるから面白いよね。

違うわ。αでもβでもいいから物を出さないのがいけない。

>>4亀レスだが
2ちゃんねるに書き込んだものってひろゆきが権利持つんじゃないだろうか。

離散フーリエ変換の定義

可換環Ｒに対して、次を満たす要素ｅが存在するとする　
・　ｅは1の原始ｎ+1乗根　（すなわち、ｅのｎ+1乗は1かつ、ｎ乗以下では1ではない）
・　ｋがｎ以下ならば、　Σ( i, 0, n) E(k) ^ i =0 ただしE(k) ＝ｅ＾ｋ

このとき、Ｇ＝｛1，Ｅ（1），・・・，Ｅ（ｎ）｝は群をなす
ＧからＲへの写像全体をＣ（ｅ）と置くと、次のようにしてR上の有限次元ベクトル空間になる
内積を、(f,g)=Σ( i, 0, n) f( Ｅ（j）) g( Ｅ（-j）)で定義すると、｛λi｝が正規直交基底を与える　（ただしλi(　Ｅ（j）　)＝Ｅ（jj））
さらにＣ（ｅ）は次で与える積でR代数になる　(f*g)( Ｅ（k）) = Σ( i, 0, n) f( Ｅ（k-i）) g( Ｅ（i）)

R[X]の要素fに対して、C(e)の要素Ｆ（ｆ）を次のように定義する　Ｆ（ｆ）( Ｅ（k）) = Σ( i, 0, n) f( Ｅ（i）) Ｅ（ki）
この対応はR代数として準同型写像を与える
これから導かれる　同型　R[X]/(X^(n+1) - 1) ≡　C(e)　を離散フーリエ変換という

たとえば、Ｒ＝Ｚ/（5）＝｛0，1，2，3，4｝において、2は1の原始4乗根であり、もう一方も条件を満たすので
Ｒ［Ｘ］/（Ｘ^4-1）　からC(2)への離散フーリエ変換が定義できる

フーリエ変換を使って高速になるのは桁数がある程度以上大きくなってからでしょう。
だからまずフーリエ変換を使わないプログラムを作るべきでは？

>>32　まあそういうことですね　どこまでフリーリエ以外が速いか知らないといけないですからね
しかしフーリエで適用できるようにデータ構造は決めないと駄目なんです

離散フーリエ変換のまとめ

可換環Ｒに対して、次を満たす要素ｅが存在するとする　
・　ｅは1の原始ｎ乗根
・　ｋがｎ未満ならば、　Σ( i, 0, n-1) E(k) ^ i =0 ただしE(k) ＝ｅ＾ｋ

このときR[X]/(x^n - 1)の自己同型写像F

f 　→　Σ( i, 0, n-1)　f( E(i) ) X^i をフーリエ変換という
逆写像は、f 　→　(1/n) Σ( i, 0, n-1)　f( E(n-i) ) X^i で与えられる

このとき、畳み込み積　　(f*g)( Ｅ（k）) = Σ( i, 0, n-1) f( Ｅ（k-i）) g( Ｅ（i）)に対して　F(fg) = F(f) * F(g)が成り立つ
X^s * X^t = X^s (s=t) or 0 (それ以外)を満たす　この性質から畳み込み積の計算は簡単にできる

34は途中で30で定義していた同型R[X]/(X^n - 1) ≡　C(e)を同一視している部分があります

まず必要なことは、素数pとZ/(p)を乗法群と見たときの最大位数の元を決定することです
最大位数が大きいほど、計算できる桁数が上がります　しかし、pが大きすぎると計算時間がかかります
どうやって最大位数の大きい元をみつけたらいいんでしょうか?

p = 2^16 - 15に対して、17はZ/(p)において位数p-1の要素だと判明したのでこれを使うことにします

素数と位数を調べるコード

#include <iostream>
using namespace std;

main(){
int n,k,N,cnt;
for(N=(1<<16)-1;N>2;N--){
for(n=2;n*n<N;n++){k=N/n; if(k*n==N)break;}
if(n*n<=N)continue;
for(n=2;n<=N/2;n++){
k=n;
for(cnt=1;k!=1;cnt++){k*=n;if(k>N)for(;k>N;k-=N);}
if(cnt+1==N){cout<<"素数 "<<N<<"に対し要素"<<n<<"の位数は"<<cnt<<endl;break;}
}}}

変更してp=2^16 + 1 と位数pの要素3を使うことにします

FFT乗法の設計の概要を書いておきます
p=2^16 + 1 とe=3とおく

f(X) = Σa_i ・X^i　(a_i　< 2^16) を多倍長整数とします

これをフーリエ変換すると、 Σf( e^i ) ・X^i　になります

g(X)も多倍長整数とするとフーリエ変換はΣg( e^i ) ・X^i　になります

fとgの畳み込み積は、 Σf( e^i ) g( e^i ) ・X^i　となり

これをフーリエ逆変換すると、Σf( e^(p-1-i) ) g( e^((p-1-i)) ) ・X^i　になります

この係数のうち、2^16であるものは桁上がりを処理して結果が出ます

簡単に言うと、多倍長整数f, gに対し、
f(1)g(1) 、　f(3)g(3)、　f(3^2)g(3^2)、・・・をmod pの範囲で求めれば良いだけです

>>29
もし本当にそうなら47氏逮捕のときにひろゆきも逮捕されてるはず

まちがえました　フーリエ逆変換は、21846・Σf( 3^(p-1-i) ) g( 3^((p-1-i)) ) ・X^i　　です
Z/(p)での演算から　1/3　＝21846

>>42
馬鹿？
金子がWinnyの実装を2ちゃんに書き込んだかよ？

高速フーリエ変換をメモリも再帰も行列計算も使わない方法を解明した
それは、2進数に対応させるんだ　例えば1024*1024のサイズのフーリエ変換ならば
0〜1023の2進数ごとに一定の手続きをすることで計算が完了する

すこしはメモリは使うけどね　全ディレクトリの探索みたいな感じで上位ディレトリのパスだけは保存しておくような感じ

たとえば、
サイズは2Nの、(Xi) = (ω^ij) (xi)というフーリエ変換は、

(X2i) = (ω^2ij) (xi+xN+i)
(X2i+1) = (ω^2ij+j) (xi+xN+iω^N)
に変換できる
これを繰り返したとき、Xの添え字、ωの指数、xの足し算の仕方を
000011101などの2進数から決定できる

(^ω^i)

計算量を考えてみたけど、フーリエ変換にN*Nかかって、これを2回やって畳み込み積がN回だから
単純な総当たり掛け算より多いような気がする　特典は桁上がりの死よりが少ないだけでは?

通常のフーリエ変換だとN*Nが高速だと、 N* log Nだね
Karatsuba法だとN^1.59だね　やはり高速フーリエ変換が最速か
しかし、より詳細に求めると、 2N* log N　+　N　だから、場合によってはN^1.59のほうが速いかもしれないね

良くわからんけど、FFTWじゃ駄目なの？

>>51
使い方分からないし、整数演算に使えるか分からないので自作してみるよ・・・
あと、自作使用した理由が、多倍長ライブラリのコンパイルが出来なかった為だったのだ
C++言語のみで、他のライブラリに依存しないで単体で動くと導入が簡単なのだ

Vector stringのようにクラスの動的確保確保ってどうやるのかわかりますか?

たとえば
Lint x;
x=2;
for(int i; i<1000; i++)x=x*x;
x.printf();

としたとき　xのサイズが自動で拡張されるようにしたいのですが

std::vector使えば(ry

それで性能上困ることがあれば、
自分でmalloc/free呼んだり、OSのメモリ確保APIを直接呼んだり、
それをメモリプールしてみたり。

サンクス
vectorは速度的に速いと思うので、内部でvector使うことにします　自分で領域を動的に管理するとバグでやすいですし

すみません

a = x + y という演算子を定義するとき、x+y演算と、代入演算を分離するとしない場合より代入分
計算コストがかかると思うのですが、tasu(&a, &x, &y)という関数と同じように定義することは出来ないですか?
これなら計算結果を直接格納できるのですが

そこらへんは最適化でうまくやってくれるのでは？

つ「式テンプレート」

>>58
サンクス　でも使い方分からなかったよ　
なんか代入にかかる計算程度をケチっても
意味がある計算には変化なさそうなので知ってるやり方でやってみる
わざと、代入が多くなるような計算すれば違いが出るとはおもうけど

なんか足し算作るだけでも面倒
符号と数値を保存してあるとして、足し算が負になるか正になるかで計算を分類しなくていけない

計算コストがかかるけど、演算結果の符号は考慮しないで最後に調べれば統一的に処理できるけど・・・

このアホはリバーシ１かな

足し算まで出来た
http://kansai2channeler.hp.infoseek.co.jp/cgi-bin/joyful/img/5338.txt

>>63
コメントとインデントつけてくれるとうれしいな

class Lint
int fugou; //正負
vector< unsigned short int > a;　//一区切りは2バイトのunsigned整数

Lint& operator=(Lint x);　// Lint の代入

Lint& operator=(int x); //整数の代入

int operator>(Lint& x);
int operator<(Lint& x); //絶対値の比較

Lint operator+(Lint& x)　//足し算　
絶対値の大きい方を先頭にして、残りの数値を引いたり足したりする

除算をどうするかが問題ではある
x / y = x * (1/y) だから1/yを求められればよい
すると、符号の他に小数(の位置)を記録する必要がある
ニュートン法で1/yを計算するのが良いらしいけど
有効桁数までもとめ続けた方が良いのだろうか

やっぱりそうだ

どうでもいいけどぉ〜符号はせめてぇ〜singじゃないのかぁ〜
あとFFTを理解するの面倒でスレを読んでないんだけど多倍長整数みたいなもんだよね？
Lint& operator=(int x)を見る限り4294967296までしか行けそうにないけど、俺が馬鹿なだけ？

＞どうでもいいけどぉ〜符号はせめてぇ〜singじゃないのかぁ〜
歌ってるのか？

ああー、signだった・・・
どうでもいいからやんわりと書こうとしただけで洒落ようとかこれっぽっちも思ってなかったのにorz

どうでもいいけどぉ〜♪符号はせめてぇ〜♪singじゃないのかぁ〜♪

(´；ω；`)

>>52のような事を言ってるようなレベルの>>1じゃ、
FFTWとかFFTEより速いもんが出来上がる訳もないわな。

ついでにその整数環FFT使った多倍長乗算は普通のFFTの物より遅いでFA出てた筈。
プロセッサの設計・製造からやるんなら話が違うかもしれんがね。

>>73
しーっ、リバーシスレのようになるから触っちゃダメ。

ふつうのFFTって何? 　Exp(iθ)とか出てくるやつ? 整数演算しようとするのに浮動小数点や複素数演算したら誤差が出る

>>74 それもそうだな。要らん事言って済まなかった。

引数はLintやLint&よりもconst参照だろ、常考。

ところでoperator +=は？
そっち作った後、それを使って残りを定義するのが普通。
More effective C++ 22項

Lint operator +(Lint const& x, Lint const& y)
{
    return Lint(x) += y;
}

あと、operator<はreturn x > *this;で十分なはず。
さらに!(*this < x)でoperator >=、!(*this > x)でoperator <=も実装できる。
ついでにoperator ==とoperator !=も作ってしまえ。

単項operator +と-も今すぐできるはずだ。

欲を言えばswapが欲しい。

>>74
やっぱりリバーシスレと同じ匂い感じた人がいたかｗ

これから修正します
あと逆数の計算方法の書いてあるサイトをのせます

1/n を求めるには

http://www.nextftp.com/swlabo/m0_pctech/hp_ultraprecision/up_820.htm
http://yowaken.dip.jp/tdiary/20070428.html

>>77
引数1つにしてくれというエラーが出ます

Lint operator +(Lint const& x, Lint const& y) { return Lint(x) += y; }

ここまで自分でなんとかしようと思う1は偉いと思う
応援するよ

あと、+=を定義して、他のメンバ関数で使おうとすると、未定義というエラーを出します
+は他で使えます　+を定義してから　+=を定義する方向で行こうと思います

俺もoperatorのオーバーライドの時の定義の仕方が曖昧なんだけど
うまくまとめてるサイトってないものかね？
http://www.nitoyon.com/vc/tips/multi_precision.htm
こんな感じのサンプル上がってるとこ見て、見よう見まねで実装してるんだけど・・・

なんかconst付けたり付けなかったりで、未定義になったりします　とりあえず動く方法でやっていきます

ここまでできた　最新版　掛け算や割り算はまだ
http://kansai2channeler.hp.infoseek.co.jp/cgi-bin/joyful/img/5341.txt

インクリメントはこれでは駄目なんですかね? 不正な構造体の演算ってでます

opPostInc() {
int i,s=a.size();
for(i=0;i<s;i++){
if(a[i]==65535)a[i]=0;
else {a[i]++;break;}}
if(i==s){a.resize(s+1);a[s]=1;}}

上はD言語の宣言でした　　普通にvoid operator ++() でした

単純な掛け算のコード書いたけど、10進表示がまだ出来ていなくてあっているか分からない

Lint operator*(Lint & x){
Lint *p,*q;
if(*this > x){p=&x; q=this;} else {p=this; q=&x;}
int s,t,u;
s=(*p).a.size(); t=(*q).a.size(); u=s*t;
Lint y; y.a.resize(u+1,0);
if(sign + x.sign==0) y.sign=-1;
int i,n; unsigned w,k=0;
for(n=0; n < u; n++){
w=k; k=0;
for(i=0; i<=min(n,s); i++){
w += (*p).a[i] * (*q).a[n-i];
k+=(w>>16); w&=LIMIT;}
y.a[n]=w;}
if(k>0)y.a[u]=k; else y.a.resize(u,0);
return y;}

a[i]==65535
ここを
a[i]==999999999
にすれば？

>>89
インクリメントの部分ですか?　それはコードの中身ではなくて初めの宣言が違っていました　
opPostInc() という名前はD言語用のインクリメントの宣言でした

2^16進数を10進表示にするにはどうしたらいいんだろう? 10で割ってその余りを出力していけばいいけど10で割るのが難しい

ダメだこりゃ。

いい方法が分かりましたよ
小数点を導入して、値を1以下になるようにシフトさせてから10倍していきます
整数になった部分が最上位の10進表示になります　それを0にして繰り返します

int x=425980;
cout<< (x&65535)<<" "<<(x%65536);

この正解は10のはずなんですけど値が違います　なぜでしょうか

すみません　間違えました
int xx=4259850;
でした

>>92は間違えました　1より小さくするときに10^(-n)で割ってシフトしなければいけなかったです　結局割り算がいりますかね?　

最初から１０進でどうぞ

>>95
10000進の掛け算と足し算でできるのでは？

>>96 ビットシフトや&が出来なくて駄目です
>>97 詳細を教えて下さい　やっぱり割り算はいると思うのですが

xを表示するなら　x%10を出力して、x=x/10として繰り返すという方法です

10進2進変換は2進多倍長の掛け算と足し算でできるのと同じように
2進10進変換が10進多倍長の掛け算と足し算でできるということです。

>>80
そのopeartor +は非メンバ関数。

swapはそれだと意味がない。代入より効率よく実装できる可能性があるからみんな定義する。
例えばstd::vectorなら、交換は内部のポインタを取り替えるだけで済む。
（代入ならメモリのコピーをしないといけない）
そしてswapという名前のメンバ関数にするのがお約束。
void swap(Lint& x) {
a.swap(x.a);
std::swap(sign, a.sign); //<algorithm>
}
そしてグローバルにもswapを定義しておく。
void swap(Lint& x, Lint& y) {
x.swap(y);
}
単項operator+と-はこういうメンバ関数でいい。
Lint opreotar +() const {
return *this;
}
Lint opreotar -() const {
Lint x = *this;
x.sign = -1;
return x;
}
っていうかコピーコンストラクタLint(Lint const& x)作れ。対称性からLint(int x)も欲しい。

>>98
こういう整数を十進法として出力するにはそれくらいしかないと思う。

コーディングテクニックと、アルゴリズムの話が混在してて、読む気にならん・・・

>>86
インクリメントはとりあえず*this + 1でいいよ。

Lint& operator ++() {
return *this += 1;
}

Lint operator ++(int) {
Lint t = *this;
return ++t;
} //あるいはreturn ++Lint(*this);とすれば1行にまとめられる。

とまあ、C++の演算子多重定義は面倒くさいので、
こんな余計なことに頭を使わないためにも、Boost.Operators推奨。
http://www.kmonos.net/alang/boost/classes/operators.html
boost::operators<Lint>とboost::operators<Lint, int>を継承し、
そこのsampleに載っているのと同じ演算子多重定義関数だけを実装すれば、残りは全て揃う。
ただし引数取るものはconst Lint&版とint版の2つを実装すること。

boost使うと、ライブラリの価値が50%は下がるな

俺はそう思わないけど、そう思う人間が多いなら、後で自分で書き直せばいい。
今はとりあえず動く、使えるようすることに注力してほしい。

>>102
<と==を定義して残りの> != <= >=をboostに任せるってやっているけど
これは問題だと思うね。
論理的には正しいが下手をすると処理時間が倍以上かかってしまいそう。

そこら辺は大抵最適化されるよバカ

>>105
>>77にあるようにオペランドの左右を入れ替えたり、結果を!で論理反転するだけ。
コンパイル時の最適化で消えたり、残ってもせいぜい数命令もあるかどうか。
だから無視できる範囲内。

それよりも>>85の現在のoperator ==（と!=）がよっぽどまずい。
!( (*this < x) || (*this > x) )と2度も比較演算を行っており、
これこそ処理時間が倍かかる。

>>107
!( (*this < x ) | ( *this > x ) )これって最適化で、( *this >= x && *this <= x )とかにしてくれるのかな？
まぁ、倍とはいかないけど大きな期待ができるほどじゃないけど・・・

boost前提のライブラリはどうかと思うけど、楽に書き直せる範囲なら
機能実装から先にやっていけばいいと思うな。
比較に処理が2倍かかろうが、どうせ後でプロファイルするんだろうし。

>>108
うわ、気付かなかった。
そういえばoperator <って符号を見ていないよな。

a<=bをa<bとa==bからどうやって構成するのかな？
素直に考えると(a<b)||(a==b)になると思うが
それだと最悪2倍程度時間がかかるだろう

あっ、そうか!(b<a)か。
勘違いした。

>>111
どんまいｗ
operator >の一つだけで一応、全部（>= < <= == !=）実装できるよ。

<<112
でもそれだと処理時間が...
自分が作るとしたらまずCompare関数を作り
残りはそれを呼び出すようにするね。
それなら処理時間は遅くならないはず。

>>113
処理時間は置いといてだから一応ね。
Compare作るくらいなら==作って、そのNOTで!=でいいんじゃないか？

>>113
流石に==は実装しないと処理時間がかかる(>>107-108)けど、
比較比較の> >= < <=はどれか1つあれば、残りは速度低下なく実装できる(>>77>>107)。

例えば、a < bが定義されていたとして
>>110で挙げられたa <= bは、!(b < a)と定義すればよい。
(a<b)||(a==b)は悪いけど率直というより愚直。

Boost.Operatorsが比較演算に関してoperator <と==だけを要求しているのは
伊達にやっているからじゃないんだよ。

a<b <=> Compare(a,b)<0
a==b <=> Compare(a,b)==0
だからCompareを作れば基本的に１つですむということ。
a<bもa==bも実際は内部で同じような処理をするわけだから
同じようなものを2つ作る必要はないでしょう。

そっちはわかる。
つい遅くなるっていうほうにだけ反応していた。

そろそろoperatorの話は置いといて>>1を待とうぜ。
>>95から出てきてないことだし。

フーリエ変換は次の値を出来るだけ高速に求めることと同じです　簡易ルーチン作ってみましたがより速くできますか?

N=2^16 、p=N+1、f(X)= a0 +　・・・　+ a(N-1) ・ X^(N-1)　に対して、f(1)、f(3)、・・・、f( 3^(N-1) ) mod p　を求める

#include <iostream>
#include <vector>
#include <time.h>
#define p 65537
unsigned int mod_p(unsigned int& x){
int y=(x&65535)-(x>>16);
if(y<0)y+=p; return y;}

main(){
int N=1<<13, i,j;
cout<<"初期データ生成中・・・\n";
vector< unsigned int > a, x; a.resize(N); x.resize(N);
for(i=0;i<N;i++)a[i]=((rand()&255)<<8)+(rand()&255);

vector< vector< unsigned int > > e;
for(i=0;i<N;i++)e.resize(N); for(i=0;i<N;i++)e[i].resize(N);
j=1; for(i=0;i<N;i++){ e[0][i]=1; e[i][0]=1; e[1][i]=j; j*=3;}
for(i=1;i<N;i++) for(j=1;j<N;j++) e[i][j] = e[1][(i*j)&(N-1)];

cout<<"計算開始・・・ ";
int cl=clock();
for(i=0;i<N;i++) for(j=0;j<N;j++){
x[i] += a[j]*e[i][j]; x[i]=mod_p(x[i]);}
cl=clock()-cl; cout<<cl<<" m秒\n";}

間違えました　Nは任意で、pは2^16+1でした

いま思ったんですが、どうやら高速フーリエ変換のアルゴリズムは適用出来ないようです
p=2^16+1に対して、p-2桁の多項式を扱う場合は良いんですが、
実用的なサイズは高々1000桁くらい　少なくとも2^15以上の桁で無いと効果はなさそうです

pの取り方と計算可能な桁数は対応してますが、例えば1000桁どおしの掛け算ならば
2000に近い素数を見つけて2000進数に変換してから高速フーリエ変換を適用する方法がありますが
前後の基数変換と、フーリエ変換にかかるコストを考えると適用出来る桁は限られそうです
単純に計算した方が速くなるケースが多そうなため

アルゴリズムの問題ではないけど、
for(i=0;i<N;i++)e.resize(N);ってこれをループにする必要ない。

>>123　ｻﾝｸｽ　そうでした
まずは、この計算を高速にするコードを作ることにします　うちのパソコンは鈍いので47秒掛かります
高速フーリエ変換が正しく実装出来れば、最大　(2^15*2^15) 　/　 (2^16 * 16 * 2 + 2^16) =496倍高速になるはずです
0.1秒以内で処理が終わるはずですが・・・
　
#include <iostream>
#include <vector>
#include <time.h>
#define rnd() (((rand()&255)<<8)+(rand()&255));

main(){
int N=1<<15, i,j,n;
vector< unsigned int > a, b, c;
a.resize(N); b.resize(N); c.resize(2*N);
for(i=0;i<N;i++){a[i]=rnd();b[i]=rnd();}

unsigned int k=0, x, cl=clock();
for(n=0; n<2*N; n++){ x=k; k=0;
for(i=0; i<N; i++){x += a[i]*b[n-i]; k+=x>>16; x&=65535;}
c[n]=x;}
cl=clock()-cl; cout<<N<<"桁どおしの掛け算の処理時間 "<<(cl+0.0)/1000<<" 秒\n";}

地の文だけじゃなく、コードにもホワイトスペース入れろよ・・・

レベル下がりすぎｗ

アルゴリズムのことはよくわかんねぇんだけどさ
x += a[i]*b[n-i];
bの大きさはN個で、n = 0 -> 2N. i = 0 -> Nなんだから
n - iのとる範囲は0 -> 2Nで範囲外アクセスじゃないの？

>>127 修正しました あとminGWとGMPいれて速度比較することにしました　コンパイルに時間が掛かります
#include <iostream>
#include <vector>
#include <time.h>
#define rnd() (((rand()&255)<<8)+(rand()&255));

main(){
int N=1<<15, i, n;
vector< unsigned int > a, b, c;
a.resize(N); b.resize(N); c.resize(2*N);
for(i=0;i<N;i++){a[i]=rnd();b[i]=rnd();}

unsigned int k=0, x, cl=clock();
for(n=0; n<=2*N-2; n++){ x=k; k=0;
for(i=max(0,n-N+1); i<=min(N-1,n); i++)
{x+=a[i]*b[n-i]; k+=x>>16; x&=65535;}
c[n]=x;}
cl=clock()-cl; cout<<N<<"桁どおしの掛け算の処理時間 "<<(cl+0.0)/1000<<" 秒\n";}

GMP速すぎかも・・・　50000！の計算が2秒しかかからない

#include<iostream>
#include<gmpxx.h>
#include<time.h>
using namespace std;

mpz_class factorial(int n){
mpz_class x= 1;
for(int i=1; i<=n; i++) x*= i;
return x;}

int main(void){
int c=clock();
mpz_class x = factorial(50000);
c=clock()-c;
// cout << x.get_str() << endl;
cout<< (c+0.0)/1000<<" sec \n";
scanf("%d",&c);}

>>128
minって関数？もしくはマクロ？
毎回計算させるよりローカルに入れたほうが早くならんかな？
最適化されてもよさそうなところだけど一応。

>>130
#include<iostream>に入っている関数だと思われ　
なんかGMPは巨大な整数には対応していないのかも・・・　65536進数の65535桁とかは扱えないかも・・　計算途中で落ちる

>>131
たしかmin関数はtemplateでしょ。
ローカルに落としたほうが早いよ。（最適化されないのなら）

65 536 * log(65 536) = 315 652.8　だから、10進で30万桁だからたいしたことはないはずだけど・・
円周率などは100億桁とか計算するはずだから
>>132 上のコードは高速化を目標にしてなくて単純なコードと、フーリエ変換のコードの速度を比較する為にうｐしました

そうそう、GMPは数がでかすぎるとやけにすっぱりと諦めるよな。
たぶんメモリを喰いつくすとかの見通しをもとにしてはいるんだろうけど。

わかっている人も多いと思うが念のため書いておこう。
VC++でもGMPつかえる。
http://fp.gladman.plus.com/computing/gmp4win.htm
>>133
その程度の桁数では大丈夫なはず
ただ、桁数が32bitで表現できる数を超えるとだめらしい。

#include<iostream>
#include<bitset>
#include<vector>
#include<time.h>
using namespace std;
#define p 65537
#define ee 3
typedef vector< unsigned int > bsu; bsu X, e; vector < bsu > x;

unsigned int mod(unsigned int z){ int y=(z&65535)-(z>>16); return y<0?p+y:y;}

void Fourier_init(){ int k=1; e.resize(p-1); for(int i=0; i<p-1; i++){ e[i]=mod(k); k*=ee; }
X.resize(1<<16); x.resize(16); for(int i=0; i<16; i++) x[i].resize( (1<<(15-i)) );}

bsu& Fourier(bsu& y){
int i, d=1, b, N=32768; bitset<32> a;
for(i=0; i<N; i++) x[0][i]=mod( y[i] );
#define A a.to_ulong()
for(;;){ if(d<15){
b=16-d-1; N=1<<b;
for(i=0; i<N; i++) x[d][i]=mod( x[d-1][i] + e[(A<<b)&65535] * x[d-1][i+N] );
d++;}

else { X[A]=mod( x[14][0] + e[A] * x[14][1] );
X[A+32768]=mod( x[14][0] + e[A+32768] * x[14][1] );
for(d=14; d>=0; d--)if(a[d])a[d]=0; else{a[d]=1;break;}
if(d<1){if(d==0)d=1; else return X;}}}}

main(){ Fourier_init();
bsu a; a.resize(32768);
for(int i=0;i<32768;i++)a[i]=rand();
int cl=clock(); Fourier(a);Fourier(a); cl=clock()-cl;cout<<(cl+0.0)/1000;}

フーリエ変換のコード書きましたよ　あっているかは不明ですが
乗法で必要な分のフーリエ変換してみると0.15秒くらいかかりました　
最大桁数32768桁の変換にかかる時間
GMPで3万桁の掛け算実行するとエラーになって分からないんですけど
だれか3万桁が出来るコードをあげてもらえませんか?

フーリエ変換は
f(X)= a0 +　・・・　+ a(N-1) ・ X^(N-1)　、N=2^16 、p=N+1
に対して、f(1)、f(3)、・・・、f( 3^(N-1) ) mod p　を求める事と同じです
3^nはあらかじめ計算しておけば、単純に計算するとN*N程度の乗法が必要ですが
{ 3^n }の周期性から、f(1)、f(3)、・・・、f( 3^(N-1) )を一斉に求めようとすると
N * log2_N程度の乗法で済みます

modとか作らなくても%でよくね？
最近のコンパイラやCPUならそっちのほうが早い気がする。
まぁmodも最適化されて%と同じになるかも知れんが。

GMPとかFFTの話に入れなくてこんなことしか言えない・・・

いろいろ駄目すぎて、お話にならない

>>139
ほんとですね　%でも同じ速度になりました　しかし最適化をうまくやらないコンパイラ用にあった方がいいと思ってます
minGWの最新版で最大の最適化でやりました

N*Nの計算量のルーチンでは9.8秒かかるところが、フーリエ変換では0.15秒ですから600倍以上の高速化になりました
しかし、答えが違いました　N*Nの方はGMPと答えが一致しています

重大なバグがある気配　なんかFourierを複数個呼び出しても速度が鈍くなってない・・・ちゃんと計算していないって事ですよね・・・

最新のソースはどれ？

まだ全然出来ていないので最新版とか無いんですよ　一応動くやつが出来たらうｐします　その他上の指摘の変更をしたやつを

どう見ても同一犯による犯行です。触るだけ無駄。

source forge とか使え。本気で作るんなら、どうせsubversionとか使うようになるだろ。
ここに張られただけのソースを見て、自分で整形までしようとは思わん

　　L　I　S　P　で　や　れ

起きっぱで眠いテンションに任せて色々書いてみた。

>>33
データ構造は結局どうなってるの？構造決まらないとFFTで積でないよ？
多倍長浮動小数だとしたら指数部を記録しないと。

>>45
そんなことは既に考えられている。というより実質ω^iの全配列を作るようなのは無駄。

>>75
加減乗除のうち加減乗は各要素が整数になることは間違いない。
なので浮動小数を使った一般的なFFTを使っても整数化誤差が例えば0.1未満など保証されていれば問題ない。
除算については除算をどう行うか（筆算のようなやり方か、ニュートン法を使うか）、どの程度の誤差を
認めるか（整数演算でキッチリ末尾までか、浮動小数演算として末尾が１違うくらいは認めるか）次第で変わる。
（>>79を書いてるということは浮動小数型？）

>>90,92-94
自分が引用してきたサイトでも見直して配列と多倍長数とがどう対応してるか再考したほうがいい。
ついでにいうと>>98や>>122の発言を見てるとB進数という言葉を勘違いしてる気がする。

>>120-122
N(要素数)が任意でp固定なんてあり得ない。ωは1の原始p乗根のはずなのでFFTに使う要素数は
p-1で固定しなきゃいけない。ということはNがどんなに小さくてもサイズp-1のFFTと同じ時間がかかる。

>>1
ついでなので。リバーシのときと同じく他人の言うことをちゃんと聴くこと。聞くではなく聴くこと。
「まずは動くモノを作る」というのなら再帰使おうが％使おうが、FFTWやGMPより遅かろうが問題ないんだから。
まずはデバグのしやすさ、見易さ優先。アルゴリズムと違うレベルでの速度を求めるのはそれから。

pが　2^16+1ならば、最大桁は2^15以下にしないと桁あふれのため計算が間違えます
あとpは固定でも桁数が2^n ごとに小さいサイズの計算で済みますよ　
例えば、桁が8桁しか無ければ8*8行列のフーリエ変換を何度かやれば済みます
pは大きくしたり、小さくしたりすると不便なんです　小さくしたら2^16進数が扱えなくなり大きくすると、32bit内に計算がおさまらなくなります
工夫すれば、別のpでも良いんですが、簡単の為に固定してます

>>149
効率は良くないですけど、浮動小数点は使わずに格納可能な限り全ての桁を記録しようと思います
たとえば、1/3なら0.3333333・・・なら65536桁まで記録して小数点の位置も記録します
整数と同じ扱いにしようと思います
データ構造は、2^16進数で、0個に一番小さい数を入れます　
もし、非常に大きい数を扱うならば、同じライブラリを2度使って65536^32786進数でやろうとおもいます

シンプルに書いたけど値が一致しない　どこが駄目なんだろう
#include <iostream>
#include <vector>
using namespace std;
typedef vector< unsigned int > Lint;
#define p 65537
#define N 65536
#define M 1000
unsigned int mod_p(unsigned int z){ int y=(z&65535)-(z>>16); return y<0?p+y:y;}
Lint e; void Fourier_init(){ e.resize(p); int k=1; for(int i=0; i<p; i++){ e[i]=k%p; k*=3;}}

Lint Fourier(Lint f){ Lint F; F.resize(N,0);
for(int n=0; n<N; n++) for(int i=0; i<f.size(); i++)
F[n] = mod_p(F[n] + mod_p( e[(n*i)%N] * f[i] )); return F;}

main(){
Fourier_init(); Lint f, g; f.resize(M,0); g.resize(M,0);
for(int i=0; i<M; i++){ f[i]=rand(); g[i]=rand(); }
// FFT乗法
Lint F=Fourier(f), G=Fourier(g); Lint a; a.resize(N);
for(int i=0; i<N; i++)a[i]=mod_p(21846 * mod_p(F[N-i] * G[N-i]));
for(int i=0; i<N; i++)if(a[i]==N){a[i]=0; a[i+1]++;}
//普通の乗法
Lint b; b.resize(N,0); unsigned int x, k=0;
for(int n=0; n<=2*M-2; n++){ x=k; k=0;
for(int i=max(0,n-M+1); i<=min(M-1,n); i++)
{x+=f[i]*g[n-i]; k+=x>>16; x&=65535;}
b[n]=x;}
for(int i=0;i<100;i++)cout<<a[i]<<" "<<b[i]<<"\n";}

よりシンプルにしてみたけど駄目だ
#include <iostream>
#include <vector>
using namespace std;
typedef vector< unsigned int > Lint;
#define p 65537
#define N 65536
#define M 2
Lint e; void Fourier_init(){ e.resize(p); int k=1; for(int i=0; i<p; i++){ e[i]=k%p; k*=3;}}

Lint Fourier(Lint f){ Lint F; F.resize(N,0);
for(int n=0; n<N; n++) for(int i=0; i<f.size(); i++)
F[n]=(F[n] + ( e[(n*i)%N] * f[i] )%p)%p; return F;}

main(){
Fourier_init(); Lint f, g; f.resize(M,0); g.resize(M,0);f[0]=f[1]=g[0]=g[1]=2;

// FFT乗法
Lint F=Fourier(f), G=Fourier(g); Lint a; a.resize(N);
for(int i=0; i<N; i++)a[i]=(21846 * (F[N-i] * G[N-i])%p)%p;
for(int i=0; i<N; i++)if(a[i]==N){a[i]=0; a[i+1]++;}
//普通の乗法
Lint b; b.resize(N,0); unsigned int x, k=0;
for(int n=0; n<=2*M-2; n++){ x=k; k=0;
for(int i=max(0,n-M+1); i<=min(M-1,n); i++)
{x+=f[i]*g[n-i]; k+=x>>16; x&=65535;}
b[n]=x;}
for(int i=0;i<10;i++)cout<<a[i]<<" "<<b[i]<<"\n";}

がんばれー

とりあえず
（１）F=Fourier(f),
（２）G=Frourier(g)
としたとき
（３）A=FxG (要素毎の積)
を求めた後に
（４）a=(1/N)Fourier^(-1)(A)
とやる必要があるな。>>153では（３）と（４）が抜けてる。

> pは大きくしたり、小さくしたりすると不便なんです　小さくしたら2^16進数が扱えなくなり大きくすると、32bit内に計算がおさまらなくなります
> 工夫すれば、別のpでも良いんですが、簡単の為に固定してます

なら精度考慮してdoubleとかを使うFFTのほうが良くない？
あと、2^16進数使うなら大体(2^16)^2*要素数を表現できるだけの精度が必要になるはずなんだけど。

なんかあちこち変じゃない？

1から5000までの数字を全部掛けたときに
末尾に並ぶ0の個数を求めるのをおながいします

1000の倍数が5個、100の倍数が50個、10の倍数が500個、555個
5000の倍数が1個、500の倍数が10個、50の倍数が100個、5の倍数が1000個
全部足して1666個と適当に予想して見る。

ってか>>1なのかな？

>>157
全く同じ趣旨の問題がVBAスレかどっかでみた。
それは兎も角、鼬害。

>>158
末尾に0が並ぶということは、並んだ数をnとすると10のn乗の倍数ということ。
つまり、2のn乗*5のn乗の倍数になる。
このnを求めればいいわけだが、2の倍数は過剰に存在しているから
5の倍数がどれだけあって、それぞれ5の何乗の倍数か判ればいいことになる。
即ち、5の倍数の個数+25の倍数の個数+125の倍数の個数+……
だと思ったんだけどね。

#include<stdio.h>
int main(void){
int i,n=1,o=0;
for(i=2;i<5001;i++){
n*=i;
while(n%10==0)
n/=10,o++;
n%=10;
}
printf("%d\n",o);
return 0;
}

1240になったけど、これであってる？

>>160
それ途中でオーバーフローしないか？

下一桁だけ計算してるから
大丈夫なんじゃない？

>>159 の言う通りに計算してみると
5000 / 5 = 1000
5000 / 25 = 200
5000 / 125 = 40
5000 / 625 = 8
5000 / 3125 = 1

1249個になるね

計算してみた。
--
4228577926605543522201064200233584405390786674626646748849782402181358052708108200690899047871706387
5370847466573006854458784860666838127363372108937727876312793903630584621606439044789869822398719297
.
.全部で16300桁以上
.
3713673959956789408847059768595145050172415177460173514309909726155093783347200000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
00000000000000000000000000
--
ということで、1249だった。

>>160のコードをn%=100にしたら1249になった。

1じゃないけと。
剰余類によるFFTとdoubleによるFFTで乗算作ってみた。

C++で作ったんだけど、どうがんばっても、doubleのほうが早い。

とにかく法mがネック。
任意のmを用いることができないので設計が大変。

2^k+1の形にしておかないと、フーリエできないし、
素数じゃないといけないし、剰余が早くないといけない。
俺は
2^32(2^32-1)+1 2^32*(2^32-4)+1 2^16(2^16-16)+1 2^16(2^16-33)+1
などを使った。

で、剰余は結構早くできるようになったんだが、それでも遅い。
剰余が入るとdoubleのFFTに比べて大体2倍になる。

つぎに最大値の問題。
乗算後の各項の係数の最大値は「項数*(基数-1)^2」になるので、
法mはそれ以上であるか、もしくは法を二つ用意して2回FFTして、
中国剰余定理から求めないといけない。これでさらに2倍かかる。

結局、doubleで精度を考慮して作ったのと剰余類でがんばって作ったのだと、
4倍くらい速度差がでる。

>>73の言うとおりだった。
俺の実装が悪いだけかも知れんので、できればFAのソースが欲しい。

何か、１０年位前のMS-DOSの時代に逆戻りしたようなスレですね。(^^)
今のCPUでは不動点小数演算の方が整数演算よりはるかに高速です。
したがって、FFTはdoubleを使うのがいいと思います。
数千桁程度まではアセンブラを使った「多倍長BASIC」が高速です。
それ以上になると・・・

Ｃ＋＋による多倍長計算
http://www.vector.co.jp/authors/VA018507/index.html
の作者です。
Wondows上でも動くgcc(minGW など)インストールすれば使えます。
Dev-C++と日本語環境をインストールすればさらに使いやすくなります。
SN Library は日本人（私^^）の書いた英語ですから分かると思います。
ちなみに、50000！は関数にしていますが、Pen4 3.0GHzで0.75secです。
普通のfor文にすると比較していませんが、GMPより遅いでしょう。

>>167
大嘘と宣伝乙。

GMPって速いんだ。

今のCPUって整数より浮動小数点数のほうが速いんだ。

不動小数点（笑）は速いのか。メモメモ

そりゃぁ、不動なら速いだろw

浮動小数点でした。（汗）
PC98にコプロセッサを積むとフラクタル図形が１時間かかっていてものが
１分以内になったことに感動しました。１０年前の話です。
「Javaによるアルゴリズム事典」も読んだ人がいないようですね。

浮動小数点の方が速いっていうのはデマだろ　そしたら円周率の計算も浮動小数点計算でやってるのか

整数演算が得意なIntel Core 2 Extreme、Intel Core 2 Duoが上位を占めてるよ

Super PI　ベンチマークランキング
http://www16.big.or.jp/~bunnywk/superpi_index.cgi

一部大きな差が出たのは、HDBENCHのメモリ項目でメモリ項目で
Athlon 64 X2が20%以上早く、浮動小数点演算も35%早いものの、
逆に円周率演算のベンチマークソフトSuperπではCore 2 Duoが20%早いという結果に。
AMDはこのデータをもとにして「ひたすら円周率の計算をしたい人ならCore 2 Duoを買ってください」と皮肉った。
http://www.watch.impress.co.jp/akiba/hotline/20060812/etc_amdevent.html

>>174
円周率の計算といえば、SuperPI, pi_css6, PiFast, QuickPiなどは軒並み、浮動小数点演算で計算してるな。
整数演算で計算してるのはGMP系(gmp-chudnovsky.cやmpfr_const_pi関数)ぐらいか？

浮動小数点ならば175の結果はどうなんだよ?

SuperPIのボトルネックはランダムライト。演算性能以前の問題。
http://www.geocities.co.jp/SiliconValley-Bay/9439/negoto/ngt006.htm

今やSuperPIは実質、ランダムライトベンチマークと思っていい。
Core 2 はキャッシュがでかくて104万桁計算時に全部キャッシュに収まるから結果がいいだけ。

>>175をみてみろ
浮動小数点が得意なAthlon64 X2 6400+ は、PentiumM やPentium 4 にも負けてるぞ

そもそも浮動小数点演算と、ランダムライトは事実なのかが不明なんだが…　作者がそういったのか

>>181
作者に聞かなくてもVTune使えばすぐに分かることだろ。使ったことないのか？

C++で作って、floatやdouble、longと差し替えられるようにしてみれ。

このラッパーつくろうぜ
GMPのバイナリ置いてある

http://www.cs.nyu.edu/exact/core/gmp/

GMPのコンパイル済みのlibを使用して
C++の演算DLL作ったら、配布して良いんですか?

FFTの理論はまだ理解できてないんだけど

FFTを使った多倍長演算とモジュラ法(中国剰余定理)を使った多倍長演算
どっちが速いの？

ちなみに変換にかかるコストを度外視すればどちらもO(n)時間で計算可能。

180桁までならCRT、181桁以上ならFFT

すごい昔に作ったな。
多倍長演算ライブラリ。
386に最適化してるんで今のCPUじゃCの方が速そうだけど。

superPIの作者やPI計算のギネスに載ってた人（T大K氏）との意見交換とかもして、
当時としてはトップクラスだったと思うよ。

速度優先だと、
A=A+B, A=B+C が違う関数になるし、
A=B^2 は掛け算とは別の関数になる。

フーリエ変換は巨大メモリにバラバラにアクセスするので、
速度を稼ぐなら適宜並び替えると良いよ。
メモリに収まりきらないような巨大なFFTを高速に行う技術を、
今用のメモリとキャッシュ用の最適化にそのまま使えると思う。

ライブラリとして公開するつもりはないけど、
アドバイスはいろいろ出来ると思う。

>>188
ソース公開しろ

>>189
ソースもバイナリも公開するつもりはありません。

>188と>190では口調が違う罠。

まぁ、>188で「公開するつもりはない」と明言しているわけで。
恐らく、>188氏はくれくれに対しては厳しくても自力で頑張る人にはヒントをくれるんじゃない?

2007/11/24の段階で、このスレがここまで続くとは誰も思っていなかっただろうな。

またそういうことを

多倍長乗算用のFFTは並べ替えをしなくていい（ビットリバースで出てくる）
というのはどう思いますか、普通のFFTの方がキャッシュのヒット率が良いと
思うのですが

なにやら大抵の計算は出来るけど亀のように遅いソフトがソフトウェア板の
「なんの役にも立たないソフトを作ったら晒すスレ」に上がってるな

5000! を計算させたら１０秒くらいかかったぞｗ

>>195
CPUのスペックにもよるだろうが、
5000!で10秒って遅いか？

数値計算ソフトだと、表示も含めて1秒以内だな。
実行したとたんに答えが出ている。

そいつを誘導して、ここで開発させよう

>>198
そんな事したら>>1は何をすれば良いんだ。

うちの関数電卓で59!(関数電卓の値の限界点)を計算させると1分くらい表示が消えたまま反応がなくなる

>>200
電池が切れるんじゃね？ｗ

カシオやシャープの関数電卓は9.99……*10^99までだと思うので、69!が計算できた記憶があるが。
前者は1秒弱の遅延があるが後者は全くなかったので、69エントリのテーブルを持っているのではないかと噂されていた。

>>167
もしも嘘じゃないなら宣伝料支払うことになるけど？

保守

保守

http://www.craig-wood.com/nick/armprime/index.html

ARM では、整数　FFT　で性能出してるよ。１もがんばれ

全然進展無いけど、
もう作るの諦めた系なんですか？

1はやればできる子.

現状 今どうなってるの

俺が昔適当に作ったプログラムでも、5000!計算するのに3秒弱だったな
(CPUはNorthwoodの3.06GHzだ)。
10秒かかるようでは遅かろう。

もう１０年位になるんだから一人位はつくれてもいいはず

いまさらだけど、アルゴリズムに著作権は無いよ。
単なる汎用多倍長演算処理のソースになら著作権はあるだろうがな。

最近使ってなかったMS-DOS用の電卓ソフトに5000!を突っ込んだらオーバーフローと
言われた。3210!なら一瞬で終わった。
ということで10秒はやっぱり遅すぎ。

アルゴリズムに著作権はなくとも特許がある場合があるんじゃないか

>>214
既に公開してる（事前提の）物に特許とな？

サブマリンは？

それ特許じゃなくて投法じゃね？

ggrks

進捗はありましたか。

著作権フリーでならいくらでも転がってるでしょ。

探すスレじゃなくて作るスレなので

著作権フリーとか曖昧な事を言われるとかえって使えない

だよね

>>222
あの時お前が間違えなかったら
こんなことにはならなかった

rejectされちゃったけどBoost.XIntを使うべし#著作権フリーよりBoostライセンスですよ

うるさいゴミ
はいはいゴミゴミ

>>216
マジレスすれば、日本は公開するからサブマリンは有り得ない。

サブマリンは非公開のまま時期を延長できることを利用してるからな

まだかな