码迷,mamicode.com
首页 > 其他好文 > 详细

SCOI2019 超矩形

时间:2020-06-19 12:12:52      阅读:53      评论:0      收藏:0      [点我收藏+]

标签:www   return   cto   oid   一个   ima   pow   sig   次方   

超矩形

求体积不超过\(n\)的所有\(k\)维长方体在每条边分别取\(b_i(0\leq b_i\leq 9)\)次方后的体积和。序列\(b\)可以被分成至多五段相同的数字。

\(T\leq 100,n\leq 10^6,\sum k\leq 10^6\)

题解

https://blog.csdn.net/sslz_fsy/article/details/104589438

首先可以暴力DP,\(f(i,n)\)表示前\(i\)维体积为\(n\)的和

\[f(i,n)=\sum_{j|n}j^{b_i}f(i-1,\frac{n}{j}) \]

写成前缀和的形式

\[F(i,n)=\sum_{j=1}^nj^{b_i}F(i-1,\lfloor\frac{n}{j}\rfloor) \]

于是可以\(O(T\sum kn^{\frac{3}{4}})\)暴力整除分块做。

考虑生成函数

\[f(x)=\prod_{i=1}^kg_{b_i}(x) \]

\[g_{b_i}(x)=\sum_{i=1}^ni^{b_i}x^i \]

乘法是狄利克雷卷积,所以可以对每一种\(b\)倍增然后合并,\(O(Tn\log n\log k)\)

注意到不同质数的幂对答案的贡献相互不影响,所以答案求的是一个关于体积的积性函数\(f(x)\)的前缀和。

考虑用Min_25筛来求出\(f(x)\)的前缀和。

只需要搞定\(p\)\(p^k\)处的取值。

\(p\)只可能在一种\(b\)中出现,对于某一种\(b\)\(f_b(p)=\text{cnt}_bp^b\)

但是\(p^k\)可以在多种\(b\)中出现,我们先求出一种的,然后将多种的卷积起来(指数做加法),\(f_b(p^k)=\binom{k+\text{cnt}_b-1}{\text{cnt}_b-1}(p^k)^b\)

预处理时间复杂度\(O(\frac{\sqrt{n}}{\ln{\sqrt{n}}}\log_p^2 n)\),可以忽略不计。

所以我们可以在\(O(T\frac{n^{\frac{3}{4}}}{\log n})\)的时间内解决这个题。

https://www.cnblogs.com/ImagineC/p/10682026.html

vector<int> operator*(CO vector<int>&a,CO vector<int>&b){
	vector<int> ans(a.size());
	for(int i=0;i<(int)a.size();++i)for(int j=0;j+i<(int)a.size();++j)
		ans[i+j]=add(ans[i+j],mul(a[i],b[j]));
	return ans;
}

CO int N=1e6+30,M=1e3+10;
int fac[N],ifac[N],pow_sum[10][N];
int n,K,q,m,ex[5],cnt[5];
int prm[M],num,prm_sum[5][M];
vector<int> poly[M];
int pos[2*M],idx,ref1[M],ref2[M],H[5][2*M];

IN int C(int n,int m){
	return mul(fac[n],mul(ifac[m],ifac[n-m]));
}
int G(int n,int x){
	if(prm[x]>=n) return 0;
	int ans=0,r=n<=m?ref1[n]:ref2[::n/n];
	for(int i=0;i<q;++i)
		ans=add(ans,mul(cnt[i],H[i][r]+mod-prm_sum[i][x]));
	for(int i=x+1;i<=num and prm[i]*prm[i]<=n;++i)
		for(int j=1,y=prm[i];y<=n;++j,y*=prm[i])
			ans=add(ans,mul(poly[i][j],(j>1)+G(n/y,i)));
	return ans;
}
void real_main(){
	read(n),read(K),read(q),m=sqrt(n);
	for(int i=0;i<q;++i) read(ex[i]),read(cnt[i]);
	fill(prm,prm+m+1,0);
	for(int i=2;i<=m;++i){
		if(!prm[i]){
			prm[++num]=i;
			for(int j=0;j<q;++j){
				int v=add(pow_sum[ex[j]][i],mod-pow_sum[ex[j]][i-1]);
				prm_sum[j][num]=add(prm_sum[j][num-1],v);
			}
		}
		for(int j=1;j<=num and i*prm[j]<=m;++j){
			prm[i*prm[j]]=1;
			if(i%prm[j]==0) break;
		}
	}
	for(int i=1;i<=n;i=n/(n/i)+1)
		pos[++idx]=n/i,n/i<=m?ref1[n/i]=idx:ref2[n/(n/i)]=idx;
	for(int i=0;i<q;++i)for(int j=1;j<=idx;++j)
		H[i][j]=pow_sum[ex[i]][pos[j]]-1;
	for(int i=1;i<=num;++i)for(int j=1;j<=idx and prm[i]*prm[i]<=pos[j];++j){
		int x=pos[j]/prm[i],r=x<=m?ref1[x]:ref2[n/x];
		for(int k=0;k<q;++k){
			int v=add(pow_sum[ex[k]][prm[i]],mod-pow_sum[ex[k]][prm[i]-1]);
			H[k][j]=add(H[k][j],mod-mul(v,H[k][r]+mod-prm_sum[k][i-1]));
		}
	}
	for(int i=1;i<=num;++i){
		int up=0;
		for(int x=1;x*prm[i]<=n;x*=prm[i]) ++up;
		poly[i].assign(up+1,0),poly[i][0]=1;
		for(int j=0;j<q;++j){
			vector<int> tmp(up+1);
			for(int k=0;k<=up;++k)
				tmp[k]=mul(C(k+cnt[j]-1,cnt[j]-1),fpow(prm[i],k*ex[j]));
			poly[i]=poly[i]*tmp;
		}
	}
	printf("%d\n",G(n,0)+1);
}
int main(){
	freopen("hyper.in","r",stdin),freopen("hyper.out","w",stdout);
	fac[0]=1;
	for(int i=1;i<N;++i) fac[i]=mul(fac[i-1],i);
	ifac[N-1]=fpow(fac[N-1],mod-2);
	for(int i=N-2;i>=0;--i) ifac[i]=mul(ifac[i+1],i+1);
	for(int k=0;k<=9;++k)for(int i=1;i<N;++i)
		pow_sum[k][i]=add(pow_sum[k][i-1],fpow(i,k));
	for(int T=read<int>();T--;) real_main();
	return 0;
}

SCOI2019 超矩形

标签:www   return   cto   oid   一个   ima   pow   sig   次方   

原文地址:https://www.cnblogs.com/autoint/p/13162102.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!