著作権フリーのC++高速汎用多倍長演算を作るスレ
■ このスレッドは過去ログ倉庫に格納されています
01531
2007/12/04(火) 09:16:45#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";}
■ このスレッドは過去ログ倉庫に格納されています