2018.01.02【CQOI2018】【BZOJ5300】【洛谷P4461】九连环（数学推理）（FFT加速高精乘）

BZOJ传送门

洛谷传送门

---

解析：

莫名其妙就BZOJ本题rk1了啊。。。

通过OEIS, 找规律 数学推理我们可以得到一个这样的递推式：
$$\{\begin{array}{r}{f}_0=0,f_1=1\\ f_n+f_{n-1}=2^n\end{array}$$

通过各种奇怪的方法我们都可以求解（比如二进制拆分），解出来是这个东西：$$f_n=\lfloor \frac{2^{n+1}}{3}\rfloor$$

然而毒瘤出题人并没有让我们取模。。。

没办法了上高精度吧。。
乘法可以压位搞一搞，不过我还是选择了压位FFT。。。

顺便优化了一下我那本来常数不小的FFT。

---

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define re register
#define gc getchar
#define pc putchar
#define cs const

namespace IO{
    namespace IOONLY{
        cs int Rlen=1<<18|1;
        char buf[Rlen],*p1,*p2;
        char obuf[Rlen],*p3=obuf;
        char ch[23];
    }
    inline char get_char(){
        using namespace IOONLY;
        return (p1==p2)&&(p2=(p1=buf)+fread(buf,1,Rlen,stdin),p1==p2)?EOF:*p1++;
    }
    inline void put_char(char c){
        using namespace IOONLY;
        *p3++=c;
        if(p3==obuf+Rlen)fwrite(obuf,1,Rlen,stdout),p3=obuf;
    }
    inline void FLUSH(){
        using namespace IOONLY;
        fwrite(obuf,1,p3-obuf,stdout);p3=obuf;
    }
    
    inline int getint(){
        re int num;
        re char c;
        while(!isdigit(c=gc()));num=c^48;
        while(isdigit(c=gc()))num=(num+(num<<2)<<1)+(c^48);
        return num;
    }
    inline void outint(int a,int minspace=1){
        using namespace IOONLY;
        do ch[++ch[0]]=a-a/10*10,a/=10; while(a);
        while(ch[0]<minspace)ch[++ch[0]]=0;
        while(ch[0])pc(ch[ch[0]--]^48);
    }
}
using namespace IO;

cs int N=50004;

cs double PI=acos(-1.0);
struct Complex{
    double x,y;
    Complex(){}
    Complex(cs double &_x,cs double &_y):x(_x),y(_y){}
    inline friend Complex operator+(cs Complex &a,cs Complex &b){
        return Complex(a.x+b.x,a.y+b.y);
    }
    inline friend Complex operator-(cs Complex &a,cs Complex &b){
        return Complex(a.x-b.x,a.y-b.y);
    }
    inline friend Complex operator*(cs Complex &a,cs Complex &b){
        return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);
    }
}omg[N],inv[N],A[N],B[N];
int r[N],lim;

inline void init(int n){
    lim=0;
    while((1<<lim)<n)++lim;
    for(int re i=0;i<n;++i){
        omg[i]=Complex(cos(2*PI/n*i),sin(2*PI/n*i));
        inv[i]=Complex(omg[i].x,-omg[i].y);
        r[i]=(r[i>>1]>>1)|((i&1)<<(lim-1));
    }
}

cs int Base=1000,Bpow=3;

inline void FFT(Complex *A,int n,Complex *omg){
    for(int re i=0;i<n;++i)if(i<r[i])swap(A[i],A[r[i]]);
    for(int re l=1;l<n;l<<=1){
        for(Complex *p=A;p!=A+n;p+=l<<1){
            for(int re i=0;i<l;++i){
                Complex tmp=omg[n/l/2*i]*p[i+l];
                p[i+l]=p[i]-tmp;
                p[i]=p[i]+tmp;
            }
        }
    }
}

struct BIGNUM{
    int num[N],len;
    BIGNUM(){}
    void operator=(int a){
        len=0;memset(num,0,sizeof num);
        if(a==0)return void(len=1);
        while(a)num[len++]=a-a/Base*Base,a/=Base;
    }
    
    void print(){
        for(int re i=len-1;~i;--i){
            outint(num[i],i==len-1?1:Bpow);
        }
        pc('\n');
    }
    
    void operator/=(cs int &x){
        re int sum=0,newlen=0;
        for(int re i=len-1;~i;--i){
            sum=sum*Base+num[i];
            if(sum<x)num[i]=0;
            else{
                if(!newlen)newlen=i+1;
                num[i]=sum/x;
                sum-=sum/x*x; 
            }
        }
        len=max(newlen,1);
    }
    
    void operator*=(cs BIGNUM &b){
        int newlen=len+b.len-1,n=1;
        while(n<newlen)n<<=1;
        init(n);
        for(int re i=0;i<n;++i){
            A[i]=Complex(i<len?num[i]:0,0);
            B[i]=Complex(i<b.len?b.num[i]:0,0);
        }
        FFT(A,n,omg);
        FFT(B,n,omg);
        for(int re i=0;i<n;++i)A[i]=A[i]*B[i];
        FFT(A,n,inv);
        for(int re i=0;i<newlen;++i)num[i]=(int)floor(A[i].x/n+0.5);
        num[len=newlen]=0;
        for(int re i=0;i<len;++i)num[i+1]+=num[i]/Base,num[i]%=Base;
        if(num[len])++len;
    }
    
}res[11],a;

int m,n[11];
signed main(){
    m=getint();
    for(int re i=1;i<=m;++i){
        res[i]=1;
        n[i]=getint()+1;
    }
    a=2;
    int re l=0;
    while(true){
        bool flag=false;
        for(int re i=1;i<=m;++i){
            if(n[i]&(1<<l)){
                res[i]*=a;
            }
            if(n[i]>(1<<l))flag=true;
        }
        if(!flag)break;
        a*=a;++l;
    }
    for(int re i=1;i<=m;++i)res[i]/=3,res[i].print();
    FLUSH();
    return 0;
}