min25筛学习笔记

min25筛

min25筛用于求一类数论函数的前缀和，适用于函数在素数处的取值可以用一个关于此素数的多项式来表示的数论函数。

处理质数部分

这部分我们需要解决$\sum _{p\subseteq prime}f(p)$，这里简单起见，假设$f(p)=p^t$

用$s_i$表示前$i$个质数之和，用$LP{F}_i$表示$i$的最小质因子，$p_i$表示第$i$个质数，设$$g(n,i)=\sum _{k=1}^n[LP{F}_k>p_i\vee k\subseteq prime]k^t$$

因此有:$$g(n,0)=\sum\limits_{k=1}^nkt$$

假设$\sqrt{n}$内有$cnt$个质数，那么

$$g(n,cnt)=\sum _{k=1}^n[k\subseteq prime]k^t$$

因此我们的目标就是求解$g(n,cnt)$，接下来考虑转移，即$g(n,i)$和$g(n,i-1)$的关系。注意到：

$$g(n,i-1)=\sum _{k=1}^n[LP{F}_k>p_i\vee k\subseteq prime]k^t+\sum _{k=1}^n[LP{F}_k=p_i\wedge k⊈prime]k^t$$

因此有转移：

• $p_i\ast p_i>n$
  $$g(n,i)=g(n,i-1)$$

• $p_i\ast p_i\le n$

$$\begin{array}{rl}g(n,i)& =g(n,i-1)-\sum _{k=1}^n[LP{F}_k=p_i\vee k⊈prime]k^t\\ & =g(n,i-1)-p_i^t\cdot \sum _{k=1}^{\lfloor \frac{n}{p_i}\rfloor }[LP{F}_k\ge p_i]k^t\\ & =g(n,i-1)-p_i^t\cdot [g(\lfloor \frac{n}{p_i}\rfloor ,i-1)-s_{i-1}]\end{array}$$

处理完整问题

设$S(n,i)=\sum _{k=1}^n[LP{F}_k>p_i]f(k)$，则我们的目标是求解$S(n,0)$

考虑如何转移：

• $p_i\ast p_i>n$

$$S(n,i-1)=g(n,i-1)-s_{i-1}$$

• $p_i\ast p_i\le n$

$$\begin{array}{rl}S(n,i-1)& =\sum _{p_j^e\le n\wedge j\ge i}f(p_j^e)\cdot ([e>1]+\sum _{k=1}^{\lfloor \frac{n}{p_j^e}\rfloor }[LP{F}_k>p_j]f(k))+g(n,cnt)-s_{i-1}\\ & =g(n,cnt)-s_{i-1}+\sum _{p_j^e\le n\wedge j\ge i}f(p_j^e)\cdot (S(\lfloor \frac{n}{p_j^e}\rfloor ,j)+[e>1])\left(1\right)\\ & =S(n,i)+\sum _{p_i^e\le n}f(p_i^e)\cdot (S(\lfloor \frac{n}{p_i^e}\rfloor ,i)+1)\left(2\right)\end{array}$$

可以直接按照$\left(1\right)$式暴力求解，也可以按照$\left(2\right)$式仿照质数部分DP递推求解。

时间复杂度分析

对于质数部分的求解，我们只关注$t\subseteq \{\lfloor \frac{n}{k}\rfloor |k\subseteq {\mathbb{N}}^{\ast }\}$处$g(t,i)$的取值，而只有当$p_i\ast p_i\le t$时才会产生转移。因此复杂度可以估计为：

$$\sum _{i=1}^{\sqrt{n}}\frac{\sqrt{i}}{\ln \sqrt{i}}+\sum _{i=1}^{\sqrt{n}}\frac{\sqrt{\frac{n}{i}}}{\ln \sqrt{\frac{n}{i}}}=O\left({\int }_1^{\sqrt{n}}\frac{\sqrt{\frac{n}{x}}}{{\log }_2\sqrt{\frac{n}{x}}}dx\right)=O\left(\frac{n^{\frac{3}{4}}}{{\log }_{2}n}\right)$$

对于第二部分的复杂度分析类似，值得一提的是，按照$\left(1\right)$式暴力求解的时间复杂度也是对的，而且这种写法更加简洁。

一些例子

• 求$\sum _{i=1}^n\sum _{i=1}^n{f}^k(\gcd (i,j))$其中$f(x)$表示$x$的次大质因子，重复的质因子计算多次，例如$f(4)=f(6)=2$规定$f(1)=0,f(p)=1$其中$p$为质数。$n,k\le 2\times {10}^9$

$$\begin{array}{rl}& \sum _{i=1}^n\sum _{i=1}^n{f}^k(\gcd (i,j))\\ & =\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }[\gcd (i,j)=1]f^k(d)\\ & =\sum _{d=1}^n{f}^k(d)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\mu (i)\cdot {\lfloor \frac{n}{d\cdot i}\rfloor }^2\\ & \mathrm{令}g(n)=\sum _{i=1}^n\mu (i)\cdot {\lfloor \frac{n}{i}\rfloor }^2，\mathrm{则}\mathrm{原}\mathrm{式}=\sum _{d=1}^n{f}^k(d)\cdot g\left(\lfloor \frac{n}{d}\rfloor \right)\\ & \mathrm{于}\mathrm{是}\mathrm{接}\mathrm{下}\mathrm{来}\mathrm{要}\mathrm{做}\mathrm{的}\mathrm{事}\mathrm{情}\mathrm{是}\mathrm{求}g\mathrm{的}\mathrm{前}\mathrm{缀}\mathrm{和}\mathrm{以}\mathrm{及}\mathrm{求}{f}^k\mathrm{的}\mathrm{前}\mathrm{缀}\mathrm{和}，g\mathrm{的}\mathrm{前}\mathrm{缀}\mathrm{和}\mathrm{可}\mathrm{以}\mathrm{用}\mathrm{数}\mathrm{论}\mathrm{分}\mathrm{块}+\mathrm{杜}\mathrm{教}\mathrm{筛}\mathrm{做}。\\ & \mathrm{然}\mathrm{后}\mathrm{考}\mathrm{虑}\mathrm{求}f\mathrm{的}\mathrm{前}\mathrm{缀}\mathrm{和}，\mathrm{设}S(n,i)=\sum _{j=1}^n[LP{F}_j>p_i\vee j\subseteq prime]f^k(j)\\ & \mathrm{若}{p}_i\cdot p_i>n,S(n,i-1)=\sum _{j=1}^n[j\subseteq prime]\\ & \mathrm{若}{p}_i\cdot p_i\le n,S(n,i-1)=S(n,i)+S(\lfloor \frac{n}{p_i}\rfloor ,i-1)-\sum _{j=1}^{\lfloor \frac{n}{p_i}\rfloor }[j\subseteq prime]+p_i\cdot ((\sum _{j=1}^{\lfloor \frac{n}{p_i}\rfloor }[j\subseteq prime])-i+1)\\ & \mathrm{对}\mathrm{于}n\mathrm{以}\mathrm{内}\mathrm{质}\mathrm{数}\mathrm{个}\mathrm{数}\mathrm{可}\mathrm{以}\mathrm{直}\mathrm{接}\mathrm{套}\mathrm{用}\mathrm{前}\mathrm{面}min25\mathrm{筛}\mathrm{的}\mathrm{质}\mathrm{数}\mathrm{部}\mathrm{分}\end{array}$$

code：

#include <cstdio>
#include <iostream>
#define int unsigned int

using namespace std;

const int M = 3e5 + 5;

int prime[M], cnt = 0, a[M], g[M << 1], f[M << 1], primeMi[M];
int n, k, v[M << 1], ss = 0, mu[M], mus[M];

inline int getIdx(int o) { return (o > M) ? (M + n / o) : o; }

int mi(int o, int p) {
    int yu = 1;
    while (p) {
        if (p & 1) yu *= o;
        o *= o;
        p >>= 1;
    }
    return yu;
}

int cal(int o) {
    if (o < M) return mu[o];
    if (mus[n / o]) return mus[n / o];
    int yu = 1;
    for (int l = 2, r; l <= o; l = r + 1) {
        r = o / (o / l);
        yu -= (r - l + 1) * cal(o / l);
    }
    return mus[n / o] = yu;
}

int find(int o) {
    int l = 1, r = ss, res = 0;
    while (l <= r) {
        int mid = (l + r) >> 1;
        if (o <= v[mid] / o) {
            res = mid;
            l = mid + 1;
        } else
            r = mid - 1;
    }
    return res;
}

signed main() {
    scanf("%u%u", &n, &k);
    mu[1] = 1;
    for (int i = 2; i < M; i++) {
        if (!a[i]) {
            prime[++cnt] = i;
            primeMi[cnt] = mi(prime[cnt], k);
            mu[i] = -1;
        }
        for (int j = 1; j <= cnt && i * prime[j] < M; j++) {
            a[i * prime[j]] = 1;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            } else {
                mu[i * prime[j]] = -mu[i];
            }
        }
    }
    mu[0] = 0;
    for (int i = 1; i < M; i++) mu[i] += mu[i - 1];
    for (int l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        v[++ss] = n / l;
        g[getIdx(n / l)] = n / l - 1;
    }
    for (int i = 1; i <= cnt; i++) {
        for (int j = 1; j <= ss; j++) {
            int t = v[j];
            if (prime[i] > t / prime[i]) break;
            g[getIdx(t)] -= g[getIdx(t / prime[i])] - i + 1;
        }
    }
    for (int j = 1; j <= ss; j++) {
        f[getIdx(v[j])] = g[getIdx(v[j])];
    }
    for (int i = cnt; i >= 1; i--) {
        int idx = find(prime[i]);
        for (int j = idx; j >= 1; j--) {
            int t = v[j];
            if (prime[i] > t / prime[i]) continue;
            int p1 = getIdx(t);
            int p2 = getIdx(t / prime[i]);
            f[p1] += f[p2] - g[p2] + primeMi[i] * (g[p2] - i + 1);
        }
    }
    int ans = 0;
    f[0] = 0;
    for (int l = 1, r; l <= n; l = r + 1) {
        int t = n / l, yu = 0;
        r = n / t;
        for (int l1 = 1, r1; l1 <= t; l1 = r1 + 1) {
            r1 = t / (t / l1);
            yu += (cal(r1) - cal(l1 - 1)) * (t / l1) * (t / l1);
        }
        ans += (f[getIdx(r)] - f[getIdx(l - 1)]) * yu;
    }
    cout << ans << endl;
    return 0;
}

• 对于正整数$n$，定义如下操作每次选择$p_i,p_j$满足$p_i|n\wedge p_j|n\wedge p_i\ne p_j$将$n$变为$\frac{n}{p_i\cdot p_j}$。求有多少n以内的正整数满足可以通过若干次操作变为$1$。$n\le {10}^{11}$

$$\begin{array}{rl}& \mathrm{首}\mathrm{先}\mathrm{可}\mathrm{以}\mathrm{通}\mathrm{过}\mathrm{若}\mathrm{干}\mathrm{次}\mathrm{操}\mathrm{作}\mathrm{变}\mathrm{为}1\mathrm{的}\mathrm{正}\mathrm{整}\mathrm{数}\mathrm{一}\mathrm{定}\mathrm{满}\mathrm{足}\mathrm{质}\mathrm{因}\mathrm{数}\mathrm{的}\mathrm{次}\mathrm{数}\mathrm{之}\mathrm{和}\mathrm{为}\mathrm{偶}\mathrm{数}\mathrm{且}\mathrm{质}\mathrm{因}\mathrm{数}\mathrm{的}\mathrm{最}\mathrm{大}\mathrm{幂}\mathrm{次}\mathrm{不}\mathrm{超}\mathrm{过}\mathrm{质}\mathrm{因}\mathrm{数}\mathrm{总}\mathrm{幂}\mathrm{次}\mathrm{的}\mathrm{一}\mathrm{半}。\\ & \mathrm{设}f(n,j,m,s)=\sum _{i=1}^n[LP{F}_i>p_j\wedge \mathrm{当}\mathrm{前}\mathrm{最}\mathrm{大}\mathrm{质}\mathrm{因}\mathrm{子}\mathrm{次}\mathrm{数}\mathrm{为}m，\mathrm{质}\mathrm{因}\mathrm{子}\mathrm{此}\mathrm{数}\mathrm{之}\mathrm{和}\mathrm{为}s\mathrm{且}\mathrm{加}\mathrm{入}i\mathrm{后}(m,s)\mathrm{能}\mathrm{到}\mathrm{达}\mathrm{合}\mathrm{法}\mathrm{状}\mathrm{态}]\\ & \mathrm{则}\mathrm{我}\mathrm{们}\mathrm{要}\mathrm{求}f(n,0,0,0)\\ & \mathrm{则}\mathrm{有}\mathrm{转}\mathrm{移}f(n,j,m,s)=[\mathrm{在}\mathrm{插}\mathrm{入}\mathrm{一}\mathrm{个}\mathrm{新}\mathrm{质}\mathrm{因}\mathrm{子}\mathrm{后}，(m,s)\mathrm{能}\mathrm{到}\mathrm{达}\mathrm{合}\mathrm{法}\mathrm{状}\mathrm{态}]\cdot ((\sum _{i=1}^n[i\subseteq prime])-j)+\\ & \sum _{i=j+1}\sum _{a=1}(f(\lfloor \frac{n}{p_i^a}\rfloor ,i,max(m,a),s+a)+[a>1\wedge \mathrm{插}\mathrm{入}{p}_i^a\mathrm{后}(m,s)\mathrm{能}\mathrm{到}\mathrm{达}\mathrm{合}\mathrm{法}\mathrm{状}\mathrm{态}])\end{array}$$

code:

#pragma GCC optimize(2)
#include <chrono>
#include <cstdio>
#include <ctime>
#include <iostream>

using namespace std;

const int M = 317005;

int prime[M], cnt = 0, a[M], ss = 0;
long long g[M << 1], n, v[M << 1];

int getIdx(long long o) { return (o > M) ? (n / o + M) : o; }

int check(int sum, int mx) { return ((mx << 1) <= sum) && ((sum & 1) == 0); }

long long S(long long o, int p, int sum, int mx) {
    if (o < prime[p + 1]) return 0;
    long long yu = (check(sum + 1, max(mx, 1))) ? (g[getIdx(o)] - p) : 0;
    for (int i = p + 1; i <= cnt; i++) {
        if (prime[i] > o / prime[i]) break;
        int t = 1;
        for (long long j = prime[i]; j <= o; j *= prime[i], t++) {
            yu += S(o / j, i, sum + t, max(mx, t));
            yu += check(sum + t, max(mx, t)) * (t > 1);
        }
    }
    return yu;
}

signed main() {
    for (int i = 2; i < M; i++) {
        if (!a[i]) prime[++cnt] = i;
        for (int j = 1; j <= cnt && i * prime[j] < M; j++) {
            a[i * prime[j]] = 1;
            if (i % prime[j] == 0) break;
        }
    }
    scanf("%lld", &n);
    for (long long l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        g[getIdx(n / l)] = n / l - 1;
        v[++ss] = n / l;
    }
    for (int i = 1; i <= cnt; i++) {
        for (int j = 1; j <= ss; j++) {
            long long t = v[j];
            if (prime[i] > t / prime[i]) break;
            g[getIdx(t)] -= g[getIdx(t / prime[i])] - i + 1;
        }
    }
    cout << S(n, 0, 0, 0) << endl;
    return 0;
}

• 求$\sum _{i=1}^n\sum _{j=1}^{m}f(i\cdot j)$,令积性函数$f(n)$满足$f(p^c)=p^{\gcd (c,k)}$，其中$p$为给定常数，$p$为素数，$c$为正整数。$n\le {10}^5,m\le {10}^{10},k\le {10}^9$

$$\begin{array}{rl}& \mathrm{令}F(n,x)=\sum _{i=1}^{n}f(i\cdot x),\mathrm{则}ans=\sum _{i=1}^{n}F(m,i)\\ & \mathrm{我}\mathrm{们}\mathrm{仿}\mathrm{照}min25\mathrm{筛}\mathrm{的}\mathrm{思}\mathrm{想}，\mathrm{每}\mathrm{次}\mathrm{提}\mathrm{取}\mathrm{出}\mathrm{一}\mathrm{个}\mathrm{质}\mathrm{数}\mathrm{的}\mathrm{贡}\mathrm{献}，\mathrm{假}\mathrm{设}\mathrm{提}\mathrm{取}\mathrm{质}\mathrm{因}\mathrm{子}p\mathrm{的}\mathrm{贡}\mathrm{献}，p\mathrm{在}x\mathrm{中}\mathrm{次}\mathrm{数}\mathrm{为}c\\ & F(n,x)=\sum _{p^e\le n}f(p^e)\cdot (F(\lfloor \frac{n}{p^e}\rfloor ,\frac{x}{p^c})-F(\lfloor \frac{n}{p^{e+1}}\rfloor ,\frac{x}{p^{c-1}}))\\ & \mathrm{这}\mathrm{样}\mathrm{看}\mathrm{起}\mathrm{来}\mathrm{状}\mathrm{态}\mathrm{数}\mathrm{有}O(n\cdot \sqrt{m}),\mathrm{但}\mathrm{如}\mathrm{果}\mathrm{每}\mathrm{次}\mathrm{都}\mathrm{取}x\mathrm{的}\mathrm{最}\mathrm{大}\mathrm{质}\mathrm{因}\mathrm{子}，\mathrm{状}\mathrm{态}\mathrm{数}\mathrm{和}\mathrm{转}\mathrm{移}\mathrm{数}\mathrm{都}\mathrm{在}\mathrm{可}\mathrm{接}\mathrm{受}\mathrm{的}\mathrm{范}\mathrm{围}\mathrm{内}\\ & \mathrm{极}\mathrm{限}\mathrm{数}\mathrm{据}\mathrm{下}\mathrm{状}\mathrm{态}\mathrm{数}\mathrm{仅}\mathrm{有}{10}^6\mathrm{级}\mathrm{别}，\mathrm{转}\mathrm{移}\mathrm{数}\mathrm{在}{10}^7\mathrm{次}\mathrm{级}\mathrm{别}。\\ & \mathrm{接}\mathrm{下}\mathrm{来}\mathrm{需}\mathrm{要}\mathrm{解}\mathrm{决}\mathrm{的}\mathrm{问}\mathrm{题}\mathrm{是}F(n,1)\\ & F(n,1)=\sum _{i=1}^{n}f(i),\mathrm{注}\mathrm{意}\mathrm{到}f\mathrm{和}id\mathrm{在}\mathrm{质}\mathrm{数}\mathrm{处}\mathrm{的}\mathrm{取}\mathrm{值}\mathrm{都}\mathrm{一}\mathrm{样}，f\mathrm{的}\mathrm{前}\mathrm{缀}\mathrm{和}\mathrm{可}\mathrm{以}\mathrm{很}\mathrm{方}\mathrm{便}\mathrm{的}\mathrm{用}PowerfulNumber\mathrm{解}\mathrm{决}。\\ & \mathrm{设}f=id\ast h,\mathrm{则}{S}_f(n)=\sum _{d\subseteq PN}h(d)\cdot \frac{\lfloor \frac{n}{d}\rfloor \cdot (\lfloor \frac{n}{d}+1\rfloor )}{2},h(p^c)=p^{\gcd (c,k)}-\sum _{t=1}^c{p}^t\cdot h(p^{c-t})\end{array}$$

#include <assert.h>

#include <algorithm>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <map>
#include <vector>
#define int long long

using namespace std;

const int mo = 1e9 + 7;
const int inv2 = (mo + 1) >> 1;
const int M = 1e5 + 5;
const int N = 2e6 + 5;

int n, m, prime[200005], cnt = 0, vis[N + 5], k, tmp = 0, pns[1000005],
                         cs[N + 5];
vector<int> h[10005];
vector<pair<int, int>> pn;
vector<int> pp[100005];
map<int, int> mp[100005];

inline int gcd(int o, int p) { return (!p) ? o : gcd(p, o % p); }

inline int getIdx(int o) { return (o <= M) ? o : (M + m / o); }

inline int Mod(int o) { return (o >= mo) ? (o - mo) : o; }

int mi(int o, int p) {
    int yu = 1;
    while (p) {
        if (p & 1) yu = yu * o % mo;
        o = o * o % mo;
        p >>= 1;
    }
    return yu;
}

vector<array<int, 3>> factor(int o) {
    vector<array<int, 3>> ans;
    for (int i = 1; i <= cnt; i++) {
        if (o % prime[i] == 0) {
            int yu = 0;
            int tt = 1;
            while (o % prime[i] == 0) {
                yu++;
                tt *= prime[i];
                o /= prime[i];
            }
            ans.push_back(array<int, 3>{yu, prime[i], tt});
        }
        if (prime[i] * prime[i] >= o) break;
    }
    if (o > 1) {
        ans.push_back(array<int, 3>{1, o, o});
    }
    return ans;
}

void init() {
    memset(pns, -1, sizeof(pns));
    memset(cs, -1, sizeof(cs));
    for (int i = 2; i <= N; i++) {
        if (!vis[i]) prime[++cnt] = i;
        for (int j = 1; j <= cnt && i * prime[j] <= N; j++) {
            vis[i * prime[j]] = 1;
            if (i % prime[j] == 0) {
                break;
            }
        }
    }
    for (int i = 1; i <= cnt; i++) {
        if (prime[i] > 100000) break;
        h[i].push_back(1);
        for (int j = prime[i], t = 1; j <= m; j *= prime[i], t++) {
            int sum = 0;
            for (int l = 0; l < t; l++) {
                sum = Mod(sum + mi(prime[i], t - l) * h[i][l] % mo);
            }
            h[i].push_back(Mod(mi(prime[i], gcd(t, k)) - sum + mo));
        }
    }
}

void pre(int o, int p, int q) {
    cs[o] = p % mo;
    for (int j = q; j <= cnt; j++) {
        if (o * prime[j] > N) break;
        for (int t = prime[j] * o, s = 1; t <= N; t *= prime[j], s++) {
            pre(t, p * mi(prime[j], gcd(s, k)) % mo, j + 1);
        }
    }
}

void getPN(int o, int p, int q) {
    pn.push_back(make_pair(o, p));
    for (int j = q; j <= cnt; j++) {
        if (prime[j] > 100000) break;
        if (prime[j] * prime[j] > m / o) break;
        for (int t = prime[j] * prime[j] * o, s = 2; t <= m;
             t *= prime[j], s++) {
            getPN(t, p * h[j][s] % mo, j + 1);
        }
    }
}

inline int s2(int o) { return o * (o + 1) % mo * inv2 % mo; }

int calH(int o) {
    int t = getIdx(o);
    if (pns[t] != -1) return pns[t];
    int ans = 0;
    for (auto [l, r] : pn) {
        if (l > o) break;
        ans = Mod(ans + r * s2(o / l % mo) % mo);
    }
    return pns[t] = ans;
}

int cal(int o, vector<array<int, 3>> p, int q) {
    if (o == 0) return 0;
    if (o == 1) return cs[q];
    if (p.empty()) return calH(o);
    if (o * q <= N) return pp[q][o];
    if (mp[q].find(o) != mp[q].end()) return mp[q][o];
    tmp++;
    int ans = 0;
    int q1 = q / p.back()[2];
    int q2 = q1 * p.back()[1];
    for (int i = 1, t = 0; i <= o; i *= p.back()[1], t++) {
        vector<array<int, 3>> p1 = p;
        vector<array<int, 3>> p2 = p;
        p1.pop_back();
        p2.pop_back();
        p2.push_back(array<int, 3>{1, p.back()[1], p.back()[1]});
        ans =
            (ans + mi(p.back()[1], gcd(k, t + p.back()[0])) *
                       (cal(o / i, p1, q1) - cal(o / i / p.back()[1], p2, q2)) %
                       mo) %
            mo;
        ans = Mod(ans + mo);
    }
    return mp[q][o] = ans;
}

signed main() {
    scanf("%lld%lld%lld", &n, &m, &k);
    init();
    pre(1, 1, 1);
    getPN(1, 1, 1);
    sort(pn.begin(), pn.end());
    for (int i = 1; i <= n; i++) {
        pp[i].push_back(0);
        for (int j = 1; i * j <= N; j++) {
            int yu = Mod(pp[i][pp[i].size() - 1] + cs[i * j]);
            pp[i].push_back(yu);
        }
    }
    int ans = 0;
    for (int i = 1; i <= n; i++) {
        ans = (ans + cal(m, factor(i), i)) % mo;
    }
    cout << ans << endl;
    return 0;
}