先定义几个符号：
[]：若方括号内为一个值，则向下取整，否则为布尔判断
集合P：素数集合。
题目分析：
题目是一个积性函数。做法之一是洲阁筛，也可以采用Min_25筛。
对于一个可以进行Min_25筛法的积性函数，它需要满足与洲阁筛相同的条件，即：
对于$f(p), p \in P$，它可以多项式表出。对于$f(p^k),p \in P$可以被快速计算出。
这道题中$f(p) = p-1$再对$2$进行修正即可。
对于1的情况我们单独考虑，现在我们对答案进行一些变换。
$$\sum_{i=2}^{n}f(i) = \sum_{2<=p^e<=n, p \in P} f(p^e)(1+\sum_{2<=i<=[\frac{n}{p^e}],i质因子大于p}f(i)) = \sum_{2<=p<=\sqrt n,p^e<=n, p \in P} f(p^e)(1+\sum_{2<=i<=[\frac{n}{p^e}],i质因子大于p}f(i))+\sum _{p \in P,\sqrt n<p<=n}f(p)$$
我们注意到中间刺眼的中文，不妨以此为根据设置dp状态。
令$$g_{n,m}=\sum_{i=2,i的质因子大于m}^{n}f(i)$$
它的递推式不难写出来。根据上面那串式子:
$$g_{n,m}=\sum_{i=m,i \in P}^{\sqrt n}[p^e \leq n]([e \neq 1]+\sum_{i=p+1,i的质因子大于p}^{[\frac{n}{p^e}]}f(i))+\sum_{i=m+1,i \in P}^{n}f(i) = \sum_{i=m,i \in P}^{\sqrt n}[p^e \leq n]([e \neq 1]+g_{[\frac{n}{p^e}],p})+\sum_{i=m+1,i \in P}^{n}f(i)$$ 之后我们可以发现后面那一部分等于$m+1$到$n$的所有素数.
朱老大的论文似乎在这个式子上有一些问题，这个式子是我加了修正的，可以对比原式找到区别。
现在我们重点考虑后面的部分，质数前缀和的求解。 我们可以用DP模拟埃拉托斯特尼筛法来完成。现在对f(p)这个多项式的每个单项式分开考虑。
令$$h_{n,m}=\sum_{i=2,i为素数或i质因子大于m}f(i)$$
由埃氏筛可知$$h_{n,m}=h_{n,m-1}-p_m^s*(h_{[\frac{n}{p}],m-1}-h_{p_m,m-1})$$起始状态将所有数看做质数，做k次幂和。
由于n以下的合数至多的因子大小是$\sqrt n$,所以对于一个n只用筛前$\sqrt n$个素数。
g数组同理。
这样我们就完成整个过程了，现在我们来思考时间复杂度上的问题。这一部分我的证明存在争议，与论文上的有出入。希望有人可以帮助我确认出入之处。
先给出两个引理吧。
素数定理:$\pi(x) \approx \frac{x}{lnx} \approx li(x)$
数论分块:对于一个数n，它对1到n整除的结果只有$2\sqrt n$个，分别是1到$\sqrt n$与$n/\sqrt n$到$n/1$
首先是dp模拟埃氏筛。对于数论分块中的每一个数，每个小于它的根号的素数都被考虑，所以时间为，
$$O(\sum_{i=2}^{\sqrt n}\frac{\sqrt i}{lni})+O(\sum_{i=1}^{\sqrt n}\frac{\sqrt \frac{n}{i}}{ln{\frac{n}{i}}})$$
对于前面的大O,将$\sqrt i$放缩到$\sqrt n$,$lni$放缩到$lnn$,因为根号的增长远大于对数，所以放缩成立。前面被证明是$O(\frac{n^{0.75}}{lnn})$的。
后面的比较复杂，将它看做积分，即$$O(\int_1^{\sqrt n}\frac{\sqrt \frac{n}{x}}{ln\frac{n}{x}}dx)$$
由于对数函数增长极慢，不妨将分母用同样的放缩。当然也可以查积分表得到与我相同的结论。这样再对单项式函数做积分就简单了，得到$$O(\frac{\sqrt n}{lnn}*[2\sqrt x]_1^{\sqrt n})=O(\frac{n^{0.75}}{lnn})$$。
埃氏筛的复杂度得证，我们再来看看Min_25筛的表现。论文中通过证明得出一个与运行效率不相符的时间复杂度，我试着用其它方法得到它的复杂度。
观察Min_25的特点，它与DP埃氏筛的区别在于多枚举了质数的次幂。这样我们不妨对于质数的次幂分开考虑。
首先是质数本身，如果只枚举它，时间与埃氏筛DP相同。实际上它代表的是二次方，否则不被考虑。
然后是质数的平方，它要小于等于n。接着是质数的三次方，它也要小于等于n，然后是四次方，五次方。我们对第x次方小于等于n做一些转化使得他们统一。
$p^x \leq n \leftrightarrow p \leq n^{1/x} \leftrightarrow p^2 \leq n^{2/x}$
这时候我们的复杂度相当于对$n^{2/2},n^{2/3},n^{2/4}$做埃氏筛DP。每个的时间在上面分析了，所以下面直接转化成积分的形式。
$$O(\int_{2}^{log_{2}n}\frac{n^{\frac{3}{2x}}}{ln{n^{\frac{2}{x}}}}dx)$$
$$=O(\frac{1}{lnn}\int_{2}^{log_{2}n}\frac{n^{\frac{3}{2x}}}{\frac{2}{x}}dx)$$
$$=O(\frac{1}{2lnn}\int_{2}^{log_{2}n}{x*n^{\frac{3}{2x}}}dx)$$
$$=O(\frac{1}{2lnn}\int_{2}^{log_{2}n}{x*a^{\frac{1}{x}}}dx),a = n^{\frac{3}{2}}$$
接下俩数学好的可以换元，令$y=\frac{1}{x}$,再令$z = a^y$，求$$\int \frac{1}{ln^3{z}dz}$$即可。
现在跳过这些步骤，得到结论是
$$O([\frac{x(x+\ln a)a^\frac{1}{x}-\ln^2 a\operatorname{li}(a^{\frac{1}{x}})}{4lnn}]_{2}^{log_{2}n})$$
这个答案是负的，证明我可能某个步骤写错了。但是代入2可以得到时间是$O(n^{\frac{3}{4}})$，递归的过程不会重复计算。
继续写下去我也写不动了。到此打止。
代码：
#include<bits/stdc++.h>
using namespace std; const int SQR = ;
const int mod = ; long long n; int prime[SQR],flag[SQR+],num;
long long h[SQR*+],sqr;
long long fh[SQR*+]; void getprime(int N){
for(int i=;i<=N;i++){
if(!flag[i]) prime[++num] = i;
for(int j=;j<=num&&i*prime[j] <= N;j++){
flag[i*prime[j]] = ;
if(i%prime[j] == ) break;
}
}
} void geth(int ft){
for(int i=;i<=num;i++){int zt = (ft > ?prime[i]:);
for(int j=;j<=sqr&&(1ll*prime[i]*prime[i]<=(n/j));j++){
long long nxt = 1ll*j*prime[i]; int plu = *sqr-j+;
if(nxt > sqr) h[plu] -= 1ll*zt*(h[n/nxt]-h[prime[i]-])%mod;
else h[plu] -= 1ll*zt*(h[*sqr-nxt+]-h[prime[i]-])%mod;
h[plu]<?h[plu]+=mod:;
}
for(int j=sqr;j>=&&(1ll*prime[i]*prime[i]<=j);j--){
h[j] -= 1ll*zt*(h[j/prime[i]]-h[prime[i]-])%mod;
h[j]<?h[j]+=mod:;
}
}
} int f(int now,int zt){
if(now <= sqr){
if(now <= || prime[zt] >= now) return ;
int ans = fh[now]-fh[prime[zt]]; if(ans < ) ans += mod;
if(1ll*prime[zt]*prime[zt] >= now) return ans;
for(int i=zt+;1ll*prime[i]*prime[i] <= now;i++){
long long dt = prime[i],k=;
while(dt <= now){
ans += ((prime[i]^k)*(f(now/dt,i)+))%mod; ans %= mod;
dt = dt*prime[i];k++;
}
ans-=(prime[i]^); if(ans < ) ans += mod;
}
return ans;
}else{
long long rd = n/(*sqr-now+);
int ans = fh[now]-fh[prime[zt]]; if(ans < ) ans += mod;
if(1ll*prime[zt]*prime[zt] >= rd) return ans;
for(int i=zt+;1ll*prime[i]*prime[i] <= rd&&i<=num;i++){
long long dt = prime[i],k=;
while(dt <= rd){
if(rd/dt <= sqr) ans += ((prime[i]^k)*(f(rd/dt,i)+))%mod;
else ans += ((prime[i]^k)*(f(*sqr-(*sqr-now+)*dt+,i)+))%mod;
ans %= mod; dt = dt*prime[i];k++;
}
ans-=(prime[i]^); if(ans < ) ans += mod;
}
return ans;
}
} void work(){
for(int i=;i<=sqr;i++){h[i] = (h[i-] + i)%mod;}
for(int i=;i<=sqr;i++){
long long fak=(n/i)%mod;h[sqr*-i+]=((+fak)*(fak-)/)%mod;
}
geth(); for(int i=;i<=*sqr;i++) fh[i] += h[i];
for(int i=;i<=sqr;i++) h[i] = i-; for(int i=;i<=sqr;i++) h[*sqr-i+]=(n/i-)%mod;
geth(); for(int i=;i<=*sqr;i++) fh[i] = fh[i]-h[i]+,fh[i] =(fh[i]+mod)%mod;
printf("%d",f(*sqr,)+);
} int main(){
scanf("%lld",&n); sqr = sqrt(n);
if(n == ){puts("");return ;}
getprime(sqr);work();
return ;
}
