Min_25筛

参考

Min_25 筛学习笔记 By zhoukangyang

Min25筛小记 By command_block

---

题意：给定一积性函数 \(f\)，其在 \(p^c\) 处的取值是一个关于 \(p^c\) 的低次多项式，求 \(\sum_{i=1}^n f(i)\)。

\(1\le n\le 10^{10}\)

---

先处理质数处的取值。

拆质数处的多项式，依次处理每个单项式。

考虑把和分为两部分，质数和非质数处的取值。

考虑一积性函数 \(g\)，其在质数处的取值和 \(f\) 对应的单项式相同。

那么 \(\sum_{i=1}^n g(i)\) 是容易求的，我们考虑去掉其中的合数，可以用每个合数的最小值因子筛掉他。

设 \(h[i,j]\) 表示 \(1...i\) 中筛去了最小质因子 \(\in\) 前 \(j\) 个质数的合数后，剩下数取值的和。

记 \(p_i\) 为第 \(i\) 个质数。考虑从 \(h[i,j-1]\) 转移，我们需要去掉最小质因子为 \(p_j\) 的数。这些数除以 \(p_j\) 后，其最小值因子都 \(\ge p_j\)。

不难发现，最小值因子为 \(p_j\) 的数的取值总和为 \(g(p_j)\cdot(h[\lfloor \frac n{p_j}\rfloor,j-1]-\sum_{k=1}^{j-1}g(k))\)

那么

\[h[i,j]=h[i,j-1]-g(p_j)\cdot (h[\lfloor \frac n{p_j} \rfloor,j-1]-\sum_{k=1}^{j-1}g(k)-g(1)) \]

注意一点，就是 \(h[i,j]\) 始终保留 \(g(1)\) 的值，因为 \(1\) 永远不会被筛掉。

再考虑合数，下文默认 \(h[i,j]\) 为已考虑的所有单项式的和。

一种想法是通过枚举他的最小值因子来求。

设 \(S[i,j]\) 表示 \(2...i\) 中筛去最小值因子 \(\in\) 前 \(j\) 个质数后剩下的数的 \(f\) 取值之和。

首先算上质数，贡献为 \(h[i,m]\)（\(m\) 为 \(\sqrt n\) 内质数个数）。枚举最小值因子为第 \(k\) 个质数，其次数为 \(c\)，贡献为 \(f(p_k^c)\cdot (S[\lfloor \frac n{p_k^c} \rfloor,k]+[c>1])\)，其中 \(c>1\) 是因为此时 \(p_k^c\) 自身也是一个合数。

\[S[i,j]=h[i,m]+f(p_k^c)\cdot (S[\lfloor \frac n{p_k^c} \rfloor,k]+[c>1]) \]

Q: 为什么质数要单独分开？

A: 因为质数值域很大，而合数的最小值因子 \(\le 10^5\)

时间复杂度不会证，\(O(\frac {n^{\frac 34}}{\log n})\)。

答案是 \(S[n,0]-\sum g(1)+1\)，注意 \(S[n,0]\) 算上的是 \(1\) 在 \(g\) 处的单项式总和，而不是 \(f(1)\)。

#include<bits/stdc++.h>
#define pr pair<ll,ll>
#define x first
#define y second
#define mkp(a,b) make_pair(a,b)
#define pb push_back
#define ll long long
using namespace std;
const ll maxn=2e6+10,mod=1e9+7;
ll n,sq,pri[maxn],m,id1[maxn],id2[maxn],sum1[maxn],sum2[maxn],w[maxn],len,g1[maxn],g2[maxn];
bool b[maxn];
void xxs()
{
	sum1[0]=sum2[0]=1;
	for(ll i=2;i<=sq;i++)
	{
		if(!b[i])
		{
			pri[++m]=i;
			sum1[m]=sum1[m-1]; sum2[m]=sum2[m-1];
			sum1[m]+=i, sum2[m]+=i*i;
			sum1[m]%=mod; sum2[m]%=mod;
		}
		for(ll j=1;j<=m&&i*pri[j]<=sq;j++)
		{
			ll k=i*pri[j]; b[k]=true;
			if(i%pri[j]==0) break;
		}
	}
}
ll id(ll x)
{
	if(x<=sq) return id1[x];
	return id2[n/x];
}
ll Sum1(ll n)
{
	n%=mod;
	return n*(n+1)%mod*(mod/2+1)%mod;
}
ll Sum2(ll n)
{
	n%=mod;
	return n*(n+1)%mod*(2*n+1)%mod*166666668%mod;
}
ll F(ll n,ll k)
{
	if(pri[k]>n) return 0;
	ll ans=((g2[id(n)]-g1[id(n)]+mod)%mod-(sum2[k]-sum1[k]+mod)%mod+mod)%mod;
	for(ll i=k+1;i<=m&&pri[i]*pri[i]<=n;i++)
		for(ll x=1,y=pri[i];y<=n;y*=pri[i],++x)
		{
			ans=(ans+y%mod*(y%mod+mod-1)%mod*(F(n/y,i)+(x>1)))%mod;
		}
	return ans;
}
int main()
{
	scanf("%lld",&n); sq=sqrt(n);
	xxs();
	for(ll i=1;i<=n;i++)
	{
		ll d=n/i, r=n/d; w[++len]=d;
		g1[len]=Sum1(d); g2[len]=Sum2(d);
		if(d<=sq) id1[d]=len;
		else id2[n/d]=len;
		i=r;
	}
	for(ll i=1;i<=m;i++)
	{
		for(ll j=1;j<=len&&pri[i]*pri[i]<=w[j];j++)
		{
			g1[j]=(g1[j]-pri[i]*(g1[id(w[j]/pri[i])]-sum1[i-1])%mod+mod)%mod;
			g2[j]=(g2[j]-pri[i]*pri[i]%mod*(g2[id(w[j]/pri[i])]-sum2[i-1])%mod+mod)%mod;
		} 
	}
	printf("%lld",(F(n,0)+1)%mod);
	return 0;
}