6508. 【GDOI2020模拟03.11】我的朋友们（多项式求逆、分治NTT）
题目描述
多项式求逆
第一次写就是这么毒瘤的题目
已知
\(A(x)B(x)\equiv 1(mod\;x^n)\)
要求
\(A(x)C(x)\equiv 1(mod\;x^{2n})\)
两式相减可得
\(A(x)(B(x)-C(x))\equiv 1(mod\;x^{n})\)
\(B(x)-C(x)\equiv 1(mod\;x^{n})\)
平方一下
\(B^2(x)+C^2(x)-2B(x)C(x)\equiv 0(mod\;x^{2n})\)，因为x0~x2n拆分后一定有一个在n以内的
两边同乘A(x)
\(A(x)B^2(x)+C(x)-2B(x)\equiv 0(mod\;x^{2n})\)
\(C(x)\equiv 2B(x)-A(x)B^2(x)(mod\;x^{2n})\)
虽然做了log n次乘法，但n的总和还是n级别的，所以还是\(O(n \log n)\)
分治FFT/NTT
一般形式（本题与此无关）：
\(f(i+j)=\sum{f(i)*g(j)}\)
cdq即可，log^2
题解
期望的一般算法：计算以某个状态为起点到终点的期望步数
设f[i]表示从i开始的答案，Pi(x)表示在滑块[i,i+L-1]中的期望多项式，其中x^k表示有k个不合法
那么可得
\(f[i]=\sum_{j=0}^{L}{[x^j]P_i(x)f[i+j]}\)
解方程可得
\(f[i]=\frac{1}{1-[x^0]P_i(x)}\sum_{j=1}^{L}{P_i(x)f[i+j]}\)
不好算，所以把数组p翻转，变成（Pi(x)变成[i-L+1,i]的期望）
\(f[i]=\frac{1}{1-[x^0]P_i(x)}\sum_{j=1}^{L}{P_i(x)f[i-j]}\)
接下来开始操♂作，定义
\(P[i]=p_ix+(1-p_i)\)（分清楚小括号和中括号），若i<1则P[i]=1
则\(P_i(x)=\prod_{i-L+1}^{i}{P[i]}\)
\(F_i(x)=\sum_{j<=i}{x^j*f(j)}\)
\(S_i(x)=\prod_{j<=i}{P[j]}\)
\(G[x][y]=\prod_{i=y-L+1}^{x}{P[i]}=S_x(x)/S_{y-L}(x)\)
（这里拓展了一下G的定义，为了方便之后的计算）
\(F[x][y]=F_{x-1}(x)*G[x][y]\)
如果能求出每个F[x][x]，那么\([x^x]F[x][x]\)的系数就是\(f[x]\)，字母不够用意会一下
\(F[1][n]=0\)，\(G[1][n]=S_1(x)/S_{n-L+1}(x)\)，直接乘会爆掉所以这里也要分治NTT
向下分治，设m=(x+y)/2
左半边：
\(G[x][m]=G[x][y]*\prod_{i=m-L+1}^{y-L}{P[i]}=G[x][y]*\prod_{i=m+1}^{y}{P[i-L]}\)
\(F[x][m]=F_{x-1}(x)*G[x][m]=F_{x-1}(x)*G[x][y]*\prod_{i=m+1}^{y}{P[i-L]}\)
\(=F[x][y]*\prod_{i=m+1}^{y}{P[i-L]}\)
有可能i-L<1，但是可以发现这样f必然为0，所以可以把P设为1
右半边：
\(G[m+1][y]=G[x][y]*\prod_{i=x+1}^{m+1}{P[i]}=G[x][y]*\prod_{i=x}^{m}{P[i+1]}\)
\(F[m+1][y]=F_{m}(x)*G[m+1][y]=F_{m}(x)*G[x][y]*\prod_{i=x}^{m}{P[i+1]}\)
为了方便接下来的操作，需要把Fm(x)的系数限制到x~y范围内
\(=F_{m}(x)*G[x][y]*\prod_{i=x}^{m}{P[i+1]}-F_{x-1}(x)*G[x][y]*\prod_{i=x}^{m}{P[i+1]}+前面那一坨\)
\(=(F_{m}(x)-F_{x-1}(x))*G[x][y]*\prod_{i=x}^{m}{P[i+1]}+F[x][y]*\prod_{i=x}^{m}{P[i+1]}\)
\(=F[x][y]*\prod_{i=x}^{m}{P[i+1]}+(F_{m}(x)-F_{x-1}(x))*G[m+1][y]\)
看似时间复杂度不变，但实际上F只需要保留以x^y结尾的2(y-x+1)项，G只用保留前(y-x+1)项即可，剩余的不会对新的区间产生影响
因为最终只需要x^x
还有一个要注意的，乘到F里的G[m+1][y]要保留(y-x+1)项
关键在于怎么实现，只可意会不可言传
一个比较方便的写法是写一个 把A(x)的一段乘上B(x)的一段后把某一段（删掉开头若干项）丢到C(x)上 的函数，要仔细考虑具体删掉多少项
极其难写+卡时间，写之前做好心理准备
code
我相信比expln之类的要难写的多（虽然还没学）
#include <bits/stdc++.h>
#define fo(a,b,c) for (a=b; a<=c; a++)
#define fd(a,b,c) for (a=b; a>=c; a--)
#define max(a,b) (a>b?a:b)
#define min(a,b) (a<b?a:b)
#define mod 998244353
#define Mod 998244351
#define ll long long
#define G 3
#define file
using namespace std;
ll aa[262144],bb[262144],bb2[262144],bb3[262144],cc[262144],A[262144],p[262144],p1[4194304],p2[4194304]; //p1=i-l p2=i+1
ll f[8388608],g[4194304],a[8388608],b[8388608],c[262144],ans[262144],P[262144]; //f=F g=G ans=f
int bg[262144],ed[262144],bg1[262144],ed1[262144],bg2[262144],ed2[262144],A2[19][262144],len2[262145],n2[262145],n3[262145],W1[19],W2[19],n,L,I,J,i,j,k,l,s,s2,N,N2,len,tot,tot1,tot2;
void look(ll *a,int x,int y)
{
int i;
fd(i,y,x) cout<<(a[i]+mod)%mod<<" ";cout<<endl;
}
ll qpower(ll a,int b)
{
ll ans=1;
while (b)
{
if (b&1)
ans=ans*a%mod;
a=a*a%mod;
b>>=1;
}
return ans;
}
void dft(ll *a,int type,int N,int len)
{
int i,j,k,l,S=N,s1=2,s2=1;
ll w,W,u,v;
fo(i,0,N-1) A[A2[len][i]]=a[i];
memcpy(a,A,N*8);
fo(i,1,len)
{
S>>=1;
w=(type==1)?W1[i]:W2[i];
fo(j,0,S-1)
{
W=1;
fo(k,0,s2-1)
{
u=a[j*s1+k];
v=a[j*s1+k+s2]*W;
a[j*s1+k]=(u+v)%mod;
a[j*s1+k+s2]=(u-v)%mod;
W=W*w%mod;
}
}
s1<<=1,s2<<=1;
}
}
void mul(ll *x,ll *y,int n,int m)
{
int i,j,k,l,N,N2,len;
len=len2[n+m];N=n2[n+m];N2=n3[n+m];
if (n<m)
fo(i,n,m) x[i]=0;
else
fo(i,m,n) y[i]=0;
fo(i,max(n,m),N-1) x[i]=y[i]=0;
memcpy(bb3,y,N*8);
dft(x,1,N,len);
dft(bb3,1,N,len);
fo(i,0,N-1) x[i]=x[i]*bb3[i]%mod;
dft(x,-1,N,len);
fo(i,0,N-1) x[i]=x[i]*N2%mod;
}
void Mul(ll *x,ll *y,ll *z,int X1,int X2,int Y1,int Y2,int Z1,int Z2,int st)
{
int i;
fo(i,X1,X2) aa[i-X1]=x[i];
fo(i,Y1,Y2) bb[i-Y1]=y[i];
mul(aa,bb,X2-X1+1,Y2-Y1+1);
fo(i,Z1,Z2) z[i]=aa[i-Z1+st];
}
void ny(ll *a,int n)
{
int i,j,k,l;
cc[0]=qpower(a[0],Mod);
for (i=2; i<=n; i*=2)
{
memcpy(bb2,a,i*8);
mul(bb2,cc,i,i);mul(bb2,cc,i,i);
fo(j,0,i+i-1)
bb2[j]=(2*cc[j]-bb2[j])%mod;
memcpy(cc,bb2,i*8);
}
memcpy(a,cc,n*8);
}
void work(int t,int x,int y)
{
int i,mid=(x+y)/2;
if (x==y)
{
bg[t]=tot+1;ed[t]=tot+2;
if (x>L)
p1[tot+1]=1-p[x-L],p1[tot+2]=p[x-L];
else
p1[tot+1]=1,p1[tot+2]=0;
if (x<n)
p2[tot+1]=1-p[x+1],p2[tot+2]=p[x+1];
else
p2[tot+1]=1,p2[tot+2]=0;
tot+=2;
return;
}
work(t*2,x,mid);
work(t*2+1,mid+1,y);
bg[t]=tot+1;ed[t]=tot+(y-x+1)+1;
Mul(p1,p1,p1,bg[t*2],ed[t*2],bg[t*2+1],ed[t*2+1],bg[t],ed[t],0);
Mul(p2,p2,p2,bg[t*2],ed[t*2],bg[t*2+1],ed[t*2+1],bg[t],ed[t],0);
tot+=(y-x+1)+1;
}
void Work(int t,int x,int y)
{
int i,mid=(x+y)/2;
if (x==y)
{
if (x>=L)
ans[x]=(f[bg1[t]+1]+1)*qpower(1-P[x],Mod)%mod;
return;
}
bg2[t*2]=tot2+1;ed2[t*2]=tot2+(mid-x+1);tot2+=mid-x+1;
Mul(g,p1,g,bg2[t],ed2[t],bg[t*2+1],ed[t*2+1],bg2[t*2],ed2[t*2],0);
bg1[t*2]=tot1+1;ed1[t*2]=tot1+min(n+1,2*(mid-x+1));tot1+=min(n+1,2*(mid-x+1));
Mul(f,p1,f,bg1[t],ed1[t],bg[t*2+1],ed[t*2+1],bg1[t*2],ed1[t*2],max(0,2*x-mid-1)-max(0,2*x-y-1));
Work(t*2,x,mid);
// ---
bg2[t*2+1]=tot2+1;ed2[t*2+1]=tot2+(y-mid);tot2+=y-mid;
Mul(g,p2,a,bg2[t],ed2[t],bg[t*2],ed[t*2],bg2[t],ed2[t],0); //Ïȱ£Áôy-x+1Ïî
fo(i,0,y-mid-1) g[bg2[t*2+1]+i]=a[bg2[t]+i];
bg1[t*2+1]=tot1+1;ed1[t*2+1]=tot1+min(n+1,2*(y-mid));tot1+=min(n+1,2*(y-mid));
Mul(f,p2,b,bg1[t],ed1[t],bg[t*2],ed[t*2],bg1[t*2+1],ed1[t*2+1],max(0,2*mid-y+1)-max(0,2*x-y-1));
fo(i,x,mid) c[i]=ans[i];
Mul(c,a,a,x,mid,bg2[t],ed2[t],bg1[t*2+1],ed1[t*2+1],max(0,2*mid-y+1)-x);
fo(i,bg1[t*2+1],ed1[t*2+1]) f[i]=(a[i]+b[i])%mod;
Work(t*2+1,mid+1,y);
}
int main()
{
freopen("friends.in","r",stdin);
#ifdef file
freopen("friends.out","w",stdout);
#endif
I=2;J=1;s=2;len2[1]=1;n2[1]=2;n3[1]=499122177;
fo(len,0,18)
{
if (len<18)
W1[len+1]=qpower(G,(mod-1)/s),W2[len+1]=qpower(W1[len+1],Mod);
s2=qpower(s,Mod);
fo(i,0,J-1)
{
if (I<=262144)
len2[I]=len+1,n2[I]=s,n3[I]=s2,++I;
j=i;k=0;
fo(l,1,len)
k=(k<<1)+(j&1),j>>=1;
A2[len][i]=k;
}
J*=2;s=s*2;
}
scanf("%d%d",&n,&L);
fd(i,n,1)
{
scanf("%d%d",&p[i],&j);
p[i]=p[i]*qpower(j,Mod)%mod;
}
P[L]=1;
fo(i,1,L) P[L]=P[L]*(1-p[i])%mod;
fo(i,L+1,n) P[i]=P[i-1]*(1-p[i])%mod*qpower(1-p[i-L],Mod)%mod;
work(1,1,n);
tot1=n+1;bg1[1]=1;ed1[1]=n+1;
tot2=n+1;bg2[1]=1;ed2[1]=n;
a[1]=p[1],a[0]=1-p[1];
fo(i,0,n)
b[i]=p1[bg[1]+i];
N=pow(2,ceil(log2(n+1)));
ny(b,N);
Mul(a,b,g,0,1,0,N-1,1,n,0);
Work(1,1,n);
printf("%lld\n",(ans[n]+mod)%mod);
fclose(stdin);
fclose(stdout);
return 0;
}
