莫比乌斯反演和杜教筛

0.数论函数

数论函数是莫比乌斯反演和杜教筛的前置知识

0.定义

• 数论函数：定义域是正整数的函数
• 积性函数：对于任意两个互质的正整数$a,b$，满足$f(ab)=f(a)f(b)$的函数
• 完全积性函数：对于任意两个正整数$a,b$，满足$f(ab)=f(a)f(b)$的函数

1.常见积性函数

1.单位函数

非常简单的定义：$e(n)=[n=1]$

2.常函数

也很简单的定义：$1(n)=1$

3.恒等函数

依旧简单的定义：$id(n)=n$

4.莫比乌斯函数

莫比乌斯函数$\mu(n)$，当n有平凡因子的时候值是$0$，否则就是$(-1)^{质因子个数}$，即：

$$\mu(n)=\begin{cases} 1\ (n=1) \\ -1^k\ (n有k个质因子且没有平方因子) \\ 0\ (n有平方因子) \end{cases}$$

莫比乌斯函数的性质：

$$\sum_{i|n}\mu(i)=[n=1]$$

5.欧拉函数

欧拉函数$\phi(n)$的定义是$[1,n]$中和n互质的数的个数

假设集合$P$包含了所有n的质因子，那么有：

$$\phi(i)=n\prod_{p\in P}\frac{p-1}{p}$$

欧拉函数的性质：

$$\sum_{i|n}\phi(i)=n$$ $$若n\not=1,\sum_{i=1}^ni[(i,n)=1]=\frac{n\phi(n)}{2}(i)$$

2.Dirichlet卷积

对于两个数论函数$f(n),g(n)$，定义Dirichlet卷积为：

$$(f*g)(n)=\sum_{d|n}f(d)g(\frac{n}{d})$$

并且Dirichlet卷积满足交换律、结合律和分配率，即：

$$f*g=g*f \\ f*(g*h)=(f*g)*h \\ f*(g+h)=f*g+f*h$$

由Dirichlet卷积我们可以推导出：

如果函数$f,g$都是积性函数，那么$f*g$也是积性函数

如果我们无法快速（在$n\log(n)$以内）求出$f\ast g$

那么我们就可以用$O(n\log((n))$的时间求出我们想要的$f\ast g $

可以这样处理：

typedef vector<int> Poly;
Poly Dirichlet(Poly f, Poly g) {
  int len = f.size();
  Poly h(len);
  for (int i = 1; i * i < len; i++) {
    for (int j = i; j * i < len; j++) {
      if (i == j) h[i * j] += f[i] * g[j];
      else h[i * j] += f[i] * g[j] + f[j] * g[i];
    }
  }
  return h;
}

1.莫比乌斯反演

1.莫比乌斯恒等式

如果两个函数$f, g $，满足

$$f(n)=\sum_{d|n}g(d)$$

那么有：

$$g(n)=\sum_{d|n}f(d)\mu(\frac{n}{d})$$

也可以表示成：

$$g(n)=\sum_{n|d}f(d)\mu(\frac{d}{n})$$

但是后面这种形式应该不太常用

证明很简单，我们其实就是要证明：

$f=g\ast1$和$g=f \ast\mu$是等价的

证明其实很简单，第一个式子两边同时卷上$\mu$或者第二个式子两边同时卷上$1$

2.使用技巧

注意在Mobius反演的时候有几种常用技巧：

1. 交换和号
2. 根号分块（下底函数分块）
3. 添加性质良好便于计算的积性函数

2.杜教筛

杜教筛，很方便，但是有其局限性

一般情况是我们需要计算一个积性函数$f$的前缀和：

$$Sum(n)=\sum_{i=1}^nf(i)$$

于是我们可以构造一个有良好特殊性质的积性函数$g$，使得$f\ast g$可以快速算出并具有一定性质

然后我们可以快速算出：

$$\sum_{i=1}^n\sum_{d|i}g(d)f(\frac{i}{d})$$

把上面的式子交换一下和号得到：

$$\sum_{d=1}^n g(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i)$$

就变成了：

$$\sum_{d=1}^n g(d)Sum(\lfloor\frac{n}{d}\rfloor)$$

也就是说：

$$Sum(n)+\sum_{d=2}^n g(d)Sum(\lfloor\frac{n}{d}\rfloor)=\sum_{i=1}^n\sum_{d|i}g(d)f(\frac{i}{d})$$

所以：

$$Sum(n)=\sum_{i=1}^n\sum_{d|i}g(d)f(\frac{i}{d})-\sum_{d=2}^n g(d)Sum(\lfloor\frac{n}{d}\rfloor)$$

我们发现上面的这个式子右边的部分，一个可以快速计算，一个可以根号分块递归成子问题

注意我们需要用线性筛预处理一部分前缀和大概$O(n^{\frac{2}{3}})$个，这样时间效率是$O(n^{\frac{2}{3}})$的

常数比较大，主要是需要map来对ll值和前缀和进行映射

---

例题：

题目大意：

定义神树集是一个二维平面上的点集，使得存在至少一条直线恰好经过两个点

求$x\in[0,n],y\in[0,n]$的神树集个数$\mod 12345678$

思路：

首先有一个结论：只有点数小于2或者所有点共线的点集不是神树集

那我们就把这一部分答案容斥掉就可以了

首先我们可以枚举两个端点的横纵坐标差然后进行统计（这里没有统计平行坐标轴的）

$$sum=2*[\sum_{x=1}^n\sum_{y=1}^n2^{(x,y)-1}(n+1-x)(n+1-y)]$$

其中$(n+1-x)(n+1-y)$是可以选择的左边端点的数量

然后我们枚举gcd

$$sum=\sum_{g=1}^n2^g\sum_{x=1}^{\lfloor\frac{n}{g}\rfloor}\sum_{y=1}^{\lfloor\frac{n}{g}\rfloor}(n-1-gx)(n+1-gy)[(x,y)=1]$$

然后把后面拆开按照g的次数分类

变成了

$$sum=\sum_{g=1}^n(n+1)^2(2^g)F(\frac{n}{g})+(n+1)(2^gg)G(\frac{n}{g})+(2^gg^2)H(\frac{n}{g})$$

其中

$$F(n)=\sum_{x=1}^n\sum_{y=1}^n[(x,y)=1]=[2\sum_{i = 1}^n\phi(i)]-1$$ $$G(n)=\sum_{x=1}^n\sum_{y=1}^n(x+y)[(x,y)=1]$$

因为我们有：

$$\sum_{i=1}^ni[(i,n)=1]=\frac{n\phi(n)}{2} (n!=1)$$

所以可以对$G(n)$进行化简，把x和y的贡献分别来看，得到：

$$G(n)=[3\sum_{i=1}^ni\phi(i)]-1$$

同理，我们也可以对$H(n)$进行化简

$$H(n)=\sum_{i=1}^ni^2\phi(i)$$

至于有关$i=1$的细节问题，请自行自习思考

然后发现$\sum_{i=1}^n\phi(i),\sum_{i=1}^ni\phi(i),\sum_{i=1}^ni^2\phi(i) $都是可以用杜教筛算出来的

然后$\sum_{i=1}^n 2^i, \sum_{i=1}^n 2^ii,\sum_{i=1}^n 2^i i^2$都是可以用做差法在$O(\log n)/O(1)$的时间内算出来的

到这里我们已经完成了一大部分

剩下的就是对还没有讨论的情况进行分析

1. 首先因为我们要算的是不合法情况，所以应该去掉被包含合法的两个点的个数

   因此我们需要算出不和坐标轴平行的两个点点集的个数

   所有点对个数是$\frac{(n+1)^2((n+1)^2-1)}{2}$，去掉与坐标轴平行的个数是$2(n+1)\frac{(n+1)n}{2}=(n+1)^2n$

2. 其次因为在与坐标轴平行的情况中，点数大于2的情况我们没有减去

   这一部分的情况实际上我们只用考虑$2(n+1)$个单独的行就可以了

   对于每一行我们用所有的方案数减去小于3的方案数

   就是$2^{n+1}-\frac{n(n+1)}{2}-(n+1)-1$

3. 最后是所有方案数中点数小于2的也不合法，一共有$(n+1)^2+1$种

然后就直接用总的方案数$2^{(n+1)^2}$种减去不合法的就可以了

注意一下$2^{(n+1)^2}$快速幂的时候指数不能直接取mod

并且还有一堆恶心的地方需要分类讨论，因为模数真的很奇葩

---

#include<bits/stdc++.h>

using namespace std;

typedef long long ll;

const ll N = 5e6 + 10;
const ll Mod = 12345678;

ll add(ll a, ll b) {
  return (a + b) % Mod;
}

ll sub(ll a, ll b) {
  return ((a - b) % Mod + Mod) % Mod;
}

ll mul(ll a, ll b) {
  if (a > Mod) a %= Mod;
  if (b > Mod) b %= Mod;
  return a * b % Mod;
}

ll fast_pow(ll a, ll b) {
  ll res = 1;
  for (; b; b >>= 1, a = mul(a, a))
    if (b & 1) res = mul(res, a);
  return res;
}

struct Node {
  ll f, g, h;
} sum[N];
map<ll, Node> mp;

ll n, tot, pri[N], phi[N];
bool vis[N];

void getprime() {
  ll limit = min(n, N - 1); 
  phi[1] = 1;
  for (ll i = 2; i <= limit; i++) {
    if (!vis[i]) {
      phi[i] = i - 1;
      pri[++tot] = i;
    }
    for (ll j = 1; j <= tot && pri[j] * i <= limit; j++) {
      vis[i * pri[j]] = 1;
      if (i % pri[j] == 0) {
        phi[i * pri[j]] = phi[i] * pri[j];
        break;
      } else {
        phi[i * pri[j]] = phi[i] * (pri[j] - 1);
      }
    }
  }
  for (ll i = 1; i < N; i++) {
    sum[i].f = add(sum[i - 1].f, phi[i]);
    sum[i].g = add(sum[i - 1].g, mul(i, phi[i]));
    sum[i].h = add(sum[i - 1].h, mul(mul(i, i), phi[i]));
  }
}

ll geta(ll x) {
  if (x & 1) {
    return mul(x, (x + 1) / 2);
  } else {
    return mul(x / 2, x + 1);
  } 
}

ll getb(ll x) {
  if (x % 6 == 0) {
    return mul(x / 6, mul(x + 1, 2 * x + 1));
  } else if (x % 6 == 1) {
    return mul(x, mul((x + 1) / 2, (2 * x + 1) / 3));
  } else if (x % 6 == 2) {
    return mul(x / 2, mul((x + 1) / 3, 2 * x + 1)); 
  } else if (x % 6 == 3) {
    return mul(x / 3, mul((x + 1) / 2, 2 * x + 1));
  } else if (x % 6 == 4) {
    return mul(x / 2, mul(x + 1, (2 * x + 1) / 3));
  } else {
    return mul(x, mul((x + 1) / 6, 2 * x + 1));
  }
}

Node solve(ll x) {
  if (x < N) return sum[x];
  if (mp.count(x)) return mp[x];
  ll pref = geta(x), preg = getb(x), preh = mul(pref, pref);
  for (ll i = 2; i <= x; i++) {
    ll j = x / (x / i);
    Node cur = solve(x / i);
    pref = sub(pref, mul(j - i + 1, cur.f));
    preg = sub(preg, mul(sub(geta(j), geta(i - 1)), cur.g));
    preh = sub(preh, mul(sub(getb(j), getb(i - 1)), cur.h));
    i = j;
  }
  return mp[x] = (Node) {pref, preg, preh};
}

ll calcf(ll x) {
  return fast_pow(2, x + 1);
}

ll calcf(ll l, ll r) {
  return sub(calcf(r), calcf(l - 1));
}

ll calcg(ll x) {
  return mul(x - 1, calcf(x)); 
}

ll calcg(ll l, ll r) {
  return sub(calcg(r), calcg(l - 1));
}

ll calch(ll x) {
  return add(sub(mul(mul(x, x), calcf(x)), 2 * calcg(x)), calcf(x));
}

ll calch(ll l, ll r) {
  return sub(calch(r), calch(l - 1));
}

 main() {
  scanf("%lld", &n);
  getprime();
  ll resf = 0, resg = 0, resh = 0;
  for (ll i = 1; i <= n; i++) {
    ll j = n / (n / i);
    Node cur = solve(n / i);
    resf = add(resf, mul(calcf(i, j), 2 * cur.f - 1));
    resg = add(resg, mul(calcg(i, j), 3 * cur.g - 1));
    resh = add(resh, mul(calch(i, j), cur.h));
    i = j;
  }
  ll ans = add(sub(mul(mul(n + 1, n + 1), resf), mul(n + 1, resg)), resh);
  ans = sub(ans, mul(mul(n + 1, n + 2), geta(n)));
  ans = add(ans, mul(mul(n + 1, n + 1), n));
  ans = add(ans, mul(n * 2 + 2, sub(fast_pow(2, n + 1), add(n + 2, geta(n)))));
  ans = sub(sub(fast_pow(fast_pow(2, n + 1), n + 1), mul(n + 1, n + 1) + 1), ans); //指数不能直接取模 
  printf("%lld", ans);
  return 0;
}