BZOJ传送门
洛谷传送门
解析:
首先呢,这个既没有 g c d gcd gcd也没有 l c m lcm lcm的式子让我们很头疼啊。。。
但是我们有一个结论,在文章的最后会我给出一个证明:
d
(
i
j
)
=
∑
k
∣
i
∑
l
∣
j
[
gcd
(
k
,
l
)
=
1
]
d(ij)=\sum_{k\mid i}\sum_{l\mid j}[\gcd(k,l)=1]
d(ij)=k∣i∑l∣j∑[gcd(k,l)=1]
通过这个结论我们可以把式子大力化简一波:
∑
i
=
1
n
∑
j
=
1
m
d
(
i
j
)
=
∑
i
=
1
n
∑
j
=
1
m
∑
k
∣
i
∑
l
∣
j
[
gcd
(
k
,
l
)
=
1
]
\begin{aligned} \sum_{i=1}^{n}\sum_{j=1}^md(ij)=&\sum_{i=1}^{n}\sum_{j=1}^m\sum_{k\mid i}\sum_{l\mid j}[\gcd(k,l)=1] \end{aligned}
i=1∑nj=1∑md(ij)=i=1∑nj=1∑mk∣i∑l∣j∑[gcd(k,l)=1]
直接反演得到: A n s = ∑ i = 1 n ∑ j = 1 m ∑ k ∣ i ∑ l ∣ j ∑ d ∣ gcd ( k , l ) μ ( d ) = ∑ d = 1 min ( n , m ) μ ( d ) ∑ i = 1 n ∑ j = 1 m ∑ k ∣ i ∑ l ∣ j [ d ∣ gcd ( k , l ) ] \begin{aligned} Ans=&\sum_{i=1}^{n}\sum_{j=1}^m\sum_{k\mid i}\sum_{l \mid j}\sum_{d\mid\gcd(k,l)}\mu(d)\\ =&\sum_{d=1}^{\min(n,m)}\mu(d)\sum_{i=1}^n\sum_{j=1}^m\sum_{k\mid i}\sum_{l \mid j}[d\mid \gcd(k,l)] \end{aligned} Ans==i=1∑nj=1∑mk∣i∑l∣j∑d∣gcd(k,l)∑μ(d)d=1∑min(n,m)μ(d)i=1∑nj=1∑mk∣i∑l∣j∑[d∣gcd(k,l)]
考虑先枚举
k
,
l
k,l
k,l
A
n
s
=
∑
d
=
1
min
(
n
,
m
)
μ
(
d
)
∑
k
=
1
n
∑
l
=
1
m
[
d
∣
gcd
(
k
,
l
)
]
⌊
n
k
⌋
⌊
m
l
⌋
=
∑
d
=
1
min
(
n
,
m
)
μ
(
d
)
∑
k
=
1
⌊
n
d
⌋
∑
l
=
1
⌊
m
d
⌋
⌊
n
d
k
⌋
⌊
m
d
l
⌋
\begin{aligned} Ans=&\sum_{d=1}^{\min(n,m)}\mu(d)\sum_{k=1}^n\sum_{l=1}^m[d\mid \gcd(k,l)]\lfloor\frac{n}{k}\rfloor\lfloor\frac{m}l\rfloor\\ =&\sum_{d=1}^{\min(n,m)}\mu(d)\sum_{k=1}^{\lfloor\frac{n}d\rfloor}\sum_{l=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{n}{dk}\rfloor\lfloor\frac{m}{dl}\rfloor \end{aligned}
Ans==d=1∑min(n,m)μ(d)k=1∑nl=1∑m[d∣gcd(k,l)]⌊kn⌋⌊lm⌋d=1∑min(n,m)μ(d)k=1∑⌊dn⌋l=1∑⌊dm⌋⌊dkn⌋⌊dlm⌋
最后分配率拆出来就是我们需要的东西: A n s = ∑ d = 1 min ( n , m ) μ ( d ) ( ∑ k = 1 ⌊ n d ⌋ ⌊ n d k ⌋ ) ( ∑ l = 1 ⌊ m d ⌋ ⌊ m d k ⌋ ) Ans=\sum_{d=1}^{\min(n,m)}\mu(d)(\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{dk}\rfloor)(\sum_{l=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{m}{dk}\rfloor) Ans=d=1∑min(n,m)μ(d)(k=1∑⌊dn⌋⌊dkn⌋)(l=1∑⌊dm⌋⌊dkm⌋)
现在考虑怎么求这个东西: g ( n ) = ∑ i = 1 n ⌊ n i ⌋ g(n)=\sum_{i=1}^{n}\lfloor\frac{n}i\rfloor g(n)=i=1∑n⌊in⌋
首先你可以整除分块预处理出来所有的 g g g,复杂度 O ( n n ) O(n\sqrt n) O(nn)
但是我的目标预处理复杂度是 O ( n ) O(n) O(n)。
这里不是推式子的事情了,转换一下思路, ⌊ n 1 ⌋ \lfloor\frac{n}1\rfloor ⌊1n⌋是 [ 1 , n ] [1,n] [1,n]中有 1 1 1这个因子的数的个数, ⌊ n 2 ⌋ \lfloor\frac{n}2\rfloor ⌊2n⌋是 [ 1 , n ] [1,n] [1,n]中有 2 2 2这个因子的数的个数。。。
所以 g g g实际上就是 d d d的前缀和。。。
线性筛,然后整除分块,这道题就做完了。
一个小证明:
d ( i j ) = ∑ k ∣ i ∑ l ∣ j [ gcd ( k , l ) = 1 ] d(ij)=\sum_{k\mid i}\sum_{l\mid j}[\gcd(k,l)=1] d(ij)=k∣i∑l∣j∑[gcd(k,l)=1]
其实就是考虑怎么建立二元组 ( k , l ) (k,l) (k,l)和 i j ij ij的约数之间的一一对应的关系
首先每一个 i j ij ij约数 q q q,考虑其中的某个质因子 p p p,使得 q = s ⋅ p t q=s\cdot p^t q=s⋅pt
考虑 i = i ′ ⋅ p a , j = j ′ ⋅ p b i=i'\cdot p^a,j=j'\cdot p^b i=i′⋅pa,j=j′⋅pb,显然有 t ≤ a + b t\leq a+b t≤a+b
那么对于每个质因子 p p p,如果 t ≤ a t\leq a t≤a,那么我们在 k k k当中分配 t t t个质因子,不给 l l l质因子 p p p,否则我们在 k k k当中不保存质因子 p p p,给 l l l分配 t − a t-a t−a个质因子,显然,这样的分配下来 gcd ( k , l ) = 1 \gcd(k,l)=1 gcd(k,l)=1,并且每一个质因子对应唯一的分解方案,由此就构成了一一对应。
(写的有点口胡,语言不是很严谨,但是正确性是毋庸置疑的)
代码:
#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define re register
#define gc get_char
#define pc putchar
#define cs const
namespace IO{
namespace IOONLY{
cs int Rlen=1<<18|1;
char buf[Rlen],*p1,*p2;
}
inline char get_char(){
using namespace IOONLY;
return (p1==p2)&&(p2=(p1=buf)+fread(buf,1,Rlen,stdin),p1==p2)?EOF:*p1++;
}
inline int getint(){
re int num;
re char c;
while(!isdigit(c=gc()));num=c^48;
while(isdigit(c=gc()))num=(num+(num<<2)<<1)+(c^48);
return num;
}
inline void outint(ll a){
static char ch[23];
if(a==0)pc('0');
while(a)ch[++ch[0]]=a-a/10*10,a/=10;
while(ch[0])pc(ch[ch[0]--]^48);
}
}
using namespace IO;
cs int P=50004;
int prime[P],pcnt,mu[P];
int d[P],minprit[P];
bool mark[P];
inline void linear_sieves(int len=P-4){
d[1]=1;mu[1]=1;
for(int re i=2;i<=len;++i){
if(!mark[i])prime[++pcnt]=i,minprit[i]=1,d[i]=2,mu[i]=-1;
for(int re j=1;i*prime[j]<=len;++j){
mark[i*prime[j]]=true;
if(i%prime[j]==0){
mu[i*prime[j]]=0;
minprit[i*prime[j]]=minprit[i]+1;
d[i*prime[j]]=d[i]/(minprit[i]+1)*(minprit[i]+2);
break;
}
mu[i*prime[j]]=-mu[i];
d[i*prime[j]]=2*d[i];
minprit[i*prime[j]]=1;
}
mu[i]+=mu[i-1];
d[i]+=d[i-1];
}
}
inline ll solve(int n,int m){
if(n>m)swap(n,m);
ll ans=0;
for(int re i=1,j;i<=n;i=j+1){
j=min(n/(n/i),m/(m/i));
ans+=(ll)d[n/i]*d[m/i]*(mu[j]-mu[i-1]);
}
return ans;
}
signed main(){
linear_sieves();
for(int re T=getint();T--;outint(solve(getint(),getint())),pc('\n'));
return 0;
}