逆元
扯一点没有多大用的东西
在数论里面，我们不把倒数叫做倒数，而叫做逆元（纯属装逼）
逆元的作用很大，先来看点easy的栗子
某些性质
a+b≡amodp+bmodp(modp) a + b ≡ a m o d p + b m o d p ( m o d p )
a−b≡amodp−bmodp(modp) a − b ≡ a m o d p − b m o d p ( m o d p )
ab≡amodp∗bmodp(modp) a b ≡ a m o d p ∗ b m o d p ( m o d p )
前面三个都是对的（不需要证明是因为我不会证）
但是
(a÷b)modp≠((amodp)÷(bmodp))modp ( a ÷ b ) m o d p ≠ ( ( a m o d p ) ÷ ( b m o d p ) ) m o d p
（÷是整除）
- 随便举个例子就知道它不成立
随便搞一下
这是不是说我们对除法下的大数模操作就GG了？naive
- 先来看一些小学级别的东东，首先
ax=1 a x = 1
- 明显 x=1a x = 1 a ，我再加一个 modp m o d p 的条件，变成了
ax≡1(modp) a x ≡ 1 ( m o d p )
那么现在 x x 就不一定等于 1a 1 a 了
- 对于 ax=1 a x = 1 ，我们把 x x 看成 1a 1 a ，而 ax≡1(modp) a x ≡ 1 ( m o d p ) 加了一个同余 p p
- 我们不把 x x 看成倒数，这时的 x x 为 a a 关于 p p 的逆元
注意：当且仅当 (a,p)=1 ( a , p ) = 1 且 p p 为质数时， x x 才存在
其他东东
又是一个栗子
∵(3×5)mod7=1 ∵ ( 3 × 5 ) m o d 7 = 1
∴3×5≡1(mod7) ∴ 3 × 5 ≡ 1 ( m o d 7 )
- 在这里， 5 5 和 13 1 3 的作用是一样的（乘3再模7都是1），所以，5是3关于7的逆元
注意一点，除了 α=5 α = 5 ，还有其他整数 α α 满足 3α≡1(mod7) 3 α ≡ 1 ( m o d 7 ) ，如 19 19 等
连百度都说逆元唯一，但是逆元不是唯一的！！！（感谢同学指出）
逆元的卵用
前面已经说了一些东东了，相信逆元很好懂
我们不妨把 x x 的逆元用 inv(x) i n v ( x ) 来表示，那么对于除法——
(a÷b)modp=(a×inv(b))modp=(amodp×inv(b)modp)modp ( a ÷ b ) m o d p = ( a × i n v ( b ) ) m o d p = ( a m o d p × i n v ( b ) m o d p ) m o d p
难搞的除法变成了乘法，水了很多
how to get 逆元？
way one
基础数论里面有两条著名的定理：欧拉定理和费马小定理
欧拉定理： aφ(n)≡1(modn)(gcd(a,n)=1) a φ ( n ) ≡ 1 ( m o d n ) ( g c d ( a , n ) = 1 )
费马小定理： ap−1≡1(modp)(gcd(a,p=1)) a p − 1 ≡ 1 ( m o d p ) ( g c d ( a , p = 1 ) )
（我两个都不会证）
明显，费马小定理是可以从欧拉定理推过来的
当 p p 为质数时 φ(p)=p−1 φ ( p ) = p − 1
我们在谈论逆元，就把费马小定理两边同除以一个p
ap−2≡1/a(modp) a p − 2 ≡ 1 / a ( m o d p )
注意到 1a 1 a 是个小数，而等号左边是个整数，别忘了我们在说的是数论倒数
想一想，我们发现应该把 1a 1 a 写成 inv(a) i n v ( a ) ——
inv(a)≡ap−2(modp) i n v ( a ) ≡ a p − 2 ( m o d p )
我们可以用快速幂求逆元，时间复杂度 O(log2p) O ( l o g 2 p )
code
long long inv(long long a,long long p)
{
long long t=1,b=p-2;
while (b)
{
if ((b%1)==1)t=t*a%p;
a=a*a%p;
b/=2;
}
return t;
}通常题目会让你求 109+7 10 9 + 7 取模，快到飞起
快速幂对于大多题目够了
way two
逆元还可以用扩展欧几里得来求
首先要知道贝祖定理，描述是
对于 a,b＞0 a , b ＞ 0 且 a,b∈N a , b ∈ N ，必然存在整数 x x 和 y y ，使得 ax+by=gcd(a,b) a x + b y = g c d ( a , b )
我们要求 a a 关于 p p 的逆元（条件是 a a 和 p p 互质且 p p 是质数）
代入贝祖定理，则
ax+py=1 a x + p y = 1
给上面那个式子做点操作，两边各 modp m o d p ，得到
axmodp+pymodp=1modp a x m o d p + p y m o d p = 1 m o d p ∴axmodp=1modp ∴ a x m o d p = 1 m o d p ∴ax≡1(modp) ∴ a x ≡ 1 ( m o d p )
所以x是a关于p的逆元（ y y 也是 p p 关于 a a 的逆元，同理）
这里用扩欧 exgcd(a,p,x,y) e x g c d ( a , p , x , y ) ，求出的 x x 即为 a a 的逆元
code
int exgcd(int a,int b,int &x,int &y)
{
if (!b)
{
x=1,y=0;
return a;
}
int ans=exgcd(b,a%b,y,x);
y-=x*(a/b);
return ans;
}
int inv(int a,int b,int p)
{
int x,y;
exgcd(a,b,x,y);
return (x%p+p)%p;
}扩展欧几里得优点是用途多，能求 gcd g c d 能求逆元能解线性方程
对于数论题是把利器，学会是必须的
way three
最后一种方法更简单更好打：递归
先来一发证明网上找的我怎么可能会证
设x=pmoda，y=p÷a 设 x = p m o d a ， y = p ÷ a
∴(x+ay)modp=pmodp=0 ∴ ( x + a y ) m o d p = p m o d p = 0
∴xmodp=−aymodp ∴ x m o d p = − a y m o d p
∴x⋅inv(a)modp=(−y)modp ∴ x · i n v ( a ) m o d p = ( − y ) m o d p
∴inv(a)=(p−y)⋅inv(x)modp ∴ i n v ( a ) = ( p − y ) · i n v ( x ) m o d p
∴inv(a)=(p−p÷a)⋅inv(pmoda)modp ∴ i n v ( a ) = ( p − p ÷ a ) · i n v ( p m o d a ) m o d p
有了这个，我们就可以用递归来求 inv(a) i n v ( a ) 了
递归出口 inv(1)=1 i n v ( 1 ) = 1 ，因为 1 1 的逆元就是 1 1
递归code
long long inv(long long a,long long p) //注意a要小于p，先把a%p一下
{
return a==1?1:(p-p/a)*inv(p%a,p)%p;
}是不是短到爆炸？
因为 pmoda p m o d a 后至多是 p2 p 2 ，所以时间复杂度是 O(log2p) O ( l o g 2 p )
这种方法还有一个好处，能用 O(n) O ( n ) 的时间预处理1~n的逆元
只不过把递归的式子改成递推而已
递推code
void init()
{
inv[1]=1;
for(int i=2;i<=n;i++)
{
inv[i]=(p-p/i)*1ll*inv[p%i]%p;
}
}个人感觉方法三更好用
易理解，码量低，时间复杂度也低，还能 O(n) O ( n ) 时间预处理
亦可赛艇
……
突然发现好像还是不怎么会扩欧嗯……
还是学习一下扩展欧几里得吧……
扩展欧几里得
普通欧几里得
普通欧几里得算法又称辗转相除法，内容为 gcd(a,b)=gcd(b,amodb) g c d ( a , b ) = g c d ( b , a m o d b )
我不会证于是证明是上网Co的
其实证明没什么必要
prove
先设 a>b a > b
a a 可以表示成 a=kb+r a = k b + r （ a,b,k,r∈N∗,r<b a , b , k , r ∈ N ∗ , r < b ）， ∴r=amodb ∴ r = a m o d b
设 d d 是 a a 和 b b 的一个公约数，那么 d|a,d|b d | a , d | b
把 r=a−kb r = a − k b 两边各除以 d d ，得到 rd=ad−kbd=m r d = a d − k b d = m
∵d|a,d|b∴m ∵ d | a , d | b ∴ m 为整数， ∴d|r ∴ d | r 也就是 d|amodb d | a m o d b
∴d ∴ d 是 b b 和 amodb a m o d b 的公因数
设 d d 是 b b 和 amodb a m o d b 的公因数中的任意一个
∴d|b,d|amodb ∴ d | b , d | a m o d b ，设 amodb=a−kb(k∈Z) a m o d b = a − k b ( k ∈ Z ) ， ∴d|a−kb ∴ d | a − k b
∴d|a ∴ d | a ， ∴d|a,d|b,d|amodb ∴ d | a , d | b , d | a m o d b ，就是说 a,b,amodb a , b , a m o d b 的所有公因数是一样的
公因数相同，最大公因数不也就相同了？ ∴gcd(a,b)=gcd(b,amodb) ∴ g c d ( a , b ) = g c d ( b , a m o d b )
code
long long gcd(long long x,long long y)
{
return !y?x:gcd(y,x%y);
}- 每 mod m o d 一次，大的那个数至多是原来一半，所以时间复杂度为 O(log2n) O ( l o g 2 n )
还有欧几里得平均迭代次数为 (12⋅ln2⋅lnN)π2+1.47,N=min(x,y) ( 12 · l n 2 · l n N ) π 2 + 1.47 , N = m i n ( x , y )
嗯……
学习扩欧
不会
ax+by=gcd(a,b)(a,b,x,y∈Z) a x + b y = g c d ( a , b ) ( a , b , x , y ∈ Z )
已经知道 a,b a , b ，求解能使该等式成立的一组 x,y x , y
扩欧就是搞这个的
不知道怎么求就滚蛋
设 ax+by=gcd(a,b) a x + b y = g c d ( a , b )
b=0 b = 0 时，明显 gcd(a,b)=a,∴x=1,y=0 g c d ( a , b ) = a , ∴ x = 1 , y = 0
但当 a,b≠0 a , b ≠ 0 时，设 ax1+by1=gcd(a,b),bx2+(amodb)y2=gcd(b,amodb) a x 1 + b y 1 = g c d ( a , b ) , b x 2 + ( a m o d b ) y 2 = g c d ( b , a m o d b )
∵gcd(a,b)=gcd(b,amodb),∴ax1+by1=bx2+(amodb)y2 ∵ g c d ( a , b ) = g c d ( b , a m o d b ) , ∴ a x 1 + b y 1 = b x 2 + ( a m o d b ) y 2
∴ax1+by1=bx2+(a−⌊ab⌋⋅b)y2=ay2+bx2+⌊ab⌋⋅by2 ∴ a x 1 + b y 1 = b x 2 + ( a − ⌊ a b ⌋ · b ) y 2 = a y 2 + b x 2 + ⌊ a b ⌋ · b y 2
根据某恒等定理 x1=y2,y1=x2−⌊ab⌋⋅y2 x 1 = y 2 , y 1 = x 2 − ⌊ a b ⌋ · y 2
所以 x1 x 1 和 y1 y 1 的值是由 x2 x 2 和 y2 y 2 来决定的
x2 x 2 和 y2 y 2 的值通过递归得到
完美
code
int exgcd(int a,int b,int &x,int &y)
{
if (!b)
{
x=1,y=0;
return a;
}
int ans=exgcd(b,a%b,y,x);
y-=x*(a/b);
return ans;
}从上面Co过来的
直接调用
exgcd(a,b,x,y)，结果是返回 gcd(a,b) g c d ( a , b )x x 和 y y 能干嘛？
扩欧有什么用？
- 三种用途
purpose one：解不定方程 ax+by=c a x + b y = c
当且仅当 gcd(a,b)|c g c d ( a , b ) | c 时，不定方程 ax+by=c a x + b y = c 有整数解
不然两边各除以 gcd(a,b) g c d ( a , b ) ，等号左边是整数，而等号右边就是个分数了
- 把 c c 乘上 gcd(a,b)c g c d ( a , b ) c ，方程化为
ax+by=gcd(a,b) a x + b y = g c d ( a , b )
直接上扩欧解出 x0,y0 x 0 , y 0
那么原方程的解就要除上 gcd(a,b)c g c d ( a , b ) c ，分别为 x1=x0cgcd(a,b),y1=y0cgcd(a,b) x 1 = x 0 c g c d ( a , b ) , y 1 = y 0 c g c d ( a , b )
而通解 x=x1+btgcd(a,b),y=y1−atgcd(a,b) x = x 1 + b t g c d ( a , b ) , y = y 1 − a t g c d ( a , b ) ，其中 t∈Z t ∈ Z
purpose two：解线性同余方程 axmodb=c a x m o d b = c
把原方程化为 ax+by=c a x + b y = c 的形式（其中 y y 通常 <0 < 0 ）
直接上扩欧解不定方程，求出一个解 x1=x0cgcd(a,b) x 1 = x 0 c g c d ( a , b )
那么通解 x=x1+btgcd(a,b) x = x 1 + b t g c d ( a , b ) ，其中 t∈Z t ∈ Z
而最小整数解为 x=(x1modbgcd(a,b)+bgcd(a,b))modbgcd(a,b) x = ( x 1 m o d b g c d ( a , b ) + b g c d ( a , b ) ) m o d b g c d ( a , b )
purpose three：求逆元
讲过了不讲
类欧几里得
类？
类欧也是一种算法
它的样子很像欧几里得，所以叫类欧
它可以解决 O(n) O ( n ) 时间的式子，时间复杂度还是普欧的 O(log2n) O ( l o g 2 n )
下面进入类欧的奇幻世界
一些东西
f(a,b,c,n)=∑i=0n⌊ai+bc⌋ f ( a , b , c , n ) = ∑ i = 0 n ⌊ a i + b c ⌋
g(a,b,c,n)=∑i=0ni⌊ai+bc⌋ g ( a , b , c , n ) = ∑ i = 0 n i ⌊ a i + b c ⌋
h(a,b,c,n)=∑i=0n⌊ai+bc⌋2 h ( a , b , c , n ) = ∑ i = 0 n ⌊ a i + b c ⌋ 2
m=⌊an+bc⌋ m = ⌊ a n + b c ⌋
随便看看
类欧可以把这三个明显 O(n) O ( n ) 用 O(log2n) O ( l o g 2 n ) 搞掉
我数学很不好，所以这坑要多久填完我不好说
下面的除法都看成整除，因为每个除法都打一次整除太蛋疼
f f
f(a,b,c,n)=∑i=0nai+bc f ( a , b , c , n ) = ∑ i = 0 n a i + b c
- 分类讨论
a=0 a = 0
f(a,b,c,n)=∑i=0n0⋅i+bc=∑i=0nbc f ( a , b , c , n ) = ∑ i = 0 n 0 · i + b c = ∑ i = 0 n b c
=bc(n+1) = b c ( n + 1 )
a≥c a ≥ c 或 b≥c b ≥ c
f(a,b,c,n)=∑i=0nai+bc=∑i=0n(aci+bc+amodc⋅i+bmodcc) f ( a , b , c , n ) = ∑ i = 0 n a i + b c = ∑ i = 0 n ( a c i + b c + a m o d c · i + b m o d c c )
=ac⋅n(n+1)2+bc(n+1)+f(amodc,bmodc,c,n) = a c · n ( n + 1 ) 2 + b c ( n + 1 ) + f ( a m o d c , b m o d c , c , n )
a≤c a ≤ c 且 b≤c b ≤ c
f(a,b,c,n)=∑i=0nai+bc= f ( a , b , c , n ) = ∑ i = 0 n a i + b c =
莫比乌斯反演
都说了定义没用的还不信
莫比乌斯反演是数论数学中很重要的内容
在许多情况下能够简化运算,可以用于解决很多组合数学的问题
这没个卵用
反正是个很奇妙的东西
并没有在网上找到莫比乌斯反演的详细解释
不墨迹了
懵逼乌斯繁衍
定义
f(n)=∑d|ng(d) f ( n ) = ∑ d | n g ( d )
莫比乌斯反演就是在已知 f f 的情况下反演求 g g
我们可以试着用 f f 来表示 g g ，则有
g(n)=∑d|nμ(nd)×f(d)=∑d|nμ(d)×f(nd) g ( n ) = ∑ d | n μ ( n d ) × f ( d ) = ∑ d | n μ ( d ) × f ( n d )
我不会证明
μ μ 是什么呢？
关于 μ μ
不必知道 μ μ 是什么东西
这只是我们自己YY出来的、觉得这里应该有的一个函数只需要知道 μ μ 就是莫比乌斯函数
推一波 f(n)=∑d|ng(d) f ( n ) = ∑ d | n g ( d ) 可以发现一些东东
白内障看不清就放大
看得出 g(n)=∑d|n???(与nd有关)×f(d) g ( n ) = ∑ d | n ? ? ? ( 与 n d 有 关 ) × f ( d )
其实 ??? ? ? ? 就是 μ(nd) μ ( n d )
容易发现 μ μ 只能取 {−1,0,1} { − 1 , 0 , 1 }
更准确地给 μ μ 来取值，那么
μ(n)=⎧⎩⎨⎪⎪1(−1)k0n=1n=∏ki=1pka2|n,a∈N,a>1(1) (1) μ ( n ) = { 1 n = 1 ( − 1 ) k n = ∏ i = 1 k p k 0 a 2 | n , a ∈ N , a > 1
上面 p p 为互质的质数，这个应该看得懂
暴力求 μ μ 肯定很慢，怎么快速求 μ μ 呢？
注意 μ μ 的性质
μ μ 的性质一
∑k|nμ(k)=[n=1] ∑ k | n μ ( k ) = [ n = 1 ]
不会证
这个东西对接下来推式子有很大的帮助
必须记住
求 μ μ
μ μ 性质二：是积性函数，但不是完全积性
思考对于数 i i ，若 i i 为质数则 μ(i)=−1 μ ( i ) = − 1
若 imodp[j]==0 i m o d p [ j ] == 0 ，则 i∗p[j] i ∗ p [ j ] 有平方因子，所以 μ(i∗p[j])=0 μ ( i ∗ p [ j ] ) = 0
否则 μ(i∗p[j])=−μ(i) μ ( i ∗ p [ j ] ) = − μ ( i )
code
线筛里面 ϕ,μ ϕ , μ 都可以求出来
void init()
{
memset(bz,1,sizeof(bz));
mu[1]=phi[1]=1,tot=0;
for (int i=2;i<MAXN;i++)
{
if (bz[i])p[tot++]=i,phi[i]=i-1,mu[i]=-1;
for (int j=0;j<tot && i*p[j]<=MAXN;j++)
{
bz[i*p[j]]=0;
if (i%p[j]==0)
{
mu[i*p[j]]=0;
phi[i*p[j]]=phi[i]*p[j];
break;
}
mu[i*p[j]]=-mu[i],phi[i*p[j]]=phi[i]*(p[j]-1);
}
}
for (int i=1;i<MAXN;i++)pre[i]=pre[i-1]+mu[i];//μ的前缀和，会有用的
}这种数学的东东还是要草稿纸什么才行
直接上反演公式套路
默认下面的 n n 小于 m m ，除都是整除
繁衍套路1
求
∑i=1n∑j=1m[gcd(i,j)=1] ∑ i = 1 n ∑ j = 1 m [ g c d ( i , j ) = 1 ]
推倒
想到性质一也就是 ∑k|nμ(k)=[n=1] ∑ k | n μ ( k ) = [ n = 1 ]
把它放进上面的式子里面
=∑i=1n∑j=1m∑k|gcd(i,j)μ(k) = ∑ i = 1 n ∑ j = 1 m ∑ k | g c d ( i , j ) μ ( k )
想一下 μ(k) μ ( k ) 出现了多少次呢？ nk∗mk n k ∗ m k 次！！
所以
=∑k=1nμ(k)nkmk = ∑ k = 1 n μ ( k ) n k m k
- 为了美观把 k k 换成 i i ，下面不再强调了，最后
=∑i=1nμ(i)nimi = ∑ i = 1 n μ ( i ) n i m i
暴力只能 O(n) O ( n ) 求了，但是我们可以分块做到根号复杂度
先看套路2
繁衍套路2
求
∑i=1n∑j=1m[gcd(i,j)=d] ∑ i = 1 n ∑ j = 1 m [ g c d ( i , j ) = d ]
推倒
和套路1差不多
根据容斥，同时除上一个 d d
=∑i=1nd∑j=1md[gcd(i,j)=1] = ∑ i = 1 n d ∑ j = 1 m d [ g c d ( i , j ) = 1 ]
=∑i=1n∑j=1m∑k|gcd(i,j)μ(k) = ∑ i = 1 n ∑ j = 1 m ∑ k | g c d ( i , j ) μ ( k )
- 然后根据我们已经推出来了的
=∑i=1ndμ(i)nidmid = ∑ i = 1 n d μ ( i ) n i d m i d
暴力时间复杂度 O(n) O ( n )
开始分块
我们容易知道的是， nd n d 最多只有 n−−√ n 个取值（不会证）
同理 md m d 最多也只有 m−−√ m 个取值，所以 nid n i d 和 mid m i d 同时不变的段数有 n−−√+m−−√ n + m 个
先看一幅图（网上盗的）
这里 n=32,m=40,k=2 n = 32 , m = 40 , k = 2 ，红色为 nkd n k d 和 mkd m k d 不变的整段
看懂了么？
既然这些段 nkd n k d , mkd m k d 都不变，意味什么？
我们取 μ μ 的前缀和，然后拿每一段的最后一位的 μ μ 的前缀和乘上 nkd∗mkd n k d ∗ m k d 不就是这一段的答案和么？
分块时间复杂度 O(n−−√) O ( n )
code
int get(int n,int m,int d)
{
if (n>m)swap(n,m);
int ans=0;
n/=d,m/=d;
for (int i=1,last=1;i<=n;i=last+1)
{
last=min(n/(n/i),m/(m/i));
ans+=(pre[last]-pre[i-1])*(n/i)*(m/i);
}
return ans;
}例题1：【BZOJ2301】 [HAOI2011]Problem b
problem
analysis
反演例题
由容斥可得， ans=ans(b,d)−ans(a−1,d)−ans(b,c−1)+ans(a−1,c−1) a n s = a n s ( b , d ) − a n s ( a − 1 , d ) − a n s ( b , c − 1 ) + a n s ( a − 1 , c − 1 )
单分块还不会？上面刚刚推出来的
于是每一个分块都是 O(n−−√+m−−√) O ( n + m ) 的时间，总 O(n(n−−√+m−−√)) O ( n ( n + m ) ) 时间复杂度
code
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define MAXN 10000001
using namespace std;
bool bz[MAXN];
int p[MAXN],phi[MAXN],mu[MAXN],pre[MAXN];
int n,a,b,c,d,k,tot;
int read()
{
int x=0,f=1;
char ch=getchar();
while (ch<'0' || '9'<ch)
{
if (ch=='-')f=-1;
ch=getchar();
}
while ('0'<=ch && ch<='9')
{
x=x*10+ch-'0';
ch=getchar();
}
return x*f;
}
void init()
{
memset(bz,1,sizeof(bz));
mu[1]=phi[1]=1,tot=0;
for (int i=2;i<MAXN;i++)
{
if (bz[i])p[tot++]=i,phi[i]=i-1,mu[i]=-1;
for (int j=0;j<tot && i*p[j]<=MAXN;j++)
{
bz[i*p[j]]=0;
if (i%p[j]==0)
{
mu[i*p[j]]=0;
phi[i*p[j]]=phi[i]*p[j];
break;
}
mu[i*p[j]]=-mu[i],phi[i*p[j]]=phi[i]*(p[j]-1);
}
}
for (int i=1;i<MAXN;i++)pre[i]=pre[i-1]+mu[i];
}
int get(int n,int m,int d)
{
if (n>m)swap(n,m);
int ans=0;
n/=d,m/=d; for (int i=1,last=1;i<=n;i=last+1) { last=min(n/(n/i),m/(m/i));
ans+=(pre[last]-pre[i-1])*(n/i)*(m/i);
}
return ans;
}
int main()
{
//freopen("readin.txt","r",stdin);
init();
n=read();
while (n--)
{
a=read(),b=read(),c=read(),d=read(),k=read();
printf("%d\n",get(b,d,k)-get(a-1,d,k)-get(b,c-1,k)+get(a-1,c-1,k));
}
return 0;
}繁衍套路3
求
∑i=1n∑j=1mgcd(i,j) ∑ i = 1 n ∑ j = 1 m g c d ( i , j )
推倒
- 枚举约数 d d
=∑d=1nd∑i=1n∑j=1m[gcd(i,j)=d] = ∑ d = 1 n d ∑ i = 1 n ∑ j = 1 m [ g c d ( i , j ) = d ]
=∑d=1nd∑i=1nd∑j=1md[gcd(i,j)=1] = ∑ d = 1 n d ∑ i = 1 n d ∑ j = 1 m d [ g c d ( i , j ) = 1 ]
- 还是再加一个 ∑ ∑ 进去再变形
=∑d=1nd∑i=1nd∑j=1md∑k|gcd(i,j)μ(k) = ∑ d = 1 n d ∑ i = 1 n d ∑ j = 1 m d ∑ k | g c d ( i , j ) μ ( k )
=∑d=1nd∑i=1ndμ(i)nidmid = ∑ d = 1 n d ∑ i = 1 n d μ ( i ) n i d m i d
