LG4351 [CERC2015]Frightful Formula

Frightful Formula

给你一个(n imes n)矩阵的第一行和第一列，其余的数通过如下公式推出：

[f_{i,j}=acdot f_{i,j-1}+bcdot f_{i-1,j}+c ]

求(f_{n,n}mod (10^6+3))。

(2≤n≤200000)

zhouzhendong的题解

利用FFT来解决此题

这题的主要思路很巧妙。

先假设(a=b=1)

我们考虑对于每一个(f_{i,j})计算它对(f_{n,n})的贡献次数。

显然贡献次数就是在网格图中从((i,j))开始只能往右或者往下走的路径条数，即(inom{n-i+n-j}{n-i})。

注意第一行和第一列的对答案的贡献式有点不同。

因为第一行的格子不能对右边的格子做贡献，所以第一行的格子((1,i))的贡献次数就是(inom{2n-i-2}{n-i})。即第一步必须往下走。

同理，第一列的格子((i,1))的贡献次数就是(inom{2n-i-2}{n-i})。

现在把(a,b)考虑进去

我们还要乘一个权值。考虑从第(i)行走到第(n)行向下走了(n-i)次，所以乘了(n-i次b)，即乘了(b^{n-i})。同理向右走要乘的权值就是(a^{n-j})。所以综上所述，我们再重新写一下贡献次数的式子。

下面的三式满足(i,j>1)。注意((1,1))是没有贡献次数的。

[(1,i) ightarrow inom{2n-i-2}{n-i}a^{n-i}b^{n-1}\ (i,1) ightarrow inom{2n-i-2}{n-i}a^{n-1}b^{n-i}\ (i,j) ightarrow inom{2n-i-j}{n-i}a^{n-j}b^{n-i} ]

考虑到((i,j))这种格子的值里面含有((i,1))和((1,i))这种的贡献，贡献不能重复计算，所以我们在计算格子((i,j))的时候就只考虑它产生的(c)对答案的贡献次数。（就是上面的那个）

所以，我们可以列出答案的式子：

[ans=sum_{i=2}^{n}f_{i,1}cdotinom{2n-i-2}{n-i}a^{n-1}b^{n-i}\ +sum_{i=2}^{n}f_{1,i}cdotinom{2n-i-2}{n-i}a^{n-i}b^{n-1}\ +sum_{i=2}^{n}sum_{j=2}^{n}ccdotinom{2n-i-j}{n-i}a^{n-i}b^{n-j} ]

各种逆元、阶乘以及幂的预处理就不说了。于是前面的两个式子都可以(O(n))搞定。

后面的那个式子，我们再给他稍稍变几个形：

[sum_{i=2}^{n}sum_{j=2}^{n}ccdotinom{2n-i-j}{n-i}a^{n-i}b^{n-j} ]

把(c)提前，并修改求和指标，枚举(n-i)和(n-j)，得：

[=csum_{i=0}^{n-2}sum_{j=0}^{n-2}frac{(i+j)!}{i!j!}a^ib^j\ =csum_{i=0}^{n-2}i!sum_{j=0}^{n-2}frac{a^{(i-j)}}{(i-j)!}cdotfrac{b^j}{j!} ]

设

[f_i=egin{cases}large{frac{a^i}{i!}}& ext{$(0leq ileq n-2)$}\0& ext{$(n-2<i)$}end{cases}\ g_i=egin{cases}large{frac{b^i}{i!}}& ext{$(0leq ileq n-2)$}\0& ext{$(n-2<i)$}end{cases} ]

则原式可以写为：

[csum_{i=0}^{n-2}(i+j)!sum_{j=0}^{n-2}f_ig_j ]

发现这个符合多项式卷积形式。于是我们再换一个写法：

[=csum_{i=0}^{2n-4}i!sum_{j=0}^{i}f_jg_{i-j} ]

于是显然可以FFT快速计算了。注意，要取模，而且卷积得到的数字较大，可以用任意模数FFT或者拆系数FFT。时间复杂度(O(nlog n))。

线性递推

我们把之前的式子继续推下去。

[csum_{i=0}^{n-2}sum_{j=0}^{n-2}frac{(i+j)!}{i!j!}a^ib^j ]

我们只需要考虑后面的式子，先不看(c)。

[sum_{i=0}^{n-2}sum_{j=0}^{n-2}inom{i+j}{i}a^ib^j ]

设

[f_n=sum_{i=0}^{n}sum_{j=0}^{n}inom{i+j}{i}a^ib^j ]

则我们要求的就是(f_{n-2})。

推导：

[f_n=sum_{i=0}^{n}sum_{j=0}^{n}inom{i+j}{i}a^ib^j\ =sum_{i=0}^{n}a^isum_{j=0}^{n}inom{i+j}{i}b^j\ =sum_{i=0}^{n-1}a^isum_{j=0}^{n-1}inom{i+j}{i}b^j+a^nsum_{j=0}^{n}inom{n+j}{n}b^j+b^nsum_{i=0}^{n}a^iinom{i+n}{i}-inom{2n}{n}a^nb^n\ =f_{n-1}+a^nsum_{i=0}^{n}inom{n+i}{n}b^i+b^nsum_{i=0}^{n}inom{n+i}{n}a^i-inom{2n}{n}a^nb^n ]

我们只需要处理(sum_{i=0}^{n}inom{n+i}{n}b^i)和(sum_{i=0}^{n}inom{n+i}{n}a^i)，就可以搞定了。

由于这两个形式一样，我这里只对a进行推导。

设

[ga_n=sum_{i=0}^{n}inom{n+i}{n}a^i ]

则

[ga_n=sum_{i=0}^ninom{n+i}{i}a^i \ = 1+sum_{i=1}^{n-1}left( inom{n-1+i}{i}+inom{n-1+i}{i-1} ight) a^i + inom{2n}{n}a^n\ = 1+ga_{n-1}-1+sum_{i=0}^{n-2}inom{n+i}{i}a^{i+1}+inom{2n}{n}a^n\ =ga_{n-1}+acdot ga_n-inom{2n-1}{n-1}a^n-inom{2n}{n}a^{n+1}+inom{2n}{n}a^n\ =ga_{n-1}+acdot ga_n+inom{2n-1}{n}a^n-inom{2n}{n}a^{n+1} ]

移项

[(1-a)ga_n=ga_{n-1}+inom{2n-1}{n}a^n-inom{2n}{n}a^{n+1} ]

• 当(a e 1)时，两边同时除以(1-a)

  [ga_n=frac{ga_{n-1}+inom{2n-1}{n}a^n-inom{2n}{n}a^{n+1}}{1-a} ]

• 当(a= 1)时，直接把a=1代入移项后的式子里面

  [(1-a)ga_n=ga_{n-1}+inom{2n-1}{n}a^n-inom{2n}{n}a^{n+1}\ Longrightarrow ga_{n-1}=inom{2n}{n}-inom{2n-1}{n}=inom{2n-1}{n-1}\ Longrightarrow ga_{n}=inom{2n+1}{n} ]

  就可以直接算(ga_n)了。

于是，

[f_n=f_{n-1}+a^ngb_n+b^nga_n-inom{2n}{n}a^nb^n ]

综合上面的两种递推式，写三个递推就可以搞定本题了。时间复杂度(O(n))。

co int N=200000+1;
int	n,a,b,c,ans;
int fac[N<<1],ifac[N<<1],pa[N],pb[N];
int f[N],ga[N],gb[N];

il int binom(int n,int m){
	return mul(fac[n],mul(ifac[m],ifac[n-m]));
}

int main(){
	read(n),read(a),read(b),read(c);
	fac[0]=pa[0]=pb[0]=1;
	for(int i=1;i<=n;++i){
		pa[i]=mul(pa[i-1],a),pb[i]=mul(pb[i-1],b);
		fac[i]=mul(fac[i-1],i);
	}
	for(int i=n+1;i<=n<<1;++i)
		fac[i]=mul(fac[i-1],i);
	ifac[n<<1]=fpow(fac[n<<1],mod-2);
	for(int i=(n<<1)-1;i>=0;--i)
		ifac[i]=mul(ifac[i+1],i+1);
	read<int>();
	for(int i=2;i<=n;++i)
		ans=add(ans,mul(mul(read<int>(),binom((n<<1)-i-2,n-i)),mul(pa[n-1],pb[n-i])));
	read<int>();
	for(int i=2;i<=n;++i)
		ans=add(ans,mul(mul(read<int>(),binom((n<<1)-i-2,n-i)),mul(pa[n-i],pb[n-1])));
	f[0]=ga[0]=gb[0]=1;
	int ia=fpow(add(1,mod-a),mod-2),ib=fpow(add(1,mod-b),mod-2);
	for(int i=1;i<=n-2;++i){
		ga[i]=a==1?binom((i<<1)+1,i):mul(add(ga[i-1],add(mul(binom((i<<1)-1,i),pa[i]),mod-mul(binom(i<<1,i),pa[i+1]))),ia);
		gb[i]=b==1?binom((i<<1)+1,i):mul(add(gb[i-1],add(mul(binom((i<<1)-1,i),pb[i]),mod-mul(binom(i<<1,i),pb[i+1]))),ib);
		f[i]=add(f[i-1],add(mul(pa[i],gb[i]),add(mul(pb[i],ga[i]),mod-mul(binom(i<<1,i),mul(pa[i],pb[i])))));
	}
	printf("%d
",add(ans,mul(c,f[n-2])));
	return 0;
}