Codechef April Challenge 2020

题目链接

(f[u]) 表示从 (u) 走到 (n) 的生成函数。

[f[u] = frac{sum_{vin e_u} f[v]*x}{1-L_ix} ]

显然最后 (f[1]) 可以写成 (frac {A}{prod 1-L_ix}) 的形式。

于是就是一个线性递推,最开始的点值就是暴力 (dp) 的答案。

现在考虑求 (A^t mod P,P = prod (x-L_i))

如果 (L) 不相等,那么显然可以代入 (L_i) ,此时 (L_i^t) 就是取模后的点值,然后可以发现答案一定是 (sum c_i L_i^t) 的形式。

如果 (L) 有相等,考虑算出来 (H(x)=A^t mod (x-k)^q) 之后,可以 ( ext{CRT}) 合并。

(H(x+k) = (x+k)^tmod x^q)。于是我们算出 (H(x+k)) 的每一个位置的贡献之后就能轻松算出答案。

那么只需要算出来 (M'* ext{inv}(M',(x-k)^q)) 就能算出贡献了,关键是后半部分怎么求。

即给 (F(x)) 求一个 (G(x)),满足:

[F(x)*G(x) = H(x) P(x)+1,P(x) = (x-k)^q ]

同理把 (x) 变换成 (x+k) 之后就是普通求逆。

注意一下卷积的均摊,这题就解决了。

复杂度 (O(nm + n^2 + Qnlog n))

代码

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod = 1e9+7;
inline int add(int a,int b){a+=b;return a>=mod?a-mod:a;}
inline int sub(int a,int b){a-=b;return a<0?a+mod:a;}
inline int mul(int a,int b){return 1ll*a*b%mod;}
inline int qpow(int a,ll b){int ret=1;for(;b;b>>=1,a=mul(a,a))if(b&1)ret=mul(ret,a);return ret;}
/* math */

const int N = 4010;
vector<int> e[N];

/* polynomial */
typedef vector<int> poly;
poly mult(const poly &a,const poly &b){
	poly c(a.size()+b.size()-1,0);
	for(size_t i=0;i<a.size();i++)for(size_t j=0;j<b.size();j++){c[i+j]=add(c[i+j],mul(a[i],b[j]));}
	return c;
}
poly inv(const poly &b,int len){
	poly c(len,0);int inv=qpow(b[0],mod-2);c[0]=1;
	for(int i=0;i<len;i++){
		c[i]=mul(c[i],inv);for(size_t j=1;j<b.size()&&i+(int)j<len;j++)c[i+j] = sub(c[i+j], mul(c[i],b[j]));
	}
	return c;
}
poly inv(const poly &b){return inv(b,b.size());}
poly div(poly a,poly b){
	int inv=qpow(b[0],mod-2);
	for(size_t i=0;i<a.size();i++){
		a[i]=mul(a[i],inv);
		for(size_t j=1;j<b.size();j++){
			a[i+j]=sub(a[i+j],mul(a[i],b[j]));
		}
	}
	a.resize(a.size()-(b.size()-1));
	return a;
}
poly Mod(poly a,poly b){
	if(a.size()<b.size())return a;
	reverse(b.begin(),b.end());
	int inv=qpow(b[0],mod-2);
	for(int i=(int)a.size()-1;i>=(int)b.size()-1;i--){
		if(a[i]){
			int d = mul(a[i],inv);
			for(size_t j=0;j<b.size();j++)a[i-j]=sub(a[i-j],mul(b[j],d));
		}
	}
	a.resize(b.size()-1);
	return a;
}

/* n^2 ver */
int fac[N],ifac[N];
inline void init(){
	fac[0]=ifac[0]=1;for(int i=1;i<N;i++)fac[i]=mul(fac[i-1],i);
	ifac[N-1]=qpow(fac[N-1],mod-2);for(int i=N-2;i;i--)ifac[i]=mul(ifac[i+1],i+1);
}

int n,m,q;
int f[N][N],l[N];
poly getpoly(){
	poly cur(1,1);
	for(int i=1;i<=n;i++){
		poly nw(2,1);
		nw[0]=sub(0,l[i]);
		cur=mult(cur,nw);
	}return cur;
}
poly Tran;

vector< pair<int,poly> >vec;
inline int bin(int a,int b){return mul(fac[a],mul(ifac[b],ifac[a-b]));}

poly Make(poly x,int k){//g(x) = f(x+k)
	poly ret(x.size(),0);
	for(size_t i=0;i<x.size();i++){
		int t = 1;for(size_t j=0;j<=i;j++,t=mul(t,k))ret[i-j]=add(ret[i-j],mul(x[i],mul(bin(i,j),t)));
	}return ret;
}
int pw[N];
poly Make(poly x,int k,int lenth){//g(x) = f(x+k)
	poly ret(lenth,0);
	pw[0]=1;for(size_t i=1;i<x.size();i++)pw[i]=mul(pw[i-1],k);
	for(size_t i=0;i<x.size();i++){
		for(size_t j=0;j<=i&&j<(int)lenth;j++)ret[j]=add(ret[j],mul(x[i],mul(bin(i,j),pw[i-j])));
	}return ret;
}
poly iMake(poly x,int k){//g(x) = f(x+k) [inverse]
	poly ret(x.size(),0);
	for(size_t i=0;i<x.size();i++){
		int t = 1;for(size_t j=0;j<=i;j++,t=mul(t,k))ret[i]=add(ret[i],mul(x[i-j],mul(bin(i,j),t)));
	}return ret;
}
void output(poly t){
	for(size_t i=0;i<t.size();i++)cout << t[i] << " ";puts("");
}
inline void CRT(poly M1,poly M,int c,int cnt){
	poly q = Make(M1,c,cnt);
	q=inv(q);q=Make(q,mod-c);
	poly cur = Mod(mult(M1,q),Tran);
	poly tp(2,0);tp[1]=1;
	poly g(cnt,0);
	for(int i=0;i<cnt;i++){
		for(size_t j=0;j<cur.size();j++){
			g[i] = add(g[i], mul(f[1][j],cur[j]));
		}
		cur = Mod(mult(cur,tp),Tran);
	}
	vec.push_back(make_pair(c,iMake(g, mod-c)));
}

void Init(){
	Tran=getpoly();
	sort(l+1,l+n+1);
	poly trans2,trans3;
	for(int u=1,v;u<=n;u=v+1){
		poly p(2,1);
		p[0]=sub(0,l[u]);
		trans3=p;
		trans2=div(Tran,p);
		v=u;while(v+1<=n&&l[v+1]==l[v])++v,trans2=div(trans2,p),trans3=mult(trans3,p);
		CRT(trans2,trans3,l[u],v-u+1);
	}
}

int main()
{
	init();
	cin >> n >> m >> q;
	for(int i=1;i<=n;i++)scanf("%d",&l[i]);
	for(int i=1;i<=m;i++){
		int u,v;scanf("%d%d",&u,&v);e[u].push_back(v);
	}
	f[n][0]=1;
	for(int i=n;i;i--){
		for(size_t ee=0;ee<e[i].size();ee++){
			int v=e[i][ee];
			for(int j=1;j<n;j++)f[i][j]=add(f[i][j],f[v][j-1]);
		}
		for(int j=1;j<n;j++)f[i][j]=add(f[i][j],mul(f[i][j-1],l[i]));
	}
	Init();
//	Tran=getpoly();
//	poly ans;
//	for(int i=0;i<n;i++)ans.push_back(f[1][i]);
//
//	for(int i=n;i<=50000;i++){
//		int t=0;
//		for(int j=1;j<=n;j++){
//			t=sub(t,mul(ans[i-j],Tran[n-j]));
//		}
//		ans.push_back(t);
//	}
    /* bruteforce */
	while(q--){
		ll t;scanf("%lld",&t);
		if(t<n){printf("%d
",f[1][t]);continue;}
		int ret=0;
		for(size_t i=0;i<vec.size();i++){
			int pw = qpow(vec[i].first,t);
			int inv = qpow(vec[i].first,mod-2);
			int Fc=1;
			for(size_t j=0;j<vec[i].second.size();j++){
				if(j)Fc=mul(Fc,(t-j+1)%mod);
				int bin=mul(Fc,ifac[j]);
				ret=add(ret,mul(mul(bin,vec[i].second[j]),pw));
				pw=mul(pw,inv);
			}
		}
		printf("%d
",ret);
	}
}
/* 
5 5 4
6 6 1 7 6
1 2
1 3
2 4
3 4
4 5
1000 2000 3000 40000
 */

原文地址:https://www.cnblogs.com/weiyanpeng/p/12694061.html