[CQOI 2018] 交错序列
\(\text{I love Wordle!!!}\)
I guessed this 5-letter word in 5/6 tries.
?????
??????
????????
?????????
??????????
Can you guess this word?
\(\frak{Description}\)
\(\rm Link.\)
\(\frak{Solution}\)
\(\text{Method 1}\)：喵喵 \(\mathtt{dp}\)
对于 \(0\) 取 \(x\) 次，\(1\) 取 \(y\) 次的序列有
\[x^ay^b=(n-y)^ay^b=\sum_{i=0}^a \binom{a}{i}\cdot n^i(-y)^{a-i}\cdot y^b=\sum_{i=0}^a \binom{a}{i}\cdot n^i\cdot (-1)^{a-i}\cdot y^{a+b-i} \]只要我们能快速计算出 所有序列的 \(y\) 的 \(b\) 到 \(a+b\) 次幂，就能以 \(\mathcal O(a)\) 的复杂度计算答案。所以设 \(dp_{i,j,0/1}\) 表示填了 \(i\) 个数字，末尾为 \(0/1\)，所有序列的 \(y\) 的 \(j\) 次方和。
那么对于末尾为零的情况（因为对序列 \(1\) 的个数没有影响）
\[dp_{i,j,0}=dp_{i-1,j,0}+dp_{i-1,j,1} \]对于末尾为 \(1\) 的情况
\[dp_{i,j,1}=\sum_{k=0}^j \binom{j}{k}\cdot dp_{i-1,k,0} \]这是因为对于所有序列，\(y\leftarrow y+1\)，由二项式定理 \((y+1)^j=\sum_{k=0}^j \binom{j}{k}\cdot y^k\)，我们可以得到这个递推式。这个玩意可以矩阵加速，做到 \(\mathcal O((2a+2b)^3\cdot \log n)\).
\(\text{Method 2}\)：我的乱搞
如果你是个文明人，千万不要点开下面这段话：
博客园竟然开始检索敏感词了，为了不辜负点开 \(\rm details\) 的小伙伴们，就放在剪贴板里了：\(\rm Link.\)
非常神奇的是，在卡了很久的常无果后，在 long long 前加上 \(\rm unsigned\) 就快到飞起了……开心??。
受到之前做的排列 \(\mathtt{dp}\) 专题的深刻影响，这题一来我就想枚举 \(y\)，然后计算 \(1\) 个数为 \(y\) 的序列方案数，那么答案就是
\[\text{Ans}=\sum_{i=0}^{(n+1)/2}\binom{n-i+1}{i}\cdot (n-i)^a\cdot i^b \]组合数的意义就是在 \(n-i\) 个零的序列里插入 \(i\) 个 \(1\)，显然有 \(n-i+1\) 个位置。用欧拉筛筛出两个函数 \(f(x)=x^a,g(x)=x^b\)，这题不就解决了吗？伞兵题！（？
写到一半突然发现可能出现 \(n\ge m\)，那么就不存在逆元，没办法算组合数。然后我就蒙了，难不成用 \(\rm lucas\) 定理？复杂度直接乘上 \(10^7\) 的 \(\log\) 啊！
不行啊，难道我就要这样屈服吗！观察到组合数是极有规律的，从 \(\binom{n+1}{0}\) 到 \(\binom{n}{1}\)，以此类推……注意到每次组合数 \(\binom{n}{m}\) 都是 \(n\leftarrow n-1, m\leftarrow m+1\)，我们可以在这上面搞些名堂。还是利用 \(\rm lucas\) 定理，它实质上是将组合数 \(\binom{n}{m}\) 中的 \(n,m\) 分别分解成两个 \(\rm mod\) 进制数，对于每一位形成一个新的组合数，然后再乘起来。根据上文的观察，我们发现 \(n,m\) 的操作和 "" 非常类似（不要脸的硬行植入广告），使用 "累计分析" 就可以很好的分析这个复杂度（不妨令 \(N=(n+1)/2\)）：对于 \(m\)，第 \(0\) 位变化 \(N\) 次，第 \(1\) 位变化 \(N/\text{mod}\) 次，依此类推……显然总变化次数的上界是 \(2N\)；对于 \(n\)，它的总变化次数实际上是 \(0\rightarrow n\) 的次数减去 \(0\rightarrow n-N\)（有点粗糙）的次数，显然这也是 \(\mathcal O(n)\) 级别的。
所以直接模拟类二进制计数器递推组合数就可以做到 \(\mathcal O(n)\) 的时间复杂度。因为卡了常，所以代码会有些长（
\(\frak{Code}\)
\(\text{Method 1}\)
#include
#define print(x,y) write(x),putchar(y)
template
inline T read(const T sample) {
T x=0; char s; bool f=0;
while((s=getchar())>'9' || s<'0')
f |= (s=='-');
while(s>='0' && s<='9')
x = (x<<1)+(x<<3)+(s^48),
s = getchar();
return f?-x:x;
}
template
inline void write(T x) {
static int writ[50], w_tp=0;
if(x<0) putchar('-'), x=-x;
do writ[++w_tp]=x-x/10*10,x/=10; while(x);
while(putchar(writ[w_tp--]^48),w_tp);
}
#include
typedef long long ll;
int n,a,b,mod,C[100][100];
struct mat {
int n,m; ll a[182][182];
mat() { memset(a,0,sizeof a); }
mat operator * (const mat& t) const {
mat r; r.n=n, r.m=t.m;
for(int i=0;i<=n;++i) for(int j=0;j<=m;++j)
if(a[i][j]) for(int k=0;k<=r.m;++k)
r.a[i][k] = r.a[i][k]+a[i][j]*t.a[j][k]%mod;
for(int i=0;i<=r.n;++i) for(int j=0;j<=r.m;++j)
r.a[i][j] = r.a[i][j]%mod;
return r;
}
} per,A;
void preWork() {
C[0][0]=1;
for(int i=1;i<=91;++i) {
C[i][0]=1;
for(int j=1;j<=i;++j)
C[i][j] = (C[i-1][j-1]+C[i-1][j])%mod;
}
per.n = per.m = (a+b<<1)+1;
for(int i=0;i<=a+b;++i) {
per.a[i][i] = per.a[i+a+b+1][i] = 1;
for(int j=i;j<=a+b;++j)
per.a[i][j+a+b+1] = C[j][i];
}
}
mat qkpow(int y) {
mat r; r.n=r.m=per.n;
for(int i=0;i<=r.n;++i) r.a[i][i]=1;
while(y) {
if(y&1) r = r*per;
per = per*per; y>>=1;
}
return r;
}
inline int qkpow(int x,int y) {
int r=1;
for(; y; y>>=1, x=1ll*x*x%mod)
if(y&1) r=1ll*r*x%mod;
return r;
}
int main() {
n=read(9), a=read(9), b=read(9), mod=read(9);
preWork(); A.n=0, A.m=per.m; A.a[0][0]=1;
A = A*qkpow(n);
int ans=0, pw=1;
for(int i=0;i<=a;++i) {
int tmp = 1ll*(A.a[0][a+b-i]+A.a[0][a+b-i+a+b+1])%mod*pw%mod*C[a][i]%mod;
pw = 1ll*pw*n%mod;
if((a-i)&1) ans = ((ans-tmp)%mod+mod)%mod;
else ans = (ans+tmp)%mod;
}
print(ans,'\n');
return 0;
}
\(\text{Method 2}\)
#include
#define print(x,y) write(x),putchar(y)
template
inline T read(const T sample) {
T x=0; char s; bool f=0;
while((s=getchar())>'9' || s<'0')
f |= (s=='-');
while(s>='0' && s<='9')
x = (x<<1)+(x<<3)+(s^48),
s = getchar();
return f?-x:x;
}
template
inline void write(T x) {
static int writ[50], w_tp=0;
if(x<0) putchar('-'), x=-x;
do writ[++w_tp]=x-x/10*10,x/=10; while(x);
while(putchar(writ[w_tp--]^48),w_tp);
}
#include
using namespace std;
typedef unsigned long long ll;
const int maxn = 1e7+5;
bool is[maxn];
ll cnm=1,fac[maxn],ifac[maxn];
int f[101],g[101],all,ban;
int A[maxn],B[maxn],p[(int)1e6+1],pc;
int n,a,b,mod;
inline int inv(int x,int y=mod-2) {
int r=1;
for(; y; y>>=1, x=1ll*x*x%mod)
if(y&1) r=1ll*r*x%mod;
return r;
}
void preWork() {
fac[0]=ifac[0]=1; int lim = min(n+1,mod-1);
for(int i=1;i<=lim;++i)
fac[i] = 1ll*fac[i-1]*i%mod;
ifac[lim] = inv(fac[lim]);
for(int i=lim-1;i;--i)
ifac[i] = 1ll*ifac[i+1]*(i+1)%mod;
}
void sieve() {
A[1]=B[1]=1;
A[0] = (!a), B[0] = (!b);
for(int i=2;i<=n;++i) {
if(!is[i]) {
p[++pc] = i;
A[i] = inv(i,a), B[i] = inv(i,b);
}
for(int j=1; j<=pc && i*p[j]<=n; ++j) {
is[i*p[j]] = 1;
A[i*p[j]] = 1ll*A[i]*A[p[j]]%mod;
B[i*p[j]] = 1ll*B[i]*B[p[j]]%mod;
if(i%p[j]==0) break;
}
}
}
inline int QuickMod(const ll& n) {
return n-mod*(n/mod);
}
inline ll plainC(const int& n,const int& m) {
return QuickMod(ifac[n-m]*QuickMod(fac[n]*ifac[m]));
}
inline int invC(const int& n,const int& m) {
return QuickMod(ifac[n]*QuickMod(fac[m]*fac[n-m]));
}
inline void adj() {
register int add1=0, add2=0;
register int pre1 = (f[1]++), pre2 = (g[1]--);
f[1] = (f[1]>=mod? (add1=1, f[1]-mod): f[1]),
g[1] = (g[1]<0? (add2=1, g[1]+mod): g[1]);
register int pd1 = (pre1<=pre2), pd2 = (f[1]<=g[1]);
if(pd1 && (!pd2)) {
++ban,
cnm = cnm*invC(pre2,pre1)%mod;
} else if((!pd1) && pd2) {
--ban,
cnm = cnm*plainC(g[1],f[1])%mod;
} else if(pd1)
cnm = cnm*plainC(g[1],f[1])%mod*invC(pre2,pre1)%mod;
for(register int i(2); i<=all && (add1 || add2); ++i) {
pre1 = f[i], pre2 = g[i];
if(add1) ++ f[i], add1=0;
if(add2) -- g[i], add2=0;
if(f[i]>=mod) f[i] = f[i]-mod, add1=1;
if(g[i]<0) g[i] = g[i]+mod, add2=1;
if(pre1<=pre2 && f[i]>g[i]) {
++ban;
cnm = cnm*invC(pre2,pre1)%mod;
} else if(pre1>pre2 && f[i]<=g[i]) {
--ban;
cnm = cnm*plainC(g[i],f[i])%mod;
} else if(pre1<=pre2) {
cnm = cnm*plainC(g[i],f[i])%mod*invC(pre2,pre1)%mod;
}
}
}
void div(int x) {
while(x) g[++all]=x%mod, x/=mod;
}
int main() {
n=read(9), a=read(9), b=read(9), mod=read(9);
preWork(); sieve(); div(n+1);
int up=(n+1>>1); register ll ans=0;
for(register int i(0);i<=up;++i) {
if(!ban) ans = ans+cnm*A[n-i]%mod*B[i]%mod;
adj();
}
print(ans%mod,'\n');
return 0;
}