LOJ#6682. 梦中的数论（min_25筛）

第一次打\(min_25\)筛，所以写一写

做法

就是穿了件肥肥的衣服~，等价于这个：

\[\frac{1}{2}(\sum\limits_{i=1}^n \sigma^2(i)-\sum\limits_{i=1}^n \sigma(i)) \]

后面那个分块一下就好了（不会有人学傻了去筛吧。。。）

看前面这个\(\sigma^2(i)\)，设\(F(i)=\sigma^2(i)\)

显然：\(F(i)\)是积性函数且\(F(p^k)=(k+1)^2,p\in prime\)

我们考虑min_25筛的公式

\(g(n,j)=\sum\limits_{i=1}^n[i\in prime || min_i(P)>P_j]F(i)\)：$$\begin{aligned}\
g(n,j)=
\begin{cases}
g(n,j-1)&P_j^2>n\
g(n,j-1)-F(P_j)[g(\frac{n}{P_j},j-1)-\sum\limits_{i=1}^{{j-1}F(P_i)]&P_j}2≤n\
\end{cases}
\end{aligned}$$

\(s(n,j)=\sum\limits_{i=1}^n[min_i(P)≥P_j]F_i\)：$$s(n,j)=g(n,|P|)-\sum\limits_{i=1}^{{j-1}F(P_i)+\sum\limits_{k=j}}{P_k^{2≤n}\sum\limits_{e=1}}{P_k^{{e+1}≤n}s(\frac{n}{P_k}e},k+1)F(P_k^e)+F(P_k)$$

\(Ans=s(n,1)+F(1)\)

• \(s(x,y)\)可以写个递归
• 而在\(s(x,y)\)中，我们只会用到\(g(n,|p|)\)，那么第二维相当于是废的，\(g\)的第二维是递推得出的
  考虑放在一个循环里做，\(g\)的本质就是埃氏筛法，类似做即可
• \(F(p)=4,p\in prime\)，则\(g(n,|p|)=4\times(\sum\limits_{i=1}^n[i\in prime])\)
  设\(g'(n,|p|)=\sum\limits_{i=1}^n[i\in prime]\)，则\(g(n,|p|)=4\times g'(n,|p|)\)，直接做\(g'\)即可
• 显然\(g'\)的第一维的个数是\(O(\sqrt n)\)级别的，那数组映射下来即可

Code

为了增加可读性，没有卡常

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const LL maxn=1e6+9,mod=998244353;
#define pb push_back
#define opt operator
LL Read(){
	LL x(0),f(1); char c=getchar();
	while(c<'0'|| c>'9'){
		if(c=='-') f=-1; c=getchar();
	}
	while(c>='0' && c<='9'){
		x=(x<<3ll)+(x<<1ll)+c-'0'; c=getchar();
	}
	return x*f;
}
LL n,m,T,ans1,ans2,tot;
LL w[maxn],pos1[maxn],pos2[maxn],f[maxn],pri[maxn],vis[maxn];
LL Pow(LL base,LL b){ LL ret(1); while(b){ if(b&1) ret=ret*base%mod; base=base*base%mod; b>>=1; } return ret; }
LL Id(LL x){ return x<maxn?pos1[x]:pos2[n/x]; }
void Fir(LL N){
	for(LL i=2;i<=N;++i){
		if(!vis[i]) pri[++tot]=i;
		for(LL j=1;j<=tot && pri[j]*i<=N;++j){
			vis[pri[j]*i]=1; if(i%pri[j]==0) break;
		}
	}
}
LL S(LL x,LL y){
	if(x<=1 || pri[y]>x) return 0;
	LL ret(0);
	ret=(f[Id(x)]-(y-1))%mod*4%mod;
	for(LL k=y;k<=tot && pri[k]*pri[k]<=x; ++k){
		LL t1(pri[k]),t2(t1*t1);
		for(LL e=1;t2<=x;++e,t1=t2,t2*=pri[k]){
			ret=(ret+S(x/t1,k+1)*(e+1)%mod*(e+1)%mod+(e+2)*(e+2)%mod)%mod;
		}
	}
	return ret;
}
int main(){
    n=Read();
    for(LL i=1,j;i<=n;i=j+1){
    	j=n/(n/i);
    	w[++w[0]]=n/i;
    	if(w[w[0]]<maxn) pos1[w[w[0]]]=w[0];
    	else pos2[j]=w[0];
    	f[w[0]]=(w[w[0]]-1+mod)%mod;
    	ans2=(ans2+(j-i+1)*w[w[0]]%mod)%mod;
	}
    Fir(sqrt(n));
    for(LL i=1;i<=tot;++i){
    	for(LL j=1;j<=w[0] && pri[i]*pri[i]<=w[j];++j){
    		f[j]=(f[j]-(f[Id(w[j]/pri[i])]-(i-1))+mod)%mod;
		}
	}
	ans1=(S(n,1)+1)%mod;
	printf("%lld\n",(ans1-ans2+mod)%mod*Pow(2,mod-2)%mod);
	return 0;
}