FFT计算多项式乘法与大数乘法
多项式乘法
![在这里插入图片描述](https://i-blog.csdnimg.cn/blog_migrate/ce46a26cc00ef41d8238e56d947445f5.png)
#include <iostream>
#include<stdio.h>
#include <complex>
const int maxN=1e6;
const double pi=3.1415926;
using namespace std;
void FFT(complex<double> *a,int n,int inv);
complex<double> a[maxN],b[maxN];
int n,n1,n2,cnt,rev[maxN],c[maxN],sign=1;
int main()
{
scanf("%d%d",&n1,&n2);
n=1;
while(n<n1+n2)n<<=1;
for(int i=0;i<n1;i++){
scanf("%lf",&a[i].real());
}
for(int i=0;i<n2;i++){
scanf("%lf",&b[i].real());
}
FFT(a,n,1);
FFT(b,n,1);
for(int i=0;i<=n;i++){
a[i]=a[i]*b[i];
}
FFT(a,n,-1);
for(int i=0;i<m;i++){
printf("%d ",(int)(a[i].real()/n+0.5));
}
return 0;
}
void FFT(complex<double> *a,int n,int inv){
int bit=0;
while((1<<bit)<n)bit++;
for(int i=0;i<n;i++){
rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
if(i<rev[i])
swap(a[i],a[rev[i]]);
}
for(int mid=1;mid<n;mid<<=1){
complex<double> wn(cos(pi/mid),inv*sin(pi/mid));
for(int i=0;i<n;i+=mid<<1){
complex<double> w(1,0);
for(int j=0;j<mid;j++){
complex<double> x=a[i+j],y=w*a[i+j+mid];
a[i+j]=x+y;a[i+j+mid]=x-y;
w*=wn;
}
}
}
}
大数乘法
数据
函数read将字符串s1,s2也就是大数转化为复数的实数部分;经过FFT转化化为c的系数,最后每位相加,注意进位;
#include <iostream>
#include <stdio.h>
#include <string.h>
#include <complex>
using namespace std;
const int maxN=1001000;
const double pi=acos(-1.0);
void FFT(complex<double> *a,int n,int inv);
int read(char s[],complex<double> *a);
int n,n1,n2,cnt,rev[maxN],c[maxN],sign=1;
void trans(int c[],int n);
void trans_1(int c[],int n);
complex<double> a[maxN],b[maxN];
int main()
{
char s1[maxN],s2[maxN];
while(~scanf("%s%s",s1,s2)){
sign=1;
n=1;
n1=read(s1,a);
n2=read(s2,b);
if(sign==-1)
printf("-");
while(n<n1+n2)n<<=1;
FFT(a,n,1);
FFT(b,n,1);
for(int i=0;i<n;i++){
a[i]=a[i]*b[i];
}
FFT(a,n,-1);
trans(c,n);
memset(a,0,(n+10)*sizeof(complex<double>));
memset(b,0,(n+10)*sizeof(complex<double>));
memset(c,0,(n+10)*sizeof(int));
}
return 0;
}
void trans(int c[],int n){
for(int i=0;i<n;i++){
c[i]+=(int)(a[i].real()/n+0.5);
if(c[i]>=10)
c[i+1]+=c[i]/10,c[i]%=10,n+=(i==n);
}
while(!c[n]&&n>0){
n--;
}
for(int i=n;i>=0;i--){
printf("%d",c[i]);
}
printf("\n");
}
void trans_1(int c[],int n){
for(int i=0;i<n;i++){
c[i]=a[i].real()/n+0.5;
}
int bit_add=0,k=0;
for(int i=0;i<n;i++){
int temp=c[i]+bit_add;
if(temp>=10){
c[k++]=temp%10;
bit_add=temp/10;
}
else{
c[k++]=temp;
bit_add=0;
}
}
while(bit_add){
c[k++]=bit_add;
bit_add/=10;
}
k--;
while(c[k]==0&&k>0)
k--;
for(int i=k;i>=0;i--){
printf("%d",c[i]);
}
printf("\n");
}
void FFT(complex<double> *a,int n,int inv){
int bit=0;
while((1<<bit)<n)bit++;
for(int i=0;i<n;i++){
rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
if(i<rev[i])
swap(a[i],a[rev[i]]);
}
for(int mid=1;mid<n;mid<<=1){
complex<double> wn(cos(pi/mid),inv*sin(pi/mid));
for(int i=0;i<n;i+=mid<<1){
complex<double> w(1,0);
for(int j=0;j<mid;j++){
complex<double> x=a[i+j],y=w*a[i+j+mid];
a[i+j]=x+y;a[i+j+mid]=x-y;
w*=wn;
}
}
}
}
int read(char s[],complex<double> *a){
int n=0;
for(int i=strlen(s)-1;i>=0;i--){
if(s[i]=='-'){
if(sign==1)
sign=-1;
else
sign=1;
}
else if(s[i]>='0'&&s[i]<='9'){
complex<double> m(s[i]-'0',0);
a[n++]=m;
}
}
return n;
}