1.树林里的树. Trees in a Wood
题目大意
求\(\frac{\sum_{i=-a}^{a}\sum_{j=-b}^{b}\lbrack\lvert gcd(i,j)\rvert=1\rbrack}{(2a+1)(2b+1)-1}(a\le 2000，b\le 2\times 10^6)\)，结果保留7位小数
分析
一道非常基础的莫反题
可以发现4个象限对答案的贡献是一样的，考虑只处理第一象限，最后将答案乘4，再加上坐标轴上的4个点，注意总点数是(2a+1)(2b+1)-1（原点不算在内）
第一象限的贡献为
\[s=\sum_{i=1}^a\sum_{j=1}^b[gcd(i,j)=1] \]可以发现是一个典型莫反的式子，按照套路
\[s=\sum_{i=1}^a\sum_{j=1}^b\sum_{d|gcd(i,j)}\mu(d) \]考虑枚举d
\[s=\sum_{d=1}^a\mu(d)\cdot\sum_{i=1}^{\lfloor\frac{a}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{b}{d}\rfloor}1 \]不难发现
\[\sum_{i=1}^{\lfloor\frac{a}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{b}{d}\rfloor}1=\lfloor\frac{a}{d}\rfloor\cdot\lfloor\frac{b}{d}\rfloor \]于是
\[s=\sum_{d=1}^a\mu(d)\cdot \lfloor\frac{a}{d}\rfloor\cdot\lfloor\frac{b}{d}\rfloor\]当\(l\le d\le r\)时，\(\lfloor\frac{a}{d}\rfloor\cdot\lfloor\frac{b}{d}\rfloor\)的值不变，考虑将\(\sum_{i=l}^r\mu(i)\)一起计算
令\(sum[d]=\sum_{i=1}^d\mu(i)\)，则\(\sum_{d=l}^r\mu(d)=sum[r]-sum[l-1]\)可以\(O[1]\)得到
考虑数论分块，因为\(\lfloor\frac{n}{d}\rfloor\)只有\(\lfloor\sqrt{n}\rfloor\)种取值，所以在预处理出\(sum[i]\)的情况下，可以\(O(\sqrt{n})\)求出\(s\)
接下来考虑如何预处理\(sum\)
由\(\mu\)函数的定义可知其为积性函数，可以用欧拉筛\(O(n)\)筛出（杜教筛更优，但是我不会），那么就可以\(O(n)\)预处理出\(\mu\)的前缀和\(sum\)
总复杂度为\(O(n+Q\cdot\sqrt{n})\)，可以通过此题~
AC·code
#include<bits/stdc++.h>
#define il inline
#define cs const
#define ri register
#define int long long
using namespace std;
namespace Q{
il int rd(){
ri int x=0,f=1;ri char c=getchar();
while(c<'0'||c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0'&&c<='9') x=x*10+(c^48),c=getchar();
return x*f;
}
cs int N=2000;
int mu[N+5],ss[N+5],cnt,sum[N+5];
bool vs[N+5];
il void init(){
mu[1]=1,sum[1]=1;
for(ri int i=2;i<=N;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1;
for(ri int j=1;j<=cnt&&ss[j]*i<=N;++j){
vs[i*ss[j]]=1;
if(i%ss[j]==0) break;
mu[i*ss[j]]=-mu[i];
}
sum[i]=sum[i-1]+mu[i];
}
}
il void done(int n,int m){
if(n>m) swap(n,m);
double as=0;
for(ri int l=1,r=0;l<=n;l=r+1){
r=min(n/(n/l),m/(m/l));
as+=(sum[r]-sum[l-1])*(n/l)*(m/l);
}
as=as*4+4;
as=as/((n*2+1)*(m*2+1)-1);
printf("%.7lf\n",as);
return;
}
}
using namespace Q;
signed main(){
init();
int a=rd(),b=rd();
while(a||b){
done(a,b);
a=rd(),b=rd();
}
return 0;
}
2.P2257 YY的GCD
题目大意
求\(\sum_{i=1}^n\sum_{j=1}^m\lbrack gcd(i,j)\in prime\rbrack (n,m\le 10^7)\)
分析
考虑枚举质数
\[s=\sum_{p\le n,p\in prime}\sum_{i=1}^{\lfloor\frac{n}{p}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{p}\rfloor}\lbrack gcd(i,j)=1\rbrack \]这样就跟上一题差不多了
\[s=\sum_{p\le n,p\in prime}\sum_{i=1}^{\lfloor\frac{n}{p}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{p}\rfloor}\sum_{d|gcd(i,j)}\mu(d) \]枚举\(d\)
\[\begin{align*} s &=\sum_{p\le n,p\in prime}\sum_{d=1}^{\lfloor\frac{n}{p}\rfloor}\mu(d)\cdot\sum_{i=1}^{\lfloor\frac{n}{p\cdot d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{p\cdot d}\rfloor}1 \\ &=\sum_{p\le n,p\in prime}\sum_{d=1}^{\lfloor\frac{n}{p}\rfloor}\mu(d)\cdot\lfloor\frac{n}{p\cdot d}\rfloor\cdot\lfloor\frac{m}{p\cdot d}\rfloor \end{align*}\]一个常用套路，令\(T=p\cdot d\)，枚举\(T\)
\[s=\sum_{T=1}^{n}\lfloor\frac{n}{T}\rfloor\cdot\lfloor\frac{m}{T}\rfloor\cdot\sum_{p|T,p\in prime}\mu(\frac{T}{p}) \]发现\(\sum_{p|T,p\in prime}\mu(\frac{T}{p})\)可以预处理
令\(f(T)=\sum_{p|T,p\in prime}\mu(\frac{T}{p})\)对于每一个质数\(p\)，枚举\(d\)，将\(f(d\cdot p)\)加上\(\mu(d)\)，第二维循环是一个调和级数，预处理总复杂度约为\(n\cdot\ln n\)
预处理出\(f(T)\)的前缀和，然后就可以用数论分块\(O(\sqrt n)\)求出\(s\)
总复杂度约为\(O(n\cdot\ln n+Q\cdot\sqrt n)\)，可过~
\(ps：f(T)\)可以\(O(n)\)线性筛出来
\[ f(i\cdot p)=\begin{cases} \mu(1) & i=1\\ \mu(i) &i\mod p=0\\ -f(i)+\mu(i) &i\mod p\ne 0 \end{cases}\]线性筛f(T)
il void init(int n){
mu[1]=1;
for(ri int i=2;i<=n;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1,f[i]=1;
for(ri int j=1;j<=tot&&ss[j]*i<=n;++j){
vs[i*ss[j]]=1;
if(i%ss[j]){
f[i*ss[j]]=-f[i]+mu[i];
mu[i*ss[j]]=-mu[i];
}
else{
f[i*ss[j]]=mu[i];
break;
}
}
sum[i]=sum[i-1]+f[i];
}
return;
}
AC·code
#include<bits/stdc++.h>
#define il inline
#define cs const
#define ri register
using namespace std;
il int rd(int &x){
x=0;ri int f=1;ri char c=getchar();
while(c<'0'||c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0'&&c<='9') x=x*10+(c^48),c=getchar();
return x*f;
}
il void wt(long long x){
if(x<0) x=-x,putchar('-');
ri unsigned short a[15]={0},tl=0;
do{a[++tl]=x%10,x/=10;}while(x);
for(ri int i=tl;i;--i) putchar(a[i]+48);
}
cs int N=1e7;
long long sum[N+5];
int ss[N+5],mu[N+5],cnt,f[N+5];
bool vs[N+5];
il void init(int n){
mu[1]=1;
for(ri unsigned int i=2;i<=n;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1;
for(ri unsigned int j=1;j<=cnt&&i*ss[j]<=n;++j){
vs[i*ss[j]]=1;
if(i%ss[j]==0) break;
mu[i*ss[j]]=-mu[i];
}
}
for(ri unsigned int j=1;j<=cnt;++j){
for(ri int i=1;ss[j]*i<=n;++i){
f[i*ss[j]]+=mu[i];
}
}
for(ri unsigned int i=1;i<=n;++i){
sum[i]=sum[i-1]+1ll*f[i];
}
}
signed main(){
int t=0,n=0,m=0;rd(t);
ri long long as=0;init(N);
while(t--){
rd(n),rd(m),as=0;
if(m<n) swap(m,n);
for(ri unsigned int l=1,r=0;l<=n;l=r+1){
r=min((n/l)?n/(n/l):n,(m/l)?m/(m/l):m);
as+=1ll*(n/l)*(m/l)*(sum[r]-sum[l-1]);
}
wt(as),putchar('\n');
}
return 0;
}
3.P2398 GCD SUM
题目大意
求\(\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)\)
分析
莫反做法
这题比上一题简单一些（那你为什么先写上一题）
考虑枚举\(gcd\)
\[\begin{align*} s &=\sum_{i=1}^n\sum_{j=1}^n\sum_{d=1}^n d\cdot\lbrack\gcd(i,j)=d\rbrack \\&=\sum_{d=1}^n d\cdot\sum_{i=1}^n\sum_{j=1}^n[gcd(i,j)=d] \end{align*}\]将后面两项同时除以\(d\)，可得
\[s=\sum_{d=1}^n d\cdot\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}[gcd(i,j)=1] \]一个标准的莫反式子
\[s=\sum_{d=1}^n d\cdot\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{g|gcd(i,j)}\mu(g) \]枚举\(g\)
\[\begin{align*} s &=\sum_{d=1}^n d\cdot\sum_{g=1}^{\lfloor\frac{n}{d}\rfloor}\mu(g)\cdot\sum_{i=1}^{\lfloor\frac{n}{d\cdot g}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d\cdot g}\rfloor}1 \\&=\sum_{d=1}^n d\cdot\sum_{g=1}^{\lfloor\frac{n}{d}\rfloor}\mu(g)\cdot\lfloor\frac{n}{d\cdot g}\rfloor^2 \end{align*}\]对于
\[\sum_{g=1}^{\lfloor\frac{n}{d}\rfloor}\mu(g)\cdot\lfloor\frac{n}{d\cdot g}\rfloor^2 \]可以线性筛预处理出\(\mu(g)\)的前缀和，再进行数论分块，复杂度\(O(\sqrt{n})\)
总复杂度\(O(n+n\cdot\sqrt{n})\) ，在该题数据范围\((n\le10^5)\)下可过
事实上可以继续往下推
考虑枚举\(T=d\cdot g\)
\[\begin{align*} s&=\sum_{T=1}^n\sum_{d|T}d\cdot\mu(\frac{T}{d})\cdot\lfloor\frac{n}{T}\rfloor^2\\ &=\sum_{T=1}^n\sum_{d|T}id(d)\cdot\mu(\frac{T}{d})\cdot\lfloor\frac{n}{T}\rfloor^2 \end{align*} \]根据狄利克雷卷积\(id*\mu=\varphi\)，得
\[s=\sum_{T=1}^n\varphi(T)\cdot\lfloor\frac{n}{T}\rfloor^2 \]可以线性筛预处理\(\varphi(d)\)前缀和，然后数论分块，总复杂度\(O(n+\sqrt{n})\)
欧拉反演
可以直接用欧拉反演做，更好推~
欧拉函数有一个性质
\[\sum_{d|n}\varphi(d)=n \]所以
\[\gcd(i,j)=\sum_{d|gcd(i,j)}\varphi(d) \]于是
\[\begin{align*} s &=\sum_{i=1}^n\sum_{j=1}^n\sum_{d|gcd(i,j)}\varphi(d) \\&=\sum_{d=1}^n\varphi(d)\cdot\sum_{i=1}^n[i|d]\cdot\sum_{j=1}^n[j|d] \\&=\sum_{d=1}^n\varphi(d)\cdot\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}1 \\&=\sum_{d=1}^n\varphi(d)\cdot\lfloor\frac{n}{d}\rfloor^2 \end{align*}\]跟上面推的式子是一样的
玄学递推
然而……还有更玄学简洁的做法……
\(f(k)\) 表示\(gcd(i,j)=k\)的\((i,j)\)对数
而\(k|gcd(i,j)\)的\((i,j)\)对数为\(\lfloor\frac{n}{k}\rfloor^2\)
所以容斥一下
\[f(k)=\lfloor\frac{n}{k}\rfloor^2-f(2\cdot k)-f(3\cdot k)-…… \]复杂度 \(n\cdot(1+\frac{1}{2}+\frac{1}{2}+……+\frac{1}{n})\) 为 \(O(n\cdot H(n))\)（\(H(n)\)为调和级数，约等于\(\ln n\)）
因其常数较小，所以比欧拉反演还快……
AC·code（莫反）
#include<bits/stdc++.h>
#define il inline
#define cs const
#define int long long
#define ri register
using namespace std;
cs int N=1e5+5;
int ss[N],vs[N],cnt,mu[N];
il void init(int n){
mu[1]=1;
for(ri int i=2;i<=n;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1;
for(ri int j=1;j<=cnt&&ss[j]*i<=n;++j){
vs[i*ss[j]]=1;
if(i%ss[j]==0) break;
mu[i*ss[j]]=-mu[i];
}
}
for(ri int i=2;i<=n;++i) mu[i]+=mu[i-1];
}
signed main(){
int n,as=0;cin>>n;init(n);
for(ri int d=1;d<=n;++d){
int m=n/d;
for(ri int l=1,r=0;l<=m;l=r+1){
r=m/(m/l);
as+=d*(mu[r]-mu[l-1])*(m/l)*(m/l);
}
}
cout<<as;
return 0;
}
AC·code（欧拉）
//题解上粘的代码
#include <cstdio>
#include <cstring>
const int MAXN = 100010;
int prime[MAXN], v[MAXN], phi[MAXN], sumPhi[MAXN], cnt, n;
void eular(int n) {//线性筛筛欧拉函数
memset(v, 0, sizeof(v));
sumPhi[1] = phi[1] = 1;
for (int i = 2; i <= n; i++) {
if (!v[i]) {
prime[++cnt] = v[i] = i;
phi[i] = i - 1;
}
for (int j = 1; j <= cnt; j++) {
if (prime[j] > v[i] || i * prime[j] > n) break;
v[i * prime[j]] = prime[j];
phi[i * prime[j]] = phi[i] * (i % prime[j] ? prime[j] - 1 : prime[j]);
}
sumPhi[i] = sumPhi[i - 1] + phi[i];//求欧拉函数的前缀和，如果整除分块的话就要用
}
}
int main() {
scanf("%d", &n);
eular(n);
long long ans = 0;
for (int l = 1, r; l <= n; l = r + 1) {//这里是整除分块写法，如果不懂可以直接for 1 to n
r = n / (n / l);
ans += (long long) (sumPhi[r] - sumPhi[l - 1]) * (n / l) * (n / l);
}
printf("%lld\n", ans);
return 0;
}
AC·code（容斥）
#include<bits/stdc++.h>
#define ri register
long long n,ans,f[100010];
int main(){
scanf("%lld",&n);
for(ri int i=n;i;--i){
f[i]=n/i*(n/i);
for(ri int j=i<<1;j<=n;j+=i)f[i]-=f[j];
ans+=f[i]*i;
}
printf("%lld",ans);
return 0;
}
4.P4450 双亲数
题目大意
求\(\sum_{a=1}^A\sum_{b=1}^B[gcd(a,b)=d]\) \((1\le A,B\le10^6,1\le d\le min(A,B))\)
分析
将式子整体除以\(d\)
\[s=\sum_{i=1}^{\lfloor\frac{A}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{B}{d}\rfloor}[gcd(i,j)=1] \]直接莫反
\[\begin{align*} s &=\sum_{i=1}^{\lfloor\frac{A}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{B}{d}\rfloor}\sum_{g|gcd(i,j)}\mu(g) \\&=\sum_{g=1}^{\lfloor\frac{A}{d}\rfloor}\mu(g)\cdot\sum_{i=1}^{\lfloor\frac{A}{d\cdot g}\rfloor}\sum_{j=1}^{\lfloor\frac{B}{d\cdot g}\rfloor}1 \\&=\sum_{g=1}^{\lfloor\frac{A}{d}\rfloor}\mu(g)\cdot\lfloor\frac{A}{d\cdot g}\rfloor\cdot\lfloor\frac{B}{d\cdot g}\rfloor \end{align*}\]线性筛加数论分块，\(O(n+\sqrt{n})\)
AC·code
#include<bits/stdc++.h>
#define il inline
#define cs const
#define ri register
using namespace std;
il int rd(){
ri int x=0,f=1;ri char c=getchar();
while(c<'0'||c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0'&&c<='9') x=x*10+(c^48),c=getchar();
return x*f;
}
il void wt(long long x){
if(x<0) x=-x,putchar('-');
ri unsigned short a[15]={0},tl=0;
do{a[++tl]=x%10,x/=10;}while(x);
for(ri int i=tl;i;--i) putchar(a[i]+48);
}
cs int N=1e6;
int ss[N+5],mu[N+5],cnt,sum[N+5];
bool vs[N+5];
il void init(int n){
mu[1]=sum[1]=1;
for(ri unsigned int i=2;i<=n;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1;
for(ri unsigned int j=1;j<=cnt&&i*ss[j]<=n;++j){
vs[i*ss[j]]=1;
if(i%ss[j]==0) break;
mu[i*ss[j]]=-mu[i];
}
sum[i]=sum[i-1]+mu[i];
}
}
signed main(){
int n=0,m=0,d=0;
ri long long as=0;init(N);
n=rd(),m=rd(),d=rd(),as=0;
n/=d,m/=d;
if(m<n) swap(m,n);
for(ri unsigned int l=1,r=0;l<=n;l=r+1){
r=min(n/(n/l),m/(m/l));
as+=1ll*(n/l)*(m/l)*(sum[r]-sum[l-1]);
}
wt(as);
return 0;
}
5.P3455 [POI2007]ZAP-Queries
题目大意
求\(\sum_{x=1}^a\sum_{y=1}^b[gcd(x,y)=d]\) \((1\le n\le5\times10^4, 1\le d\le a,b\le5\times10^4)\)
分析
跟上一题一样（双倍经验）
将式子整体除以\(d\)莫反可得
\[s=\sum_{g=1}^{\lfloor\frac{a}{d}\rfloor}\mu(g)\cdot\lfloor\frac{a}{d\cdot g}\rfloor\cdot\lfloor\frac{b}{d\cdot g}\rfloor \]线性筛加数论分块，\(O(n+Q\cdot\sqrt{n})\)
AC·code
#include<bits/stdc++.h>
#define il inline
#define cs const
#define ri register
using namespace std;
il int rd(){
ri int x=0,f=1;ri char c=getchar();
while(c<'0'||c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0'&&c<='9') x=x*10+(c^48),c=getchar();
return x*f;
}
il void wt(long long x){
if(x<0) x=-x,putchar('-');
ri unsigned short a[15]={0},tl=0;
do{a[++tl]=x%10,x/=10;}while(x);
for(ri int i=tl;i;--i) putchar(a[i]+48);
}
cs int N=5e4;
int ss[N+5],mu[N+5],cnt,sum[N+5];
bool vs[N+5];
il void init(int n){
mu[1]=sum[1]=1;
for(ri unsigned int i=2;i<=n;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1;
for(ri unsigned int j=1;j<=cnt&&i*ss[j]<=n;++j){
vs[i*ss[j]]=1;
if(i%ss[j]==0) break;
mu[i*ss[j]]=-mu[i];
}
sum[i]=sum[i-1]+mu[i];
}
}
signed main(){
int t=rd(),n=0,m=0,d=0;
ri long long as=0;init(N);
while(t--){
n=rd(),m=rd(),d=rd(),as=0;
if(m<n) swap(m,n);
for(ri unsigned int l=1,r=0;l<=n;l=r+1){
r=min(n/(n/l),m/(m/l));
as+=1ll*(n/l/d)*(m/l/d)*(sum[r]-sum[l-1]);
}
wt(as),putchar('\n');
}
return 0;
}
6.P2522 [HAOI2011]Problem b
题目大意
求\(\sum_{x=a}^b\sum_{y=c}^d[gcd(x,y)=k]\) \((1\le n,k\le5\times10^4,1\le a\le b\le5\times10^4,1\le c\le d\le5\times10^4)\)
分析
(三倍经验)
发现此题比上一题多了下界要求，而上一题已经求了下界均为\(1\)的特殊情况
考虑下界不都为\(1\)，
可以用类似二维前缀和的容斥思想
设\(s(i,j)\)表示上界分别为\(i,j\)，下界均为\(1\)
那么
\[s=s(b,d)-s(a-1,d)-s(b,c-1)+s(a-1,c-1) \]就很好求了
复杂度还是\(O(n+Q\cdot\sqrt{n})\)，不过常数大一些
AC·code
#include<bits/stdc++.h>
#define il inline
#define cs const
#define ri register
using namespace std;
il int rd(){
ri int x=0,f=1;ri char c=getchar();
while(c<'0'||c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0'&&c<='9') x=x*10+(c^48),c=getchar();
return x*f;
}
il void wt(long long x){
if(x<0) x=-x,putchar('-');
ri unsigned short a[15]={0},tl=0;
do{a[++tl]=x%10,x/=10;}while(x);
for(ri int i=tl;i;--i) putchar(a[i]+48);
}
cs int N=5e4;
int ss[N+5],mu[N+5],cnt,sum[N+5];
bool vs[N+5];
il void init(int n){
mu[1]=sum[1]=1;
for(ri unsigned int i=2;i<=n;++i){
if(!vs[i]) ss[++cnt]=i,mu[i]=-1;
for(ri unsigned int j=1;j<=cnt&&i*ss[j]<=n;++j){
vs[i*ss[j]]=1;
if(i%ss[j]==0) break;
mu[i*ss[j]]=-mu[i];
}
sum[i]=sum[i-1]+mu[i];
}
}
il int done(int n,int m,int d){
if(n>m) swap(n,m);
ri int as=0;
for(ri unsigned int l=1,r=0;l<=n;l=r+1){
r=min(n/(n/l),m/(m/l));
as+=1ll*(n/l/d)*(m/l/d)*(sum[r]-sum[l-1]);
}
return as;
}
signed main(){
int t=rd(),a,b,c,d,k,as=0;init(N);
while(t--){
a=rd(),b=rd(),c=rd(),d=rd(),k=rd();
as=done(b,d,k)-done(a-1,d,k)-done(b,c-1,k)+done(a-1,c-1,k);
wt(as),putchar('\n');
}
return 0;
}