标签：frac limits loj sum varphi 反演 int rg 6179
题目描述
求 \(\sum\limits_{i=1}^n\sum\limits_{j=1}^m \varphi(ij)(mod\ 998244353)\)
\(T\) 组询问
\(1 \leq n,m,T \leq 10^5\)
分析
令 \(n<m\)
首先，我们把 \(\varphi(ij)\) 拆成 \(\varphi(i)\varphi(j)\frac{gcd(i,j)}{\varphi(gcd(i,j))}\)
考虑求单个欧拉函数的方法
\(\varphi(i)=i\prod\limits_{d \in prime,d|i}\frac{d-1}{d}\)
\(\varphi(i)\varphi(j)\) 比 \(\varphi(ij)\) 多乘了 \(i,j\) 共有的质因子
所以我们要把这些共有的质因子做出的贡献消除
即乘上一个 \(\frac{gcd(i,j)}{\varphi(gcd(i,j))}\)
那么式子就变成了
\(\sum\limits_{i=1}^n\sum\limits_{j=1}^m\varphi(i)\varphi(j)\frac{gcd(i,j)}{\varphi(gcd(i,j))}\)
按照莫比乌斯反演的套路，枚举 \(gcd(i,j)\)
\(\sum\limits_{d=1}^n\frac{d}{\varphi(d)}\sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{m/d}\varphi(id)\varphi(jd)[gcd(i,j)=1]\)
\(\sum\limits_{d=1}^n\frac{d}{\varphi(d)}\sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{m/d}\varphi(id)\varphi(jd)\sum\limits_{k|gcd(i,j)}\mu(k)\)
\(\sum\limits_{d=1}^n\frac{d}{\varphi(d)}\sum\limits_{k=1}^{n/d}\mu(k)\sum\limits_{i=1}^{n/dk}\sum\limits_{j=1}^{m/dk}\varphi(idk)\varphi(jdk)\)
\(\sum\limits_{T=1}^n\sum\limits_{d|T}\frac{d}{\varphi(d) }\mu(\frac{T}{d})\sum\limits_{i=1}^{n/T}\sum\limits_{j=1}^{m/T}\varphi(iT)\varphi(jT)\)
设 \(g(T)=\sum\limits_{d|T}\frac{d}{\varphi(d) }\mu(\frac{T}{d}),s(n,m)=\sum\limits_{i=1}^n\varphi(mi)\)
原式变成
\(\sum\limits_{T=1}^ng(T)s(n/T,T)s(n/T,T)\)
\(g\) 数组和 \(s\) 数组可以 \(nlogn\) 预处理出来
然后就可以 \(O(n)\) 计算每一次的答案
总复杂度就是 \(O(nT)\)，还是不够优秀
我们可以人为地设一个值 \(top\)
当 \(T<\frac{m}{top}\) 时暴力计算
否则 \(\frac{n}{T},\frac{m}{T}\) 一定是小于 \(top\) 的
那么我们就可以开一个数组把小于 \(top\) 的这一部分预处理出来
查询的时候就可以直接用
这道题的内存限制是 \(64M\)
所以不能直接开数组
而要开一个内存池或者用 \(vector\)
代码
#include<cstdio>
#include<cmath>
#include<vector>
#define rg register
const int maxn=1e5+5,maxm=55,mod=998244353,tp=52;
inline int addmod(rg int now1,rg int now2){
return now1+=now2,now1>=mod?now1-mod:now1;
}
inline int delmod(rg int now1,rg int now2){
return now1-=now2,now1<0?now1+mod:now1;
}
inline int mulmod(rg long long now1,rg int now2){
return now1*=now2,now1>=mod?now1%mod:now1;
}
int pri[maxn],mmax,phi[maxn],mu[maxn],g[maxn],*f[maxm][maxm],buf[maxn*160],*o=buf,*s[maxn];
int ksm(rg int ds,rg int zs){
rg int nans=1;
while(zs){
if(zs&1) nans=mulmod(nans,ds);
ds=mulmod(ds,ds);
zs>>=1;
}
return nans;
}
bool not_pri[maxn];
void pre(){
not_pri[0]=not_pri[1]=1;
phi[1]=mu[1]=1;
for(rg int i=2;i<=mmax;i++){
if(!not_pri[i]){
pri[++pri[0]]=i;
phi[i]=i-1;
mu[i]=mod-1;
}
for(rg int j=1;j<=pri[0] && 1LL*i*pri[j]<=mmax;j++){
not_pri[i*pri[j]]=1;
if(i%pri[j]==0){
phi[i*pri[j]]=mulmod(pri[j],phi[i]);
break;
} else {
phi[i*pri[j]]=mulmod(phi[pri[j]],phi[i]);
mu[i*pri[j]]=mulmod(mu[pri[j]],mu[i]);
}
}
}
rg int cs;
for(rg int i=1;i<=mmax;i++){
cs=mulmod(i,ksm(phi[i],mod-2));
for(rg int j=i,now=1;j<=mmax;j+=i,now++){
g[j]=addmod(g[j],mulmod(cs,mu[now]));
}
}
for(rg int i=1;i<=mmax;i++){
cs=mmax/i;
s[i]=o;
o+=(cs+1);
for(rg int j=1;j<=cs;j++){
s[i][j]=addmod(s[i][j-1],phi[i*j]);
}
}
for(rg int i=1;i<=tp;i++){
for(rg int j=1;j<=tp;j++){
f[i][j]=o;
f[i][j][0]=0;
cs=std::min(mmax/i,mmax/j);
o+=(cs+1);
for(rg int k=1;k<=cs;k++){
f[i][j][k]=(addmod(f[i][j][k-1],mulmod(g[k],mulmod(s[k][i],s[k][j]))));
}
}
}
}
int t,n,m;
int main(){
mmax=1e5;
pre();
scanf("%d",&t);
while(t--){
scanf("%d%d",&n,&m);
if(n>m) std::swap(n,m);
rg int nans=0,mmax=m/tp,orz1,orz2;
for(rg int i=1;i<=mmax;i++){
nans=addmod(nans,mulmod(g[i],mulmod(s[i][n/i],s[i][m/i])));
}
for(rg int l=mmax+1,r;l<=n;l=r+1){
r=std::min((n/(n/l)),m/(m/l)),orz1=n/l,orz2=m/l;
nans=addmod(nans,delmod(f[orz1][orz2][r],f[orz1][orz2][l-1]));
}
printf("%d\n",nans);
}
return 0;
}
标签：frac,limits,loj,sum,varphi,反演,int,rg,6179 来源： https://www.cnblogs.com/liuchanglc/p/14261292.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。