FFT

 

例题1:力

可以把题目给的式子转化为卷积的形式,然后通过FFT可以求得(推公式过程待补)

FFT
 //#include<bits/stdc++.h>
 #include<cstdio>
 #include<cmath>
 #include<iostream>
 #include<algorithm>
 #include<vector>
 #include<map>
 #include<stack>
 #include<queue>
 #include<cstring>
 #include<string>
 #include<set>
 #include<complex>
 //#include<unordered_map>
 using namespace std;
 #define rep(i,a,b) for(int i=a;i<=b;i++)
 #define dep(i,b,a) for(int i=b;i>=a;i--)
 #define m_p make_pair
 #define fi first
 #define se second
 #define pb push_back
 #define sz(x) (int)(x).size()
 #define pi acos(-1)
 #define Io ios::sync_with_stdio(false);cin.tie(0)
 #define io ios::sync_with_stdio(false)
 #define IOS ios::sync_with_stdio(false),cin.tie(0),cout.tie(0)
 typedef  long long ll;
 typedef pair<int,int> pii;
 typedef pair<ll,ll>pll;
 typedef complex<double> Comp;
 const int inf=0x3f3f3f3f;
 const ll INF=0x3f3f3f3f3f3f3f3f;
 const double  eps=1e-9;
 inline int gcd(int a,int b){return b?gcd(b,a%b):a;}
 inline int lcm(int a,int b){return a*b/gcd(a,b);}
 namespace  IO  
 {
     const int len=4e7;char buf[len];int sz,p;
     void begin(){p=0;sz=fread(buf,1,len,stdin);}
     inline bool read(ll &x)
     {
         if (p==sz)return 0;int f=1,d=0;char s=buf[p++];
         while(s<'0'||s>'9'&&p<sz){if(s=='-') f=-1;s=buf[p++];}
         while(s>='0'&&s<='9'&&p<sz){d=d*10+s-'0';s=buf[p++];}
         x=f*d; return p!=sz;
     }
 }
 inline int rd()
 {
 int x=0,f=1;
 char ch=getchar();
 while(ch<'0'||ch>'9') { if(ch=='-') f=-1; ch=getchar();}
 while(ch>='0'&&ch<='9') { x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
 return x*f;
 }
 const int maxn=5e5+50;
 struct Complex
{
    double x, y;
    Complex(double x=0,double y=0):x(x),y(y){}

} a[maxn], b[maxn], c[maxn];
Complex operator + (Complex a, Complex b) { return Complex(a.x + b.x , a.y + b.y);}
Complex operator - (Complex a, Complex b) { return Complex(a.x - b.x , a.y - b.y);}
Complex operator * (Complex a, Complex b) { return Complex(a.x * b.x - a.y * b.y , a.x * b.y + a.y * b.x);} 
int  rev [maxn];
int lim = 1, l;
void fft(Complex f[],double tag)
{
    for (int i = 0; i < lim;++i) 
        if(i<rev[i])   swap(f[i], f[rev[i]]);
    for (int mid = 1; mid <lim;mid<<=1)
    {
        Complex wn(cos(pi / mid), tag * sin(pi / mid));
        for (int R = mid << 1, j = 0; j < lim;j+=R)
        {
            Complex w(1, 0);
            for (int k = 0; k < mid;k++,w=w*wn)
            {
                Complex x = f[j + k], y = w * f[j + mid + k];
                f[j + k] = x + y;
                f[j + mid + k] = x - y;
            }
        }
     }
}
 void run()
 {
     int n = rd();
     for (int i = 1;i<=n;++i)
     {
        scanf("%lf", &a[i].x);
        b[i].x = 1.0 / i / i;
        c[n - i + 1].x = a[i].x;
     }
     while(lim<=2*n)
         lim <<= 1, l++;
     for (int i = 0; i < lim;++i)
         rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (l - 1)));
     fft(a, 1), fft(b, 1), fft(c, 1);
     for (int i = 0; i <= lim;++i)
         a[i] = a[i] * b[i], c[i] = b[i] * c[i];
     fft(a, -1), fft(c, -1);
     for (int i = 1; i <= n;++i)
         a[i].x /= lim, c[i].x /= lim;
     for (int i = 1;i<=n;++i)
         printf("%.3lf\n", a[i].x - c[n-i+1].x);
 }
 int main()
 {
    run();
    return 0;
 }
 /* 
 start time:
 over time:
 */
 
View Code

 

例题2:大数乘法

FFT
 //#include<bits/stdc++.h>
 #include<cstdio>
 #include<cmath>
 #include<iostream>
 #include<algorithm>
 #include<vector>
 #include<map>
 #include<stack>
 #include<queue>
 #include<cstring>
 #include<string>
 #include<set>
 #include<complex>
 //#include<unordered_map>
 using namespace std;
 #define rep(i,a,b) for(int i=a;i<=b;i++)
 #define dep(i,b,a) for(int i=b;i>=a;i--)
 #define m_p make_pair
 #define fi first
 #define se second
 #define pb push_back
 #define sz(x) (int)(x).size()
 #define pi acos(-1)
 #define Io ios::sync_with_stdio(false);cin.tie(0)
 #define io ios::sync_with_stdio(false)
 #define IOS ios::sync_with_stdio(false),cin.tie(0),cout.tie(0)
 typedef  long long ll;
 typedef pair<int,int> pii;
 typedef pair<ll,ll>pll;
 typedef complex<double> Comp;
 const int inf=0x3f3f3f3f;
 const ll INF=0x3f3f3f3f3f3f3f3f;
 const double  eps=1e-9;
 inline int gcd(int a,int b){return b?gcd(b,a%b):a;}
 inline int lcm(int a,int b){return a*b/gcd(a,b);}
 namespace  IO  
 {
     const int len=4e7;char buf[len];int sz,p;
     void begin(){p=0;sz=fread(buf,1,len,stdin);}
     inline bool read(ll &x)
     {
         if (p==sz)return 0;int f=1,d=0;char s=buf[p++];
         while(s<'0'||s>'9'&&p<sz){if(s=='-') f=-1;s=buf[p++];}
         while(s>='0'&&s<='9'&&p<sz){d=d*10+s-'0';s=buf[p++];}
         x=f*d; return p!=sz;
     }
 }
 inline int rd()
 {
 int x=0,f=1;
 char ch=getchar();
 while(ch<'0'||ch>'9') { if(ch=='-') f=-1; ch=getchar();}
 while(ch>='0'&&ch<='9') { x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
 return x*f;
 }
 const int maxn=5e6+50;
 struct Complex
{
    double x, y;
    Complex(double x=0,double y=0):x(x),y(y){}

} a[maxn], b[maxn];
Complex operator + (Complex a, Complex b) { return Complex(a.x + b.x , a.y + b.y);}
Complex operator - (Complex a, Complex b) { return Complex(a.x - b.x , a.y - b.y);}
Complex operator * (Complex a, Complex b) { return Complex(a.x * b.x - a.y * b.y , a.x * b.y + a.y * b.x);} 
int  rev [maxn];
int lim = 1, l;
void fft(Complex f[],double tag)
{
    for (int i = 0; i < lim;++i) 
        if(i<rev[i])   swap(f[i], f[rev[i]]);
    for (int mid = 1; mid <lim;mid<<=1)
    {
        Complex wn(cos(pi / mid), tag * sin(pi / mid));
        for (int R = mid << 1, j = 0; j < lim;j+=R)
        {
            Complex w(1, 0);
            for (int k = 0; k < mid;k++,w=w*wn)
            {
                Complex x = f[j + k], y = w * f[j + mid + k];
                f[j + k] = x + y;
                f[j + mid + k] = x - y;
            }
        }
     }
}
char s1[maxn],s2[maxn];
int ans[maxn];
void run()
{
    scanf("%s", s1 );
    scanf("%s", s2 );
    int len1 = strlen(s1), len2 = strlen(s2);
    int tot1 = 0, tot2 = 0;
    for (int i = len1-1; i >= 0;--i)
        a[tot1++].x = s1[i] - '0';
    for (int i = len2-1; i >= 0;--i)
        b[tot2++].x = s2[i] - '0';
    while (lim <= len1 + len2)
            lim <<= 1, l++;
    for (int i = 0; i < lim; ++i)
        rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (l - 1)));
    fft(a, 1), fft(b, 1);
    for (int i = 0; i <= lim; ++i)
        a[i] = a[i] * b[i];
    fft(a, -1);
    for (int i = 0; i <= lim;++i)
    {
        ans[i] += int(a[i].x / lim + 0.5);
        if(ans[i]>=10)
        {
            ans[i + 1] = ans[i] / 10;
            ans[i] %= 10;
            lim += (i == lim);
        }
    }
    while(!ans[lim]&&lim>=1)
        lim--;
    lim++;
    while (--lim>=0) printf("%d",ans[lim]);
    puts("");
}
int main()
{
    run();
    return 0;
}
 /* 
 start time:
 over time:
 */
 
View Code

 

上一篇:CSUST 看直播 题解(二分+dp)


下一篇:2021-07-28