精度损失到自闭，菜鸡求助

tommy0221 @ 2020-03-29 20:55:05

#include<bits/stdc++.h>
using namespace std;
#define rint register int
typedef long double ld;
typedef long long LL;
inline int rd(){
   int x=0,f=1;
   char ch=getchar();
   while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}
   while(isdigit(ch)) x=x*10+(ch^48),ch=getchar();
   return x*f;
}
const int N=400010;
const int M=32768;
const ld pi=acos(-1);
struct cp{
    ld x,y;
    cp(){}
    cp(ld _x,ld _y){x=_x,y=_y;}
};
inline cp operator + (const cp &a,const cp &b){return cp(a.x+b.x,a.y+b.y);}
inline cp operator - (const cp &a,const cp &b){return cp(a.x-b.x,a.y-b.y);}
inline cp operator * (const cp &a,const cp &b){return cp(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
int n,m,p,lim,lg,rev[N],a[N],b[N],c[N];
cp A1[N],B1[N],A2[N],B2[N],A[N],B[N],C[N];
inline void FFT(cp *a,int lim,int o) {
    for(rint i=0;i<lim;++i)if(i>rev[i])swap(a[i],a[rev[i]]);
    for(rint i=1;i<lim;i<<=1) {
        cp wn=cp(cos(pi/i),o*sin(pi/i));
        for(rint j=0;j<lim;j+=(i<<1)) {
            cp w0=cp(1,0);
            for(rint k=0;k<i;++k,w0=w0*wn) {
                cp X=a[j+k],Y=w0*a[i+j+k];
                a[j+k]=X+Y;
                a[i+j+k]=X-Y;
            }
        }
    }
}
inline void init(int n) {
    for(lim=1,lg=0;lim<=n;++lg,lim<<=1);
    for(rint i=0;i<lim;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg-1));
}
inline void MTT(int mod,int *a,int *b,int *ans) {
    for(rint i=0;i<lim;++i) {
        A1[i].x=a[i]/M,B1[i].x=a[i]%M;
        A2[i].x=b[i]/M,B2[i].x=b[i]%M;
    }
    FFT(A1,lim,1),FFT(B1,lim,1),FFT(A2,lim,1),FFT(B2,lim,1);
    for(rint i=0;i<lim;++i) {
        A[i]=A1[i]*A2[i];
        B[i]=A1[i]*B2[i]+A2[i]*B1[i];
        C[i]=B1[i]*B2[i];
    }
    FFT(A,lim,-1),FFT(B,lim,-1),FFT(C,lim,-1);
    for(rint i=0;i<lim;++i) {
        ans[i]=(LL)round(A[i].x)%mod;
        ans[i]=(1ll*ans[i]*M%mod+(LL)round(B[i].x)%mod)%mod;
        ans[i]=(1ll*ans[i]*M%mod+(LL)round(C[i].x)%mod)%mod;
    }
}
signed main() {
    n=rd(),m=rd(),p=rd();
    for(rint i=0;i<=n;++i)a[i]=rd()%p;
    for(rint i=0;i<=m;++i)b[i]=rd()%p;
    init(n+m),MTT(p,a,b,c);
    for(rint i=0;i<=n+m;++i)printf("%d%c",c[i]," \n"[i==n+m]);
    return 0;
}

by Smile_Cindy @ 2020-03-29 20:59:17

@世外明月 我的代码给你看看，虽然我也不知道为什么

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;
const long double PI=acos(-1);
const int MAX_N=1<<18;
const int base=1<<15;
int n,m,Mod; 
struct Complex
{
    long double x,y;
    Complex(){}
    Complex(long double _x,long double _y)
    {
        x=_x;
        y=_y;
    }
};
Complex operator + (Complex w1,Complex w2)
{
    return Complex(w1.x+w2.x,w1.y+w2.y);
}
Complex operator - (Complex w1,Complex w2)
{
    return Complex(w1.x-w2.x,w1.y-w2.y);
}
Complex operator * (Complex w1,Complex w2)
{
    return Complex(w1.x*w2.x-w1.y*w2.y,w1.x*w2.y+w1.y*w2.x);
}
int rev[MAX_N];
void FFT(Complex A[],int n,int flag)
{
    for(int i=0;i<n;i++)if(rev[i]>i)swap(A[i],A[rev[i]]);
    for(int i=1;i<n;i<<=1)
    {
        Complex w1=Complex(cos(2*PI*flag/(i<<1)),sin(2*PI*flag/(i<<1)));
        for(int j=0;j<n;j+=(i<<1))
        {
            Complex w2=Complex(1,0);
            for(int k=0;k<i;k++)
            {
                Complex t1=A[j+k],t2=A[i+j+k]*w2;
                A[j+k]=t1+t2;
                A[i+j+k]=t1-t2;
                w2=w2*w1;
            }
        }
    }
    if(flag==-1)for(int i=0;i<n;i++)A[i].x/=n;
}
int len=1,L=0;
void init()
{
    while(len<=n+m)
    {
        len<<=1;
        L++;
    }
    for(int i=1;i<len;i++)rev[i]=rev[i>>1]>>1|(i&1)<<(L-1);
}
Complex A1[MAX_N],A2[MAX_N],B1[MAX_N],B2[MAX_N];
Complex A[MAX_N],B[MAX_N],C[MAX_N];
int F[MAX_N],G[MAX_N];
int ans[MAX_N];
signed main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    cin>>n>>m>>Mod;
    for(int i=0;i<=n;i++)cin>>F[i];
    for(int i=0;i<=m;i++)cin>>G[i];
    init();
    for(int i=0;i<=n;i++)A1[i].x=F[i]/base,B1[i].x=F[i]%base;
    for(int i=0;i<=m;i++)A2[i].x=G[i]/base,B2[i].x=G[i]%base;
    FFT(A1,len,1);
    FFT(A2,len,1);
    FFT(B1,len,1);
    FFT(B2,len,1);
    for(int i=0;i<len;i++)A[i]=A1[i]*A2[i];
    for(int i=0;i<len;i++)B[i]=A1[i]*B2[i]+A2[i]*B1[i];
    for(int i=0;i<len;i++)C[i]=B1[i]*B2[i];
    FFT(A,len,-1);
    FFT(B,len,-1);
    FFT(C,len,-1);
    for(int i=0;i<len;i++)
    {
        ans[i]=(LL)(A[i].x+0.5)%Mod;
        ans[i]=(1LL*ans[i]*base+(LL)(B[i].x+0.5))%Mod;
        ans[i]=(1LL*ans[i]*base+(LL)(C[i].x+0.5))%Mod;
    }
    for(int i=0;i<=n+m;i++)cout<<ans[i]<<" \n"[i==n+m];
    return 0;
}

by Smile_Cindy @ 2020-03-29 20:59:31

难得碰到一个和我一样7次DFT的。

by tommy0221 @ 2020-03-29 20:59:39

@Alpha 谢谢，我看看

by tommy0221 @ 2020-03-29 21:00:16

网上随便找了一篇qwq

by tommy0221 @ 2020-03-29 21:19:38

啊，IDFT的时候忘记/lim了，但是还是炸

by tommy0221 @ 2020-03-29 21:22:17

@Alpha AC了，谢谢！！！其实问题是我过了样例以后一直在测大数据没测样例，被我瞎改结果改错了，IDFT的时候忘记/lim了，我个sb