UOJ #450「集训队作业2018」复读机

UOJ #450

---

题意

有$ k$台复读机，每时每刻有且只有一台复读机进行复读

求$ n$时刻后每台复读机的复读次数都是$ d$的倍数的方案数

$ 1\leq d \leq 3,k \leq 5·10^5,n \leq 10^9$ 

当$ d=3$时$ k \leq 10^3$

---

题解

$ d=1$的略过

对复读机构建生成函数

发现这是指数生成函数

即我们要计算的是$$(\sum_{i=0}^n[d|i]\frac{x^i}{i!})^k[x^n]=(\sum_{i=0}^n[d|i]e^x)^k [x^n]$$

直接单位根反演，将后式化简为$$(\frac{1}{d}\sum_{i=0}^{d-1}e^{w_d^ix})^k[x^n]$$

当$ d=2$时我们实际要求的就是$$(\frac{e^x+e^{-x}}{2})^k[x^n]$$

直接二项式展开即可

当$ d=3$时用类似的做法三项式展开即可

时间复杂度$ O(k^{d-1}·\log)$

---

代码 

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<queue>
#include<vector>
#define p 19491001
#define rt register int
#define ll long long
using namespace std;
inline ll read(){
    ll x=0;char zf=1;char ch=getchar();
    while(ch!='-'&&!isdigit(ch))ch=getchar();
    if(ch=='-')zf=-1,ch=getchar();
    while(isdigit(ch))x=x*10+ch-'0',ch=getchar();return x*zf;
}
void write(ll y){if(y<0)putchar('-'),y=-y;if(y>9)write(y/10);putchar(y%10+48);}
void writeln(const ll y){write(y);putchar('\n');}
int k,m,n,x,y,z,cnt,ans,d;
int inv[500010],jc[500010],njc[500010];
int ksm(int x,int y=p-2){
    x=(x%p+p)%p;
    int ans=1;
    for(rt i=y;i;i>>=1,x=1ll*x*x%p)if(i&1)ans=1ll*ans*x%p;
    return ans;
}
inline int C(int x,int y){
    return 1ll*jc[x]*njc[y]%p*njc[x-y]%p;
}
namespace subtask1{
    void main(){
        write(ksm(k,n));
        return;
    }
}
namespace subtask2{
    void main(){
        int ans=0;
        for(rt i=0;i<=k;i++)(ans+=1ll*C(k,i)*ksm(i+i-k,n)%p)%=p;
        ans=1ll*ans*ksm(ksm(2,k))%p;
        writeln(ans);
        return;
    }
}
namespace subtask3{
    void main(){
        int ans=0;
        int w[3];w[0]=1;w[1]=ksm(7,(p-1)/3);w[2]=1ll*w[1]*w[1]%p;
        for(rt i=0;i<=k;i++)
        for(rt j=0;i+j<=k;j++){
            int xs=(1ll*w[0]*i%p+1ll*w[1]*j%p+1ll*w[2]*(k-i-j)%p)%p;
            xs=ksm(xs,n);
            (ans+=1ll*C(k,i)*C(k-i,j)%p*xs%p)%=p;
        }
        ans=1ll*ans*ksm(ksm(3,k))%p;
        writeln(ans);
        return;
    }
}
int main(){
    cin>>n>>k>>d;
    for(rt i=0;i<=1;i++)jc[i]=njc[i]=inv[i]=1;
    for(rt i=2;i<=k;i++){
        jc[i]=1ll*jc[i-1]*i%p;
        inv[i]=1ll*inv[p%i]*(p-p/i)%p;
        njc[i]=1ll*njc[i-1]*inv[i]%p;
    } 
    if(d==1)subtask1::main();
    if(d==2)subtask2::main();
    if(d==3)subtask3::main();
    return 0;
}