【数论】卷积反演大集合

不知道为啥脑抽要学数论,骂声一片中发现数论还没入门(悲)。

1. 狄利克雷卷积与数论函数#

1.1 数论函数#

定义:数论函数为值域为整数的函数。

简单数论函数:

  • I(n)恒等函数,恒等为 1
  • e(n)元函数,卷积中的单位元,若 n=1e(n)=1。否则为 e(n)=0
  • id(n)单位函数id(n)=n
  • idx(n)幂函数idx(n)=nx

在数论函数的基础上,有以下。

定义:完全积性函数f(1)=1,且对于任意的整数 abf(ab)=f(a)f(b)

e(n),I(n),id(n) 是完全积性函数。

定义:积性函数f(1)=1,且对于任意互质的 abf(ab)=f(a)f(b)

欧拉函数与莫比乌斯函数皆为积性函数。证明显然

1.2 狄利克雷卷积#

定义:(fg)(n)=d|nf(d)g(nd)

(fg)(n)=xy=nf(x)g(y)。(优美的对称结构)

交换律:fg=gf

结合律:(fg)h=f(gh)

会用就行。

2. 莫比乌斯反演#

2.1 莫比乌斯函数#

函数 fF 满足如下关系:

F(n)=d|nf(d)

由狄利克雷卷积得:F(n)=d|n(1f(d))

于是乎就有:F=If

若已知 fF,则需反演。

变形:f=FI1

这里的 I1 即为莫比乌斯函数,记为 μ

积性函数的逆还是积性函数,所以 μ 为积性函数。

μ(n)={1n=10n 含有平方因子(1)kk 为 n 的本质不同质因子个数

2.2 莫比乌斯反演#

  • 嵌套式反演公式:由 μI=ed|nμ(d)=[n=1]

好了,你现在已经数论入门了!!

3. 欧拉反演#

3.1 欧拉函数#

φ(n) 表示小于 n 的正整数中与 n 互质的数的个数。

计算方法:将总数减去与 n 不互质的数的个数。将 n 质因数分解后,记 p1,p2,p3,,pkn 的质因子,则有:

φ(n)=ni=1nnpi+1i<jnnpipj1i<j<knnpipjpk++(1)k+1np1p2p3pk

因数分解可得:

φ(n)=n(11p1)(11p2)(11p3)(11pk)

计算欧拉函数的步骤中运用到了容斥原理,可见容斥原理的应用范围非常广。

欧拉函数的性质:

  1. gcd(a,b)=1,则有 φ(a)·φ(b)=φ(ab),即欧拉函数为积性函数

    证明:

    Aa 的质因子集合,集合元素记为 p1,p2,p3,,pkBb 的质因子集合,集合元素记为 P1,P2,P3,,Pm

    gcd(a,b)=1AB=.

    φ(a)=a(11p1)(11p2)(11p3)(11pk),φ(b)=b(11P1)(11P2)(11P3)(11Pm)

    φ(ab)=ab(11p1)(11p2)(11p3)(11pk)(11P1)(11P2)(11P3)(11Pm)

    φ(a)·φ(b)=φ(ab)

  2. p 为质数,φ(p)=p1.

    证明:

    p 为质数, 小于 p 的正整数中与 p 互质的数的个数为 p1.

  3. p 为质数,φ(p)=pkpk1.

    证明:

    p 为质数,p 不互质的数只能为 p 的倍数,共有 pk1 个。

    φ(p)=pkpk1.

不用怀疑了,全从数学大礼包上贺的

3.2 欧拉反演#

  • 嵌套式欧拉反演:由 φI=idd|nφ(d)=n

4. 反演习题#

P3455 [POI2007] ZAP-Queries#

Problem#

给定 a,b,c,求 i=1aj=1b[(i,j)=c].

Solve#

莫反好题。

转化形式:

i=1acj=1bc[(i,j)=1]

代入嵌入式莫比乌斯反演:

i=1acj=1bcd|(i,j)μ(d)

d|(i,j)d|i,d|j 可知:

i=1acj=1bcd|i,d|jμ(d)

交换枚举顺序:

d=1min(ac,bc)μ(d)d|iacd|jbc1

化简:

d=1min(ac,bc)μ(d)acdbcd

μ(d) 可以用前缀和处理,acdbcd 整除分块即可。

Code#

点击查看代码
#include <bits/stdc++.h>
#define int long long
#define H 19260817
#define rint register int
#define For(i,l,r) for(rint i=l;i<=r;++i)
#define FOR(i,r,l) for(rint i=r;i>=l;--i)
#define MOD 1000003
#define mod 1000000007

using namespace std;

inline int read() {
  rint x=0,f=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();}
  while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+(ch^48);ch=getchar();}
  return x*f;
}

void print(int x){
  if(x<0){putchar('-');x=-x;}
  if(x>9){print(x/10);putchar(x%10+'0');}
  else putchar(x+'0');
  return;
}

const int N = 5e4 + 10;

int T, a, b, d, prime[N], tot, mu[N], sum[N];

bool vis[N];

void Euler() {
  vis[1] = mu[1] = 1;
  for (int i = 2; i <= N-10; i++) {
    if(!vis[i]) prime[++tot] = i, mu[i] = -1;
    for (int j = 1; i * prime[j] <= N-10 && j <= tot; j++) {
      vis[i * prime[j]] = 1;
      mu[i * prime[j]] = (i % prime[j] == 0 ? 0 : -mu[i]);
      if(i % prime[j] == 0) break;
    } 
  }
}

signed main() {
  T = read();
  Euler();
  For(i,1,N) {
    sum[i] = sum[i-1] + mu[i];
  }
  while(T--) {
    a = read(), b = read(), d = read();
    int ans = 0, n = a / d, m = b / d;
    for (int l = 1, r = 0; l <= min(n, m); l = r + 1) {
      r = min(n / (n / l), m / (m / l));
      ans += (sum[r] - sum[l-1]) * (n / l) * (m / l);
    }
    cout << ans << '\n';
  }
  return 0;
}

P2522 [HAOI2011] Problem b#

Problem#

给定 a,b,c,d,k,求 i=abi=cd[(i,j)=k].

Solve#

根据上一题得到启发:

i=1nj=1m[(i,j)=k]=d=1min(nk,mk)μ(d)nkdmkd

若上题求得是一个 n×m(i,j) 矩阵,则此题求以 (a,c),(b,d) 为顶点的子矩阵。

在上一题的基础上二维前缀和求解即可。

Code#

点击查看代码
#include <bits/stdc++.h>
#define int long long
#define H 19260817
#define rint register int
#define For(i,l,r) for(rint i=l;i<=r;++i)
#define FOR(i,r,l) for(rint i=r;i>=l;--i)
#define MOD 1000003
#define mod 1000000007

using namespace std;

inline int read() {
  rint x=0,f=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();}
  while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+(ch^48);ch=getchar();}
  return x*f;
}

void print(int x){
  if(x<0){putchar('-');x=-x;}
  if(x>9){print(x/10);putchar(x%10+'0');}
  else putchar(x+'0');
  return;
}

const int N = 5e4 + 10;

int T, a, b, c, d, k, prime[N], tot, mu[N], sum[N];

bool vis[N];

void Euler() {
  vis[1] = mu[1] = 1;
  for (int i = 2; i <= N-10; i++) {
    if(!vis[i]) prime[++tot] = i, mu[i] = -1;
    for (int j = 1; i * prime[j] <= N-10 && j <= tot; j++) {
      vis[i * prime[j]] = 1;
      mu[i * prime[j]] = (i % prime[j] == 0 ? 0 : -mu[i]);
      if(i % prime[j] == 0) break;
    } 
  }
}

int calc(int n, int m) {
  int ans = 0;
  n /= k, m /= k;
  for (int l = 1, r = 0; l <= min(n, m); l = r + 1) {
    r = min(n / (n / l), m / (m / l));
    ans += (sum[r] - sum[l-1]) * (n / l) * (m / l);
  }
  return ans;
}

signed main() {
  T = read();
  Euler();
  For(i,1,N) {
    sum[i] = sum[i-1] + mu[i];
  }
  while(T--) {
    a = read(), b = read(), c = read(), d = read(), k = read();
    cout << calc(b, d) - calc(a-1, d) - calc(b, c-1) + calc(a-1, c-1) << '\n';
  }
  return 0;
}

P2158 [SDOI2008] 仪仗队#

Problem#

给定一个 n×n 的点矩阵,求位于左下角能看见多少个不被遮挡的点。

Solve#

(x,y) 能被看见,则满足 (i,j)=1

若存在点 (kx,ky) 则它一定会被 (x,y) 挡住。

所求即为 i=0n1j=0n1[(i,j)=1].

莫比乌斯反演求解。

注意细节:计算时按 2+i=1n1j=1n1[(i,j)=1] 计算。n=1 需要特判。

Code#

点击查看代码
#include <bits/stdc++.h>
#define int long long
#define H 19260817
#define rint register int
#define For(i,l,r) for(rint i=l;i<=r;++i)
#define FOR(i,r,l) for(rint i=r;i>=l;--i)
#define MOD 1000003
#define mod 1000000007

using namespace std;

inline int read() {
  rint x=0,f=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();}
  while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+(ch^48);ch=getchar();}
  return x*f;
}

void print(int x){
  if(x<0){putchar('-');x=-x;}
  if(x>9){print(x/10);putchar(x%10+'0');}
  else putchar(x+'0');
  return;
}

const int N = 4e4 + 10;

int n, prime[N], tot, mu[N], sum[N];

bool vis[N];

void Eular() {
  vis[1] = mu[1] = 1;
  for (int i = 2; i <= 4e4; i++) {
    if(!vis[i]) prime[++tot] = i, mu[i] = -1;
    for (int j = 1; i * prime[j] <= 4e4 && j <= tot; j++) {
      vis[i * prime[j]] = 1;
      mu[i * prime[j]] = (i % prime[j] == 0 ? 0 : -mu[i]); 
      if(i % prime[j] == 0) break;
    }
  }
}

int calc(int n) {
  int ans = 0; n--; 
  for (int l = 1, r; l <= n; l = r + 1) {
    r = n / (n / l);
    ans += (sum[r] - sum[l-1]) * (n / l) * (n / l); 
  }
  return ans + 2;
}

signed main() {
  n = read();
  Eular();
  For(i,1,n) sum[i] = sum[i-1] + mu[i];
  cout << (n == 1 ? 0 : calc(n)) << '\n';
  return 0;
}

P2398 GCD SUM#

Problem#

给定正整数 n,求 i=1nj=1ngcd(i,j).

Solve#

欧拉反演模板题。

原式:

i=1nj=1n(i,j)

嵌入式欧拉反演:

i=1nj=1nd|(i,j)φ(d)

d|(i,j)d|i,d|j 可知:

i=1nj=1nd|i,d|jφ(d)

转化形式:

d=1nφ(d)d|ind|jn1

化简:

d=1nφ(d)(nd)2

前缀和 + 整除分块即可求解。注意欧拉筛欧拉函数时,是以质数为讨论根据。

Code#

点击查看代码
#include <bits/stdc++.h>
#define int long long
#define H 19260817
#define rint register int
#define For(i,l,r) for(rint i=l;i<=r;++i)
#define FOR(i,r,l) for(rint i=r;i>=l;--i)
#define MOD 1000003
#define mod 1000000007

using namespace std;

inline int read() {
  rint x=0,f=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();}
  while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+(ch^48);ch=getchar();}
  return x*f;
}

void print(int x){
  if(x<0){putchar('-');x=-x;}
  if(x>9){print(x/10);putchar(x%10+'0');}
  else putchar(x+'0');
  return;
}

const int N = 1e5 + 10;

int n, prime[N], tot, phi[N], sum[N];

bool vis[N];

void Eular() {
  vis[1] = phi[1] = 1;
  for (int i = 1; i <= N - 1; i++) {
    if(!vis[i]) prime[++tot] = i, phi[i] = i-1;
    for (int j = 1; j <= tot && i * prime[j] <= N - 1; j++) {
      vis[i * prime[j]] = 1;
      phi[i * prime[j]] = (i % prime[j] == 0 ? phi[i] * prime[j] : phi[i] * (prime[j]-1));
      if(i % prime[j] == 0) break;
    }
  }
}

int calc(int n) {
  int ans = 0;
  for (int l = 1, r; l <= n; l = r + 1) {
    r = n / (n / l);
    ans += (sum[r] - sum[l-1]) * (n / l) * (n / l);
  }
  return ans;
}

signed main() {
  n = read();
  Eular();
  For(i,1,n) sum[i] = sum[i-1] + phi[i];
  cout << calc(n) << '\n';
  return 0;
}

作者:Daniel-yao

出处:https://www.cnblogs.com/Daniel-yao/p/18026382

版权:本作品采用「署名-非商业性使用-相同方式共享 4.0 国际」许可协议进行许可。

posted @   Daniel_yzy  阅读(31)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 地球OL攻略 —— 某应届生求职总结
· 周边上新:园子的第一款马克杯温暖上架
· Open-Sora 2.0 重磅开源!
· 提示词工程——AI应用必不可少的技术
· .NET周刊【3月第1期 2025-03-02】
more_horiz
keyboard_arrow_up light_mode palette
选择主题
menu
点击右上角即可分享
微信分享提示