hihocoder #1388 : Periodic Signal NTT&FFT

Posted 2020-08-10 youmi

传送门：hihocoder #1388 : Periodic Signal 

先来几个大牛传送门：  (模板) NTT long long 版

解法一：因为我们知道FFT会精度不够，所以坚持用NTT，但是模数不够大，然后就一直GG，看来我们的搜索姿势也有问题，居然没有搜到上面大神的板子，真的是GG

http://www.cnblogs.com/WABoss/p/5903927.html

/**************************************************************
    Problem:
    User: youmi
    Language: C++
    Result: Accepted
    Time:
    Memory:
****************************************************************/
//#pragma comment(linker, "/STACK:1024000000,1024000000")
//#include<bits/stdc++.h>
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <map>
#include <stack>
#include <set>
#include <sstream>
#include <cmath>
#include <queue>
#include <deque>
#include <string>
#include <vector>
#define zeros(a) memset(a,0,sizeof(a))
#define ones(a) memset(a,-1,sizeof(a))
#define sc(a) scanf("%d",&a)
#define sc2(a,b) scanf("%d%d",&a,&b)
#define sc3(a,b,c) scanf("%d%d%d",&a,&b,&c)
#define scs(a) scanf("%s",a)
#define sclld(a) scanf("%I64d",&a)
#define pt(a) printf("%d\\n",a)
#define ptlld(a) printf("%I64d\\n",a)
#define rep(i,from,to) for(int i=from;i<=to;i++)
#define irep(i,to,from) for(int i=to;i>=from;i--)
#define Max(a,b) ((a)>(b)?(a):(b))
#define Min(a,b) ((a)<(b)?(a):(b))
#define lson (step<<1)
#define rson (lson+1)
#define eps 1e-6
#define oo 0x3fffffff
#define TEST cout<<"*************************"<<endl
const double pi=4*atan(1.0);

using namespace std;
typedef long long ll;

int n;
const ll P = 50000000001507329LL; //190734863287 * 2 ^ 18 + 1
//const ll P = 1004535809LL; //479 * 2 ^ 21 + 1
//const ll P = 1004535809; // 119 * 2 ^ 23 + 1
const int N = 1 << 18;
const int G = 3;
int len;
ll A[N],B[N];
long long a[N],b[N],wn[30];

ll mul(ll x, ll y) {
    return (x * y - (ll)(x / (long double)P * y + 1e-3) * P + P) % P;
}

ll qpow(ll x, ll k, ll p) {
    ll ret = 1;
    while(k) {
        if(k & 1) ret = mul(ret, x);
        k >>= 1;
        x = mul(x, x);
    }
    return ret;
}

void getwn()
{
    for(int i = 1; i <= 18; ++i)
    {
        int t = 1 << i;
        wn[i] = qpow(G, (P - 1) / t, P);
    }
}

void change(ll *y, int len)
{
    for(int i = 1, j = len / 2; i < len - 1; ++i)
    {
        if(i < j) swap(y[i], y[j]);
        int k = len / 2;
        while(j >= k)
        {
            j -= k;
            k /= 2;
        }
        j += k;
    }
}

void NTT(ll *y, int len, int on)
{
    change(y, len);
    int id = 0;
    for(int h = 2; h <= len; h <<= 1)
    {
        ++id;
        for(int j = 0; j < len; j += h)
        {
            ll w = 1;
            for(int k = j; k < j + h / 2; ++k)
            {
                ll u = y[k];
                ll t = mul(y[k+h/2], w);
                y[k] = u + t;
                if(y[k] >= P) y[k] -= P;
                y[k+h/2] = u - t + P;
                if(y[k+h/2] >= P) y[k+h/2] -= P;
                w = mul(w, wn[id]);
            }
        }
    }
    if(on == -1)
    {
        for(int i = 1; i < len / 2; ++i) swap(y[i], y[len-i]);
        ll inv = qpow(len, P - 2, P);
        for(int i = 0; i < len; ++i)
            y[i] = mul(y[i], inv);
    }
}
void work()///卷积，点乘，插值
{
    NTT(a,len,1);
    NTT(b,len,1);
    for(int i=0;i<len;i++)
        a[i]=mul(a[i],b[i]);
    NTT(a,len,-1);
}
ll solve()
{
    zeros(a);
    zeros(b);
    rep(i,0,n-1)
        a[i]=A[i];
    rep(i,0,n-1)
        b[i]=B[i];
    reverse(b,b+n);
    work();
    ll ans=0;
    rep(i,0,n-1)
        a[i]+=a[i+n];
    rep(i,0,n-1)
        ans=max(ans,2*a[i]);
    return ans;
}
int main()
{
    #ifndef ONLINE_JUDGE
    freopen("in.txt","r",stdin);
    #endif
    int T_T;
    scanf("%d",&T_T);
    getwn();
    for(int kase=1;kase<=T_T;kase++)
    {
        sc(n);
        len=1;
        while(len<=2*n)
            len<<=1;
        rep(i,0,n-1)
            cin>>A[i];
        rep(i,0,n-1)
            cin>>B[i];
        ll temp=0;
        rep(i,0,n-1)
            temp+=A[i]*A[i];
        rep(i,0,n-1)
            temp+=B[i]*B[i];
        ll ans=solve();
        ans=temp-ans;
        cout<<(ans)<<endl;
    }
}

View Code

 

解法二：这个解法确实很漂亮，比赛的时候一直徘徊找一个更大的 模数，然后就GG了，http://www.cnblogs.com/smartweed/p/5903838.html

解法三：其实这种解法我们也尝试了，队友说NTT搞了那么久，说明暴力应该可以，不过最后只剩几分钟来不及找到合适的循环次数，http://www.cnblogs.com/cshg/p/5905398.html