zoukankan      html  css  js  c++  java
  • hihocoder #1388 : Periodic Signal NTT&FFT

    传送门: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
    ",a)
    #define ptlld(a) printf("%I64d
    ",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

    不为失败找借口,只为成功找方法
  • 相关阅读:
    zw版【转发·台湾nvp系列Delphi例程】HALCON MoveRectangle
    zw版【转发·台湾nvp系列Delphi例程】HALCON SetIcon1
    zw版【转发·台湾nvp系列Delphi例程】HALCON SetIcon2
    zw版【转发·台湾nvp系列Delphi例程】HALCON DispArc
    zw版【转发·台湾nvp系列Delphi例程】HALCON HSerializedItem
    zw版【转发·台湾nvp系列Delphi例程】HALCON union1
    zw版【转发·台湾nvp系列Delphi例程】HALCON SqrtImage
    zw版【转发·台湾nvp系列Delphi例程】HALCON SubImage
    zw版【转发·台湾nvp系列Delphi例程】HALCON SetWindowExtent
    Raphael绘制圆圈环绕方法
  • 原文地址:https://www.cnblogs.com/youmi/p/5905679.html
Copyright © 2011-2022 走看看