给定整数区间 \([1,m]\) ,从中可重复的选择 \(n\) 个数,形成一个数列 \(\{a_n\}\) 。问:所有满足 \(\gcd(a_1,...,a_n)\le q\) 并且 \(\text{lcm}(a_1,...,a_n)\ge p\) 的数列的乘积和。
官方题解其实已经很明了了,我这里再做个翻译。题目要求的是 \(\gcd\le q\) 且 \(\text{lcm} \ge p\) 的所有数列的乘积和。根据 \(\gcd|\text{lcm}\) 这一性质可以枚举 \((\gcd,\text{lcm})\) 数对来计算,但此题的 \(\text{lcm}\) 会很大,无法直接枚举所有大于 \(p\) 的 \(\text{lcm}\) ,所以需要转换思路,大的枚举不完,但小的好枚举,于是我们做这么一个容斥:
\[(\gcd\le q,\text{lcm}\ge p)=S-(\gcd>q)-(\gcd\le q,\text{lcm}<p) \]其中 \((a,b)\) 表示满足 \(a\) 和 \(b\) 条件下的所有数列的乘积和, \(S\) 代表全体数列的乘积和。
全体数列的乘积和可以这么表示:每一个数的取值范围都是 \([1,m]\) ,利用分配律,可知
\[S=(1+2+...+m)(1+2+...+m)...(1+2+...+m)=(1+2+...+m)^n \]然后再求 \((\gcd>q)\) 。这是个经典问题,很容易由容斥算出来。
然后就是 \((\gcd\le q,\text{lcm}<p)\) 。考虑上述的分配律,我们也可以把满足 \(\gcd=g,\text{lcm}=l\) 的所有数列的乘积和表示成一些数乘积的和。考虑到唯一分解定理,我们使用素数来进行表示。对于两个数 \(a,b\) ,把他们写成唯一分解的形式:
\[a=p_1^{k_1}p_2^{k_2}...p_s^{k_s}\\ b=p_1^{t_1}p_2^{t_2}...p_s^{t_s}\\ \gcd(a,b)=p_1^{\min(k_1,t_1)}p_2^{\min(k_2,t_2)}...p_s^{\min(k_s,t_s)}\\ \text{lcm}(a,b)=p_1^{\max(k_1,t_1)}p_2^{\max(k_2,t_2)}...p_s^{\max(k_s,t_s)} \]由此可知,一对 \(\gcd,\text{lcm}\) 便确定了素数幂次的上界和下界,利用上述的分配律计算贡献即可。
这份代码没有优化,1887ms,过的非常危险。
#include <bits/stdc++.h> #define rep(i,a,b) for(int i=a;i<b;i++) using namespace std; using ll = long long; constexpr int N = 2e5+5, P=998244353; ll f[N], mind[N]; ll qpow(ll a, ll b){ ll res=1; for(;b;b>>=1,a=a*a%P)if(b&1)res=res*a%P; return res; } int main() { #ifndef ONLINE_JUDGE freopen("1.in", "r", stdin); #endif // !ONLINE_JUDGE cin.tie(nullptr)->ios::sync_with_stdio(false); rep(i,2,N)if(!mind[i])for(int j=i;j<N;j+=i)mind[j]=i; ll n,m,p,q; cin>>n>>m>>p>>q; ll inv2=qpow(2,P-2); ll ans=qpow((m*(m+1)%P*inv2)%P,n); for(int i=m;i>q;i--){ ll cnt=m/i; f[i]=qpow(i,n)*qpow((cnt*(cnt+1)%P*inv2%P),n)%P; for(int j=i+i;j<=m;j+=i)f[i]=(f[i]-f[j]+P)%P; ans=(ans-f[i]+P)%P; }//cout<<ans<<endl; rep(l,1,m+1){ f[l]=1; for(int tmp=l,pr;tmp!=1;){ pr=mind[tmp]; ll x=0,y=1; for(;tmp%pr==0;x=(x+y)%P,y=y*pr%P)tmp/=pr; ll temp=(qpow(x+y,n)-qpow(x,n)-qpow(x+y-1,n)+qpow(x-1,n))%P; f[l]=(f[l]*temp)%P; } } rep(i,1,q+1)for(int j=i;j<p;j+=i)ans=(ans-qpow(i,n)*f[j/i]%P+P)%P; cout<<ans; return 0; }