hdu 5730 Shell Necklace [分治fft | 多项式求逆]

hdu 5730 Shell Necklace

题意：求递推式(f_n = sum_{i=1}^n a_i f_{n-i})，模313

---

多么优秀的模板题

可以用分治fft，也可以多项式求逆

分治fft

注意过程中把r-l+1当做次数界就可以了，因为其中一个向量是[l,mid]，我们只需要[mid+1,r]的结果。

多项式求逆

变成了

[A(x) = frac{f_0}{1-B(x)} ]

的形式

要用拆系数fft，直接把之前的代码复制上就可以啦

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;
const int N = (1<<18) + 5, mo = 313;
const double PI = acos(-1.0);
inline int read() {
    char c=getchar(); int x=0,f=1;
    while(c<'0'||c>'9') {if(c=='-')f=-1;c=getchar();}
    while(c>='0'&&c<='9') {x=x*10+c-'0';c=getchar();}
    return x*f;
}

struct meow{
    double x, y;
    meow(double a=0, double b=0):x(a), y(b){}
};
meow operator +(meow a, meow b) {return meow(a.x+b.x, a.y+b.y);}
meow operator -(meow a, meow b) {return meow(a.x-b.x, a.y-b.y);}
meow operator *(meow a, meow b) {return meow(a.x*b.x-a.y*b.y, a.x*b.y+a.y*b.x);}
meow conj(meow a) {return meow(a.x, -a.y);}
typedef meow cd;

namespace fft {
    int maxlen = 1<<17, rev[N];
    cd omega[N], omegaInv[N];
    void init() {
        for(int i=0; i<maxlen; i++) {
            omega[i] = cd(cos(2*PI/maxlen*i), sin(2*PI/maxlen*i));
            omegaInv[i] = conj(omega[i]);
        }
    }

    void dft(cd *a, int n, int flag) {
        cd *w = flag == 1 ? omega : omegaInv;
        int k = 0; while((1<<k) < n) k++;
        for(int i=0; i<n; i++) {
            rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
            if(i < rev[i]) swap(a[i], a[rev[i]]);
        }
        for(int l=2; l<=n; l<<=1) {
            int m = l>>1;
            for(cd *p = a; p != a+n; p += l) 
                for(int k=0; k<m; k++) {
                    cd t = w[maxlen/l*k] * p[k+m];
                    p[k+m] = p[k] - t;
                    p[k] = p[k] + t;
                }
        }
        if(flag == -1) for(int i=0; i<n; i++) a[i].x /= n;
    }
}

int n, q[N], f[N];

cd a[N], b[N];
void cdq(int l, int r) {
    using fft::dft;
    if(l == r) return;
    int mid = (l+r)>>1;
    cdq(l, mid);
    int n = 1; while(n < r-l+1) n <<= 1;
    for(int i=0; i<n; i++) a[i] = b[i] = cd();
    for(int i=l; i<=mid; i++) a[i-l].x = f[i];
    for(int i=0; i<=r-l; i++) b[i].x = q[i];
    dft(a, n, 1); dft(b, n, 1);
    for(int i=0; i<n; i++) a[i] = a[i] * b[i];
    dft(a, n, -1);
    for(int i=mid+1; i<=r; i++) f[i] = (f[i] + (int) floor(a[i-l].x + 0.5)) %mo;
    cdq(mid+1, r);
}

int main() {
    freopen("in", "r", stdin);
    fft::init();
    while( (n=read()) ) {
        //f[0] = 1;
        for(int i=1; i<=n; i++) f[i] = q[i] = read() % mo;
        cdq(0, n);
        printf("%d
", f[n] %mo);
    }
}

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;
const int N = (1<<18) + 5, mo = 313, P = 313;
const double PI = acos(-1.0);
inline int read() {
    char c=getchar(); int x=0,f=1;
    while(c<'0'||c>'9') {if(c=='-')f=-1;c=getchar();}
    while(c>='0'&&c<='9') {x=x*10+c-'0';c=getchar();}
    return x*f;
}

inline int Pow(int a, int b) {
	int ans = 1;
	for(; b; b>>=1, a=a*a%P)
		if(b&1) ans=ans*a%P;
	return ans;
}

struct meow{
    double x, y;
    meow(double a=0, double b=0):x(a), y(b){}
};
meow operator +(meow a, meow b) {return meow(a.x+b.x, a.y+b.y);}
meow operator -(meow a, meow b) {return meow(a.x-b.x, a.y-b.y);}
meow operator *(meow a, meow b) {return meow(a.x*b.x-a.y*b.y, a.x*b.y+a.y*b.x);}
meow conj(meow a) {return meow(a.x, -a.y);}
typedef meow cd;

namespace fft {
    int maxlen = 1<<18, rev[N];
    cd omega[N], omegaInv[N];
    void init() {
        for(int i=0; i<maxlen; i++) {
            omega[i] = cd(cos(2*PI/maxlen*i), sin(2*PI/maxlen*i));
            omegaInv[i] = conj(omega[i]);
        }
    }

    void dft(cd *a, int n, int flag) {
        cd *w = flag == 1 ? omega : omegaInv;
        int k = 0; while((1<<k) < n) k++;
        for(int i=0; i<n; i++) {
            rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
            if(i < rev[i]) swap(a[i], a[rev[i]]);
        }
        for(int l=2; l<=n; l<<=1) {
            int m = l>>1;
            for(cd *p = a; p != a+n; p += l) 
                for(int k=0; k<m; k++) {
                    cd t = w[maxlen/l*k] * p[k+m];
                    p[k+m] = p[k] - t;
                    p[k] = p[k] + t;
                }
        }
        if(flag == -1) for(int i=0; i<n; i++) a[i].x /= n;
    }

    cd a[N], b[N], c[N], d[N]; int z[N];
    void inverse(int *x, int *y, int l) {
        if(l == 1) {y[0] = 1; return;}
        inverse(x, y, l>>1);
		int n = l<<1;
        
        for(int i=0; i<l; i++) a[i] = cd(y[i]>>15), b[i] = cd(y[i]&32767);
        for(int i=l; i<n; i++) a[i] = b[i] = cd();
        dft(a, n, 1); dft(b, n, 1);
        for(int i=0; i<n; i++) {
            cd _a = a[i], _b = b[i];
            a[i] = _a * _a;
            b[i] = _a * _b + _a * _b;
            c[i] = _b * _b;
        }
        dft(a, n, -1); dft(b, n, -1); dft(c, n, -1);
        for(int i=0; i<l; i++)
            z[i] = ( (ll(a[i].x + 0.5) %mo << 30) %mo + (ll(b[i].x + 0.5) %mo << 15) %mo + ll(c[i].x + 0.5) %mo) %mo;

        for(int i=0; i<l; i++) 
            a[i] = cd(x[i]>>15), b[i] = cd(x[i]&32767), c[i] = cd(z[i]>>15), d[i] = cd(z[i]&32767); 
        for(int i=l; i<n; i++) a[i] = b[i] = c[i] = d[i] = cd();
        dft(a, n, 1); dft(b, n, 1); dft(c, n, 1); dft(d, n, 1);
        for(int i=0; i<n; i++) {
            cd _a = a[i], _b = b[i], _c = c[i], _d = d[i];
            a[i] = _a * _c;
            b[i] = _a * _d + _b * _c;
            c[i] = _b * _d;
        }
        dft(a, n, -1); dft(b, n, -1); dft(c, n, -1);
        for(int i=0; i<l; i++) {
            ll t = ( (ll(a[i].x + 0.5) %mo << 30) %mo + (ll(b[i].x + 0.5) %mo << 15) %mo + ll(c[i].x + 0.5) %mo) %mo;
            y[i] = (y[i] * 2 %mo - t + mo) %mo;
        }
    }
}

int n, a[N], b[N];
int main() {
    freopen("in", "r", stdin);
    fft::init();
    while( (n=read()) ) {
		memset(b, 0, sizeof(b)); memset(a, 0, sizeof(a));
        for(int i=1; i<=n; i++) b[i] = mo - read() % mo;
		b[0] = (b[0] + 1) %mo;
		int len = 1; while(len <= n) len <<= 1;
		fft::inverse(b, a, len);
        printf("%d
", a[n] %mo);
    }
}