【学习笔记】特殊数论函数求和

前言

已经不知道第几次修改了。修过很多锅，也陆续新增了少许内容。数学的魅力就在于，它是没有被穷尽的。我也不相信它会有被穷尽的一天。

Comment. 用 $\mathbb{P}$ 表示质数集，用 $■$ 表示证毕，用 Comment 表示注释。

$\text{min25}$ 筛 $/$ 洲阁筛

$\text{min25}$ 筛或洲阁筛的主要目的是求积性函数 $f(x)$ 的前缀和。——其实也可以拓展到部分非积性函数 [2]。

Comment. 其实它也不能叫做筛法吧，但是其思想确实是筛法的。

我们首先用直观的方法描述之。该算法的总体思路是：

• 第一步，求出质数 $p$ 对应的 $f(p)$ 的前缀和。
• 第二步，求出 $f(x)$ 的前缀和。

为什么会这样？因为质数的幂——相较于质数而言——数量很少。

小质数

一般的积性函数是复杂的。所以我们需要用 完全积性函数 去 拟合 原函数：只需要找到完全积性函数 $\hslash (x)$ 使得 $\hslash (p)=f(p)(p\in \mathbb{P})$ 不变。因为 $\text{step one}$ 只需要求出 ${\sum }_{p\in \mathbb{P}}f(p)$ 。

Comment. 其实也可以不是完全积性的，可见后文 $\text{zzt}$ 求和法。

用 $g(n,j)$ 表示 $[2,n]$ 去除前 $j$ 个质数的倍数之后，剩余数字的 $f(x)$ 之和。注意：要求 $g(n,j)$ 时刻包含 $[2,n]$ 之间所有质数的 $f(x)$，也就是说，去除的是质数的 “真倍数”（不包含质数本身）。

形式化地，记 { p } \{p\} {p} 为从小到大的质数（角标从 $1$ 开始），则
$$g(n,j)=\sum _{i=1}^{j}f(p_i)+\sum _{i=1}^n[p_1\nmid i\wedge p_2\nmid i\wedge \cdots \wedge p_j\nmid i]f(i)$$

我们的目标是求出 $g(n,+\infty )$，也就是只留下了质数。

考虑 $\mathtt{dp}$ 求解。转移过程就是去掉 $k{p}_j(k>1)$ 。显然 $k$ 不应该是前 $(j-1)$ 个质数的倍数，因为那是已经被移除的数。又因为 $f(x)$ 是完全积性函数，有 $f(k{p}_j)=f(p_j)f(k)$，于是有转移式
$$g(n,j)=g(n,j-1)-f(p_j)\left[g\left(\lfloor \frac{n}{p_j}\rfloor ,j-1\right)-\mathrm{sumf}(j-1)\right](p_j⩽\sqrt{n})$$

这里 $\mathrm{sumf}(j-1)={\sum }_{i=1}^{j-1}f(p_i)$，因为 $g$ 里面没有挖掉质数本身，减去它以修正。

特别地，如果 $n<p_j^2$，显然没有更多的数字需要被筛掉，此时 $g(i,j)=g(i,j-1)$ 。若在 $j$ 这一维上做滚动数组，则无需修改这些位置。同时可知，需要的质数最大是 $\sqrt{n}$，所以直接线性筛即可求出 $\mathrm{sumf}$ 了。

我们要对 $g(i,0)$ 赋初值。所以又涉及一个问题是，拟合函数的前缀和要易于求解。一般来说，我们会选择单项式作为拟合函数。不会真的有毒瘤题目需要用可以杜教筛的数论函数来拟合吧。

Comment. 多项式不总是完全积性的。将其拆解为多个单项式，分别求 $g$ 再相加即可。

$n$ 的范围很大。注意到求解 $g(n,j)$ 时，只递归到 $n^{\prime }=\lfloor \frac{n}{v}\rfloor$ 的值；最后我们会看到，若求解 $s$ 长度的前缀和，则我们只需求出 $n=\lfloor \frac{s}{v}\rfloor (v\in {\mathbb{N}}^{+})$ 的答案，这样的 $n$ 只有 $\mathcal{O}(\sqrt{s})$ 个。

具体复杂度计算见后文。

大质数

注意下面的方法中 $f(1)$ 总是会被算漏。记得最后加上。

$\text{min25}$ 筛

我们前面只求出了质数的求和是 $g(s,+\infty )$，所以现在我们只需要考虑合数了。

好消息：合数的最小质因数是不超过 $\sqrt{s}$ 的！于是咱可以枚举这玩意儿。记 $x$ 的最小质因数为 $\gamma (x)(x⩾2)$ 。

记 $h(n,j)={\sum }_{i=2}^n[\gamma (i)⩾p_j]f(i)$ 。欲求即 $h(s,1)$ 。

枚举最小质因子与其幂次，可以写出转移式
$$\begin{gathered} h(n,j)=g(n,+\infty )-\mathrm{sumf}(j-1) \\ +\sum _{i⩾j}\sum _{k=1}^{\lfloor {\log }_{p_i}n\rfloor }f(p_i^k)\left[h(\lfloor \frac{n}{p_i^k}\rfloor ,i+1)+[k\ne 1]\right] \end{gathered}$$

用 $g(n,+\infty )-\mathrm{sumf}(j-1)$ 计算了 $[2,n]$ 中的质数 $p_c(c⩾j)$，用 $f(p_i^k)\cdot h(\dots )$ 计算了含至少两个质因子的合数，用 $[k\ne 1]f(p_i^k)$ 计算了质数的幂，不重不漏。

与 $\text{step one}$ 同理，只在 $p_i⩽\sqrt{n}$ 时转移。有趣的是，我们可以直接 递归，而且 不进行记忆化。详细复杂度计算见后文。

洲阁筛

令 $h(n,j)$ 始终包含 $[2,s]$ 全体质数。这是利于转移的要求。
$$h(n,j)=h(n,j+1)+\sum _{k=1}^{\lfloor {\log }_{p_j}n\rfloor }f(p_j^k)\left[h(\lfloor \frac{i}{p_j^k}\rfloor ,k+1)-\mathrm{sumf}(j)+[k\ne 1]\right]$$

转移条件仍为 $p_j⩽\sqrt{n}$ 。初值为 $h(n,+\infty )=g(n,+\infty )$，相当于在 $g$ 的基础上接着递推。

没错，二者明明只是递归和递推的区别，但它的名字就变成了洲阁筛。小编也很好奇。

代码实现

嗯，直接上代码。以板题 $f(p^k)=p^k\cdot (p^k-1)$ 为例：需要用 $x^2$ 和 $x$ 分别拟合，于是将 $g$ 设为 $\text{pair}$ 类型。

这里的代码是递归的。在 补充信息 中的题目有递推版代码。

#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cctype>
using namespace std;
# define rep(i,a,b) for(int i=(a); i<=(b); ++i)
# define drep(i,a,b) for(int i=(a); i>=(b); --i)
typedef long long llong;
inline int readint(){
	int a = 0, c = getchar(), f = 1;
	for(; !isdigit(c); c=getchar())
		if(c == '-') f = -f;
	for(; isdigit(c); c=getchar())
		a = (a<<3)+(a<<1)+(c^48);
	return a*f;
}

const int MOD = 1e9+7, SQRTN = 100005;
int primes[SQRTN], primes_size;
bool isPrime[SQRTN];
struct Node{
	int one, two;
	Node() = default;
	Node(int x):one(x),two(int(llong(x)*x%MOD)){}
	Node(int _o,int _t):one(_o),two(_t){}
	Node operator * (const Node &t) const {
		return Node(int(llong(one)*t.one%MOD),
			int(llong(two)*t.two%MOD));
	}
	Node operator - (const Node &t) const {
		return Node((one-t.one+MOD)%MOD,(two-t.two+MOD)%MOD);
	}
	int val() const { return (two-one+MOD)%MOD; }
};
Node sumf[SQRTN];
void sieve(int n){
	memset(isPrime+2,true,n-1);
	for(int i=2; i<=n; ++i){
		if(isPrime[i]){
			primes[++ primes_size] = i;
			Node &v = sumf[primes_size] = sumf[primes_size-1];
			if((v.one += i) >= MOD) v.one -= MOD;
			v.two = int((v.two+llong(i)*i)%MOD);
		}
		for(int j=1; j<=primes_size&&primes[j]<=n/i; ++j){
			isPrime[i*primes[j]] = false;
			if(i%primes[j] == 0) break;
		}
	}
}

const int inv2 = (MOD+1)>>1, inv3 = (MOD+1)/3;
int haxi[2][SQRTN]; // value or divisor
inline int& index_(const llong &x,const llong &n){
	return (x < SQRTN) ? haxi[0][x] : haxi[1][n/x];
}
llong w[SQRTN<<1]; ///< positions to get sum
Node g[SQRTN<<1]; ///< 2 times SQRTN
void step_one(const llong &n){
	int tot = 0; ///< allocate index
	for(llong i=1; i<=n; i=n/(n/i)+1){
		w[++ tot] = n/i; index_(n/i,n) = tot;
		g[tot].one = int(((w[tot]%MOD+1)*(w[tot]%MOD)>>1)%MOD)-1;
		g[tot].two = int((w[tot]<<1|1)%MOD*(w[tot]%MOD+1)
			%MOD*(w[tot]%MOD)%MOD*inv2%MOD*inv3%MOD)-1;
	}
	for(int j=1; j<=primes_size; ++j){
		if(primes[j] > n/primes[j]) break;
		for(int i=1; i<=tot; ++i){
			if(primes[j] > w[i]/primes[j]) break;
			g[i] = g[i]-Node(primes[j])*(
				g[index_(w[i]/primes[j],n)]-sumf[j-1]);
		}
	}
}

inline llong func(const llong &x){
	return (x-1)%MOD*(x%MOD)%MOD; // definition
}
int step_two(const llong &x,int i,const llong &n){
	if(primes[i] > x) return 0;
	int res = g[index_(x,n)].val()-sumf[i-1].val();
	for(; i<=primes_size&&primes[i]<=x/primes[i]; ++i){
		llong t = primes[i], fk = x/t;
		for(; t<=fk; t*=primes[i])
			res = int((func(t*primes[i])+func(t)
				*step_two(x/t,i+1,n)+res)%MOD);
	}
	return (res >= 0) ? res : (res+MOD);
}

int main(){
	sieve(SQRTN-1); // just do it
	llong n; scanf("%lld",&n); step_one(n);
	printf("%d\n",step_two(n,1,n)+1);
	return 0;
}

时间复杂度

小质数

$\text{step one}$ 的复杂度比较好算，因为转移是 $\mathcal{O}(1)$ 的，只算转移次数。每个 $w$ 要用到不超过 $\sqrt{w}$ 的质数，所以复杂度约为

$$\sum _{i=1}^{\sqrt{n}}\frac{\sqrt{\frac{n}{i}}}{\ln \sqrt{\frac{n}{i}}}\approx {\int }_1^{\sqrt{n}}\frac{\sqrt{n}}{\sqrt{x}\ln \left(\frac{n}{x}\right)}\text{d}x$$

$2021/8/1\mathtt{update}$：由于被 $\mathsf{HandInDevil}$ 臭骂了一顿，稍微讲讲用积分计算复杂度的事儿。它实际上是将数论函数的求和，转化为了连续函数的积分。这是否能作为渐进复杂度呢？

一般而言，原数论函数是单调的，且不存在奇点。在此情形下，考虑用 ${\int }_{\delta }^{\delta +1}g(x)\text{d}x$ 拟合 $f(\delta )(\delta \in \mathbb{Z})$，则结果是偏小（或偏大）的；而用 ${\int }_{\delta -1}^{\delta }g(x)\text{d}x$ 去拟合 $f(\delta )$ 又会偏大（或偏小）。由于有 $0$ 是奇点的风险，上界（或下界）可能是 ${\int }_1^{n}g(x)\text{d}x+f(1)$，下界（或上界）是 ${\int }_2^{n+1}g(x)\text{d}x$ 。二者往往是等阶的，那么就 “夹逼” 出了确界。

现在回到正题上。由于我已知道了结果，我告诉你
$$\int \frac{1}{\sqrt{x}\ln (\frac{n}{x})}\text{d}x=\mathcal{O}\left[\frac{\sqrt{x}}{\ln (\frac{n}{x})}\right]$$

因为 $\ln x$ 型函数作分母比较特殊：右式的导数实际上是 $\frac{1}{2\sqrt{x}\ln (\frac{n}{x})}+\frac{1}{\sqrt{x}{\ln }^2(\frac{n}{x})}$，但是在大 $\mathcal{O}$ 表示法下就是左式。代回原式可知复杂度为
$$\mathcal{O}\left(\frac{n^{\frac{3}{4}}}{\ln n}\right)$$

$\text{min25}$ 筛

复杂度即 $h$ 计算中的求和次数和。这和 $[k\ne 1]f(p_i^k)$ 被计入的次数相同。

注意到递归的本质是搜索，可以设这个 $f(p_i^k)$ 实际上来自 $f(l)$ 。设 $l=m{p}_i$，则 $p_i$ 恰为 $m$ 的最大质因子，因为 $k\ne 1$ 。另一方面，设 $\mathrm{big}(x)$ 为 $x$ 的最大质因子，则 $x\mathrm{big}(x)⩽s$ 必须在 $\mathrm{big}(x)$ 处统计贡献。因此它的复杂度就是 ${\sum }_{i=2}^s[i\mathrm{big}(i)]⩽s$ 。

Lemma. 对于实数 $\alpha \in (0,1)$，令 Q ( n ) = { i ⩽ n : big ⁡ ( i ) ⩽ i α } Q(n)=\{i\leqslant n:\operatorname{big}(i)\leqslant i^{\alpha}\} Q(n)={i⩽n:big(i)⩽iα}，则 $|Q(n)|\sim n\rho ({\alpha }^{-1})$，其中 $\rho$ 是 $\text{Dickman function}$ [3]。

若令 M ( n ) = { i : i big ⁡ ( i ) ⩽ n } M(n)=\{i:i\operatorname{big}(i)\leqslant n\} M(n)={i:ibig(i)⩽n}，则 $\forall a\in (0,1),|M(n)|=\Omega (n^{\alpha })$ 。因为 $\text{Lemma}$ 告诉我们 P = { i ⩽ n α : big ⁡ ( i ) ⩽ n 1 − α } P=\{i\leqslant n^\alpha:\operatorname{big}(i)\leqslant n^{1-\alpha}\} P={i⩽nα:big(i)⩽n1−α} 的大小是 $\Omega (n^{\alpha })$，而 $P\subseteq M(n)$ 。

另一方面，取 $t=\lceil \log \log n\rceil$，则 { x ⩽ n : big ⁡ ( x ) ⩽ p t } \{x\leqslant n:\operatorname{big}(x)\leqslant p_t\} {x⩽n:big(x)⩽pt​} 大小不超过 $(\log n{)}^t=o(\frac{n}{\log \log n})$，因为每个质因数的次数不超过 $\log n$ 。而其他的 $i\mathrm{big}(i)⩽n$ 的数满足 $i⩽\frac{n}{p_t}⩽\frac{n}{\log \log n}$，因此 $|M(n)|=\mathcal{O}(\frac{n}{\log \log n})$ 。

于是乎，我们得到 $|M(n)|=\Theta (n^{1-ϵ})$ 。也就是说，它的渐进复杂度是错误的（从严格的数学意义上）。

Comment. 这里 $n^{1-ϵ}$ 指优于 $n$ 又劣于 $n^{\alpha }(\alpha <1)$，应该容易理解。

但是打表发现，$n⩽10^{13}$ 时，对于 $p⩽n^{1/4}$，满足 $\mathrm{big}(i)=p$ 的 $i$ 的个数是 $\sqrt{n}$ 级别的，而一共只有 $\mathcal{O}(\frac{n^{1/4}}{\ln n})$ 个 $p$，因此其提供的贡献是 $\frac{n^{3/4}}{\ln n}$ 的 [3]。对于 $\mathrm{big}(i)⩾n^{1/4}$ 不难发现其提供贡献是 $\mathcal{O}(\frac{n^{3/4}}{\ln n})$ 的。

一句话：渐进复杂度是错误的，但是它跑得足够快。

洲阁筛

这就很简单了。质数密度 $\pi (n)\simeq \frac{n}{\ln n}$，很容易看出 ${\sum }_{j=2}^{\lfloor {\log }_{2}n\rfloor }\pi (\sqrt[j]{w})=o(\pi (n))$，因此复杂度为
$$\sum _{w=1}^{\sqrt{n}}\sum _{j=1}^{\log w}\pi \left(\sqrt[j]{n/w}\right)=\mathcal{O}\left(\sum _{w=1}^{\sqrt{n}}\pi (w)\right)=\mathcal{O}\left(\frac{n^{3/4}}{\ln n}\right)$$

这就是所谓 “质数的幂的数量是很少的——相较于质数而言”。

$\mathsf{zzt}$ 求和法

源自 [4]，也可以参阅 [5] 获得更多信息。

若原函数在素数处的取值可以被易于计算前缀和的数论函数拟合，则可以用该方法。后面会说到，这种方法有 $\text{min25}$ 筛的影子。

前置知识：狄利克雷生成函数和杜教筛。

考虑 $\mathrm{Dirichlet}$ 双曲线法计算 $f\ast g$ 在所有 $\lfloor \frac{n}{v}\rfloor$ 位置的前缀和。若 $f,g$ 的非零项密度为 $\mathcal{O}(\frac{1}{\ln n})$，则计算 $f\ast g$ 的前 $x$ 项的和的复杂度是 $\frac{\sqrt{x}}{\ln n}$，故总复杂度是
$$\sum _{x=1}^{\lfloor \frac{n}{S}\rfloor }\frac{\sqrt{\frac{n}{x}}}{\ln n}=\mathcal{O}\left(\frac{n^{2/3}}{\sqrt{S}\ln n}\right)$$

而前 $S$ 项是狄利克雷卷积求出的，则复杂度为 $\mathcal{O}(\frac{S\log S}{{\ln }^{2}n})$ 。仍取 $S=n^{2/3}$，但复杂度已经是 $\mathcal{O}(\frac{n^{2/3}}{\ln n})$ 了，而且不基于线性筛（不基于积性）。

考虑狄利克雷生成函数，设 $F_p$ 为质数 $p$ 上的狄利克雷生成函数，则
$$F(x)=\left[\prod _{p⩽n^{1/6}}{F}_p(x)\right]\left[\prod _{n^{1/6}<p⩽n^{1/2}}{F}_p(x)\right]\left[\prod _{n^{1/2}<p}{F}_p(x)\right]$$

先求第二部分。

定理：恰含 $k$ 个质因子，且最小质因子至少为 $n^{\alpha }(\alpha <0.5)$ 的 $n$ 以内的数的个数是 $\mathcal{O}[({\alpha }^{-1}-1{)}^{k}n{\ln }^{-1}n]$ 。简记为 ${\varphi }_k(n^{\alpha },n)$ 。

证明：当 $k=1$ 时显然成立。归纳法，枚举最小质因子 $x\in [2,\sqrt{n}]$，则 ${\varphi }_k(x,\frac{n}{x})=\mathcal{O}(\frac{({\alpha }^{-1}-1{)}^{k}n}{x\log x})$，于是估算
$${\varphi }_{k+1}(n^{\alpha },n)\approx {\int }_{n^{\alpha }}^{\sqrt{n}}\frac{\text{d}x}{\ln x}\frac{({\alpha }^{-1}-1{)}^k\cdot n}{x\ln x}=\mathcal{O}[({\alpha }^{-1}-1{)}^{k+1}n{\ln }^{-1}n]$$

即证。$■$

为了叙述方便，设第二部分对应的积性函数是 $f_2(x)$，并记 ${\mathbb{P}}_2=\mathbb{P}\cap (n^{1/6},n^{1/2}]$ 。

考虑忽略质因子顺序，每次任选某个质数，将其指数加一。设非积性的函数 $\text{inc}(p)=[p\in {\mathbb{P}}_2]f_2(p)$，我们从 $\text{inc}(x)$ 开始，递推求出 ${\text{inc}}^{k+1}={\text{inc}}^k\ast \text{inc}(k⩾1)$ 。

由定理一知 ${\text{inc}}^k(x)$ 非零项密度是 $\mathcal{O}(\frac{1}{\ln n})$ 。因此这样乘法的复杂度是 $\mathcal{O}(\frac{n^{2/3}}{\ln n})$ 。注意 $k⩽\lfloor \frac{1}{\alpha }\rfloor$ 是常数，因此可以全部求出。

将它们累加，得到 $g(x)={\sum }_{k⩾1}{g}_k(x)$ 。设 $n=\prod p_i^{t_i}$，不难发现，其值是二项式系数 $g(n)=\frac{(\sum t_i)!}{\prod (t_i!)}\prod \text{inc}(p_i{)}^{t_i}$ 。

于是构造积性函数 $\text{fit}(p^t)=\text{inc}(p{)}^t(t!{)}^{-1}(p\in {\mathbb{P}}_2)$，我们求出 $\frac{g(n)}{(\sum t_i)!}$ 为 $\text{fit}(x)$ 的前缀和。注意到 $\text{fit}(p)=\text{inc}(p)=f_2(p)$，因此用 $\text{powerful number}$ 的方法，可以 $\mathcal{O}(\sqrt{n})$ 从 $\text{fit}$ 前缀和过渡到 $f_2$ 的前缀和。

先考虑怎么算第一部分。直接依次加入质数即可。每次加入时 $\mathcal{O}(\sqrt{n}{\log }_{p}n)$ 暴力更新，复杂度
$$\begin{gathered} \sum _{p\in {\mathbb{P}}_1}\sqrt{n}{\log }_{p}n\sim {\int }_1^{n^{1/6}}\frac{\text{d}p}{\ln n}\sqrt{n}{\log }_{p}n \\ \sim \sqrt{n}\mathrm{li}(n^{1/6})=\mathcal{O}\left(\frac{n^{2/3}}{\ln n}\right) \end{gathered}$$

怎么算第三部分。用易于求解前缀和的数论函数来拟合之，用类似 $\text{min25}$ 的方式，筛掉小质数。但是我们直接求出小质数对应的积性函数，则可以用上面的方法（第二部分 $+$ 第一部分）做到 $\mathcal{O}(\frac{n^{2/3}}{\log n})$ 。

Comment. 若对 $p<\sqrt{n}$ 直接依次加入，则复杂度是 $\mathcal{O}(\frac{n^{3/4}}{\log n})$ 的。这似乎就是 $\text{min25}$ 呢（枚举每个小质数并更新 $\mathcal{O}(\sqrt{n})$ 个前缀和值），但是这就不再要求是积性函数了。估计跟第二部分原理相同。

然后现在的问题是，知 $f_3\ast a_1=a_2$，反推 $f_3$ 。由于 $f_3$ 在 $[1,\sqrt{n}]$ 内只有 $f_3(1)=1$，因此算 $w$ 前缀和时，需要枚举的部分只有 $\mathcal{O}(\frac{w}{\sqrt{n}})$ 长，易知其复杂度为 $\mathcal{O}(\sqrt{n}\log n)$ 。

算出来 $f_3$ 后，我们将其与 $f_1\ast f_2$ 合并。仍然利用 $f_3$ 在 $[1,\sqrt{n}]$ 内只有 $f_3(1)=1$ 非零的性质，合并还是 $\mathcal{O}(\sqrt{n}\log n)$ 的。

所以我们得到了总复杂度为 $\mathcal{O}(\frac{n^{2/3}}{\log n})$ 的，素数处取值是多项式的积性函数的求和法。但是，据原作者说，这个做法并没有效率提升 😂

$\text{SBT}$ 割线法

由于我能力较弱，这里就只做实用主义的探讨。

算法最初的目标是，求解 ${\sigma }_0(n)={\sum }_{d\mid n}1$ 约数个数函数的前缀和。

方案是，用折线去拟合曲线 $xy=n$，使得曲线下方（含边界）的整点就是折线下方（不含边界）的全体整点。

仍然用 $\mathrm{Dirichlet}$ 双曲线法，则 $x,y⩽\sqrt{n}$ 部分 $\mathcal{O}(1)$ 即得，因此只需计算 $y⩽\sqrt{n}$ 的部分。

从 $x=\lceil \frac{n+1}{y}\rceil ,y=1+\lfloor \sqrt{n}\rfloor$ 开始，每次找正整数 $p,q$ 满足 $(x+p)(y-q)>n$，让 $\frac{q}{p}$ 最大的前提下最小化 $p$ 。特别地，若 $y=1$ 则终止该过程。不难发现这样得到的折线符合要求且唯一存在。

怎么找？利用 $\text{Stern-Brocot Tree}$ 即可。其简单讲解在此。

一个重要的点在于，判断某子树内是否还有可能会用到的斜率。也就是说，目前得到两个分数 $\frac{a}{b},\frac{c}{d}$，已知 $\frac{a}{b}$ 可行而 $\frac{c}{d}$ 不可行，问是否存在互质的正整数 $p,q$ 使得 $\frac{pa+qc}{pb+qd}$ 可行。

代码给出的判断是：若 $(x+d{)}^2\cdot a⩾n\cdot b$ 就停止。尚不清楚其原理。

这个做法被证明是 $\tilde{\mathcal{O}}(n^{1/3})$ 的，最多是 $\mathcal{O}(n^{1/3}{\log }^{3}n)$ 的，据称可被证明是 $\mathcal{O}(n^{1/3}\log n)$ 。

补充信息

这道题是非积性函数，也可以做。内有递推版代码。

据说 $\text{min25}$ 有了新版，不知道啥情况，这里只 $\mathrm{mark}$ 一下。

$\mathrm{Reference}$

[1] jokerwyt. Min 25筛入门[EB/OL]. blog.csdn.net/j….

[2] GuessYCB. 关于min_25筛的一些理解[EB/OL]. cnblogs.com/G….

[3] 朱震霆. 一些特殊的数论函数求和问题[EB/OL]. private-link.

[4] whzzt. 积性函数求和问题的一种筛法[EB/OL]. blog.csdn.net/w….

[5] negiizhao. zzt求和法的简化版. negiizhao.blog.uoj.ac/blog/7165.