Loj#6183. 看无可看

Loj#6183. 看无可看

题目描述
首先用特征根求出通项公式\(A_n=p\cdot 3^n+q\cdot(-1)^n\)。通过给定的\(f_0,f_1\)可以解出\(p,q\)。

然后我们要求的就是\(\sum_{|s'|=k}\Pi_{x\in s'}a_x\)。这就是个背包。

考虑它的生成函数就是\(\Pi(1+a_ix)\)。用分治\(FFT\)求解。

代码：

#include<bits/stdc++.h>
#define ll long long
#define N 100005

using namespace std;
inline int Get() {int x=0,f=1;char ch=getchar();while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}while('0'<=ch&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}

const ll mod=99991;
int Mod(int a) {return a<mod?a:a-mod;}
struct matrix {
    ll a[2][2];
    void Init() {memset(a,0,sizeof(a));}
}g,f;
matrix operator *(const matrix &x,const matrix &y) {
    static matrix tem;
    tem.Init();
    tem.a[0][0]=(x.a[0][0]*y.a[0][0]+x.a[0][1]*y.a[1][0])%mod;
    tem.a[0][1]=(x.a[0][0]*y.a[0][1]+x.a[0][1]*y.a[1][1])%mod;
    tem.a[1][0]=(x.a[1][0]*y.a[0][0]+x.a[1][1]*y.a[1][0])%mod;
    tem.a[1][1]=(x.a[1][0]*y.a[0][1]+x.a[1][1]*y.a[1][1])%mod;
    return tem;
}
matrix operator +(const matrix &x,const matrix &y) {
    static matrix tem;
    tem.Init();
    tem.a[0][0]=Mod(x.a[0][0]+y.a[0][0]);
    tem.a[0][1]=Mod(x.a[0][1]+y.a[0][1]);
    tem.a[1][0]=Mod(x.a[1][0]+y.a[1][0]);
    tem.a[1][1]=Mod(x.a[1][1]+y.a[1][1]);
    return tem;
}

matrix ksm(matrix g,int x) {
    static matrix ans;
    ans.Init();
    for(int i=0;i<2;i++) ans.a[i][i]=1;
    for(;x;x>>=1,g=g*g) if(x&1) ans=ans*g;
    return ans;
}

int n,k;
int a[N];
bool vis[N];
matrix G[N];
ll ans;
struct Com {
    double a,r;
    Com() {a=0,r=0;}
    Com(double x,double y) {a=x,r=y;}
};

Com operator +(const Com &x,const Com &y) {return Com(x.a+y.a,x.r+y.r);}
Com operator -(const Com &x,const Com &y) {return Com(x.a-y.a,x.r-y.r);}
Com operator *(const Com &x,const Com &y) {return Com(x.a*y.a-x.r*y.r,x.a*y.r+x.r*y.a);}
Com operator /(const Com &x,const double &y) {return Com(x.a/y,x.r/y);}
const double pi=acos(-1);
Com W[20][N<<2],M[20][N<<2];

void FFT(Com *a,int d,int flag) {
    static int rev[N<<2];
    int n=1<<d;
    for(int i=0;i<n;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<d-1);
    for(int i=0;i<n;i++) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int s=1;s<=d;s++) {
        int len=1<<s,mid=len>>1;
        for(int i=0;i<n;i+=len) {
            for(int j=0;j<mid;j++) {
                Com u=a[i+j],v=(flag==1)?a[i+j+mid]*W[d][j*n/len]:a[i+j+mid]*M[d][j*n/len];
                a[i+j]=u+v;
                a[i+j+mid]=u-v;
            }
        }
    }
    if(flag==-1) for(int i=0;i<n;i++) a[i]=a[i]/n;
}

ll ksm(ll t,ll x) {
    ll ans=1;
    for(;x;x>>=1,t=t*t%mod)
        if(x&1) ans=ans*t%mod;
    return ans;
}

ll F[N<<2];
Com L[N<<2],R[N<<2];
Com A[N<<2];
Com tem[N<<2];
ll t;

void solve(int l,int r) {
    if(l==r) {
        F[l]=ksm(t,a[l]);
        return ;
    }
    int mid=l+r>>1;
    solve(l,mid),solve(mid+1,r);
    int d=ceil(log2(r-l+2));
    for(int k=0;k<1<<d;k++)
        L[k]=R[k]=Com(0,0);
    for(int i=l;i<=mid;i++) {
        L[i-l+1]=Com(F[i],0);
    }
    for(int i=mid+1;i<=r;i++) {
        R[i-mid]=Com(F[i],0);
    }
    L[0]=R[0]=Com(1,0);
    int len=1<<d;
    FFT(L,d,1);
    FFT(R,d,1);
    for(int k=0;k<1<<d;k++) {
        tem[k]=(L[k]*R[k]);
    }
    FFT(tem,d,-1);
    for(int k=l;k<=r;k++) {
        F[k]=((ll)(tem[k-l+1].a+0.5))%mod;
    }
}

int main() {
    
    n=Get(),k=Get();
    g.a[0][0]=0,g.a[0][1]=3;
    g.a[1][0]=1,g.a[1][1]=2;
    int d=ceil(log2(n*2));
    for(int i=1;i<=d;i++)
        for(int j=0;j<1<<i;j++)
            W[i][j]=Com(cos(2*j*pi/(1<<i)),sin(2*j*pi/(1<<i))),
            M[i][j]=Com(cos(-2*j*pi/(1<<i)),sin(-2*j*pi/(1<<i)));
    for(int i=1;i<=n;i++) {
        a[i]=Get();
        G[i]=ksm(g,a[i]);
    }
    f.a[0][0]=Get(),f.a[0][1]=Get();
    ll q=(f.a[0][0]+f.a[0][1])*ksm(4,mod-2)%mod,p=(3*f.a[0][0]-f.a[0][1]+mod)*ksm(4,mod-2)%mod;

    ll ans=0;
    t=3;
    solve(1,n);
    ans=q*F[k]%mod;
    t=mod-1;
    solve(1,n);
    (ans+=p*F[k])%=mod;
    cout<<ans;
    return 0;
}

转载于:https://www.cnblogs.com/hchhch233/p/10458324.html