超矩形
求体积不超过(n)的所有(k)维长方体在每条边分别取(b_i(0leq b_ileq 9))次方后的体积和。序列(b)可以被分成至多五段相同的数字。
(Tleq 100,nleq 10^6,sum kleq 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,lfloorfrac{n}{j}
floor)
]
于是可以(O(Tsum 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(Tnlog nlog k))。
注意到不同质数的幂对答案的贡献相互不影响，所以答案求的是一个关于体积的积性函数(f(x))的前缀和。
考虑用Min_25筛来求出(f(x))的前缀和。
只需要搞定(p)及(p^k)处的取值。
(p)只可能在一种(b)中出现，对于某一种(b)，(f_b(p)= ext{cnt}_bp^b)。
但是(p^k)可以在多种(b)中出现，我们先求出一种的，然后将多种的卷积起来（指数做加法），(f_b(p^k)=inom{k+ ext{cnt}_b-1}{ ext{cnt}_b-1}(p^k)^b)。
预处理时间复杂度(O(frac{sqrt{n}}{ln{sqrt{n}}}log_p^2 n))，可以忽略不计。
所以我们可以在(O(Tfrac{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
",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;
}