IOI2021集训队作业 191NF Frightful Formula
一个矩阵，给出第一行第一列，后面每个位置的数由此得到：\(F_{i,j}=F_{i,j-1}*a+F_{i-1,j}*b+c\)
求\(F_{n,n}\)。
答案对\(10^6+3\)取模。
\(n\le 2*10^5\)
如果\(c=0\)，考虑组合意义：对于从\((x,1)\)，贡献相当于\((x,1)\)到\((n,n)\)，每次可以向右或向下，向右有\(a\)的贡献，向下有\(b\)的贡献。第一步必须要向右。简单组合数即可计算。
如果\(c>0\)，类似地考虑除第一行第一列外的点，相当于它一开始有个\(c\)的贡献，向右或向下走到\((n,n)\)。这一部分显然可以列出式子：\(\sum_{x=0}^{n-2}\sum_{y=0}^{n-2}\binom{x+y}{x}a^xb^y\)。
这里其实可以直接卷积了。由于答案的模数比较小，所以可以直接FFT或大质数NTT（模\(29*2^{53}+1\)，原根为\(3\)）。
其实也不用。按照\(x+y\)的值分类计算，\(x+y\le n-2\)的贡献直接\(\sum_{s=0}^{n-2}(a+b)^s\)计算。
剩下的部分，考虑放到原来的题目中：先算出\(x+y=n-2\)的贡献，然后算下一个。乘\((a+b)\)相当于向右或向下走一格，那这个时候把走出矩阵的不合法贡献减掉即可。
NTT做法
using namespace std;
#include <cstdio>
#include <cstring>
#include <algorithm>
#define N 524288
#define ll long long
const int mo=1000003;
ll qpow(ll x,ll y=mo-2){
ll r=1;
for (;y;y>>=1,x=x*x%mo)
if (y&1)
r=r*x%mo;
return r;
}
ll fac[N*2],ifac[N*2];
void initC(int n){
fac[0]=1;
for (int i=1;i<=n;++i)
fac[i]=fac[i-1]*i%mo;
ifac[n]=qpow(fac[n]);
for (int i=n-1;i>=0;--i)
ifac[i]=ifac[i+1]*(i+1)%mo;
}
ll C(int m,int n){return fac[m]*ifac[n]%mo*ifac[m-n]%mo;}
int n,a,b,c;
ll A[N],B[N],s[N];
int p[N],q[N];
namespace NTT{
const ll mo=(29ll<<57)+1;
ll multi(ll x,ll y){
// return x*y%mo;
//x%=mo,y%=mo
x=(x+mo)%mo;
y=(y+mo)%mo;
ll d=(long double)x*y/mo;
d=x*y-d*mo;
if (d<0) d+=mo;
if (d>=mo) d-=mo;
// printf("%lld*%lld%%%lld=%lld\n",x,y,mo,d);
// printf("%lf\n",((double)x*y-d)/mo);
return d;
}
int re[N],nN;
ll qpow(ll x,ll y=mo-2){
ll r=1;
for (;y;y>>=1,x=multi(x,x))
if (y&1)
r=multi(r,x);
return r;
}
void setlen(int n){
int bit=0;
for (nN=1;nN<=n;++bit,nN<<=1);
re[0]=0;
for (int i=1;i<nN;++i)
re[i]=re[i>>1]>>1|(i&1)<<bit-1;
}
void dft(ll A[],int flag){
for (int i=0;i<nN;++i)
if (i<re[i])
swap(A[i],A[re[i]]);
static ll wnk[N];
for (int i=1;i<nN;i<<=1){
ll wn=qpow(3,flag==1?(mo-1)/(2*i):mo-1-(mo-1)/(2*i));
wnk[0]=1;
for (int k=1;k<i;++k)
wnk[k]=multi(wnk[k-1],wn);
for (int j=0;j<nN;j+=i<<1)
for (int k=0;k<i;++k){
ll x=A[j+k],y=multi(A[j+k+i],wnk[k]);
A[j+k]=(x+y)%mo;
A[j+k+i]=(x-y)%mo;
}
}
for (int i=0;i<nN;++i)
A[i]=(A[i]+mo)%mo;
if (flag==-1){
ll invn=qpow(nN);
for (int i=0;i<nN;++i)
A[i]=multi(A[i],invn);
}
}
void multi(ll c[],ll a[],ll b[]){
// for (int i=0;i<=n;++i)
// for (int j=0;j<=n;++j)
// (c[i+j]+=multi(a[i],b[j]))%=mo;
setlen(n*2);
static ll A[N],B[N];
for (int i=0;i<=n;++i)
A[i]=a[i],B[i]=b[i];
dft(A,1),dft(B,1);
for (int i=0;i<nN;++i)
A[i]=multi(A[i],B[i]);
dft(A,-1);
for (int i=0;i<=n*2;++i)
c[i]=A[i];
}
}
int main(){
// freopen("in.txt","r",stdin);
// freopen("out.txt","w",stdout);
scanf("%d%d%d%d",&n,&a,&b,&c);
A[0]=B[0]=1;
for (int i=1;i<=n;++i){
A[i]=A[i-1]*a%mo;
B[i]=B[i-1]*b%mo;
}
for (int i=1;i<=n;++i) scanf("%d",&p[i]);
for (int i=1;i<=n;++i) scanf("%d",&q[i]);
initC(n*2);
ll ans=0,sum=0;
for (int i=2;i<=n;++i) (ans+=A[n-1]*B[n-i]%mo*p[i]%mo*C(n-2+n-i,n-2))%=mo;
for (int i=2;i<=n;++i) (ans+=B[n-1]*A[n-i]%mo*q[i]%mo*C(n-2+n-i,n-2))%=mo;
n-=2;
for (int i=0;i<=n;++i){
A[i]=A[i]*ifac[i]%mo;
B[i]=B[i]*ifac[i]%mo;
}
// for (int i=0;i<=n;++i) printf("%lld ",A[i]);
// printf("\n");
// for (int i=0;i<=n;++i) printf("%lld ",B[i]);
// printf("\n");
NTT::multi(s,A,B);
// for (int i=0;i<=n*2;++i) printf("%lld ",s[i]);
// printf("\n");
for (int i=0;i<=n*2;++i)
(sum+=s[i]%mo*fac[i])%=mo;
(ans+=sum*c)%=mo;
printf("%lld\n",(long long)ans);
return 0;
}
