51Nod1220 约数之和 【杜教筛】

题目描述：

求$$\sum _{i=1}^n\sum _{j=1}^n{\sigma }_1(i\ast j)$$
$n\le 10^9，mod10^9+7$

题目分析：

有个结论：$${\sigma }_k(i\ast j)=\sum _{x|i}\sum _{y|j}[(x,y)==1]x^k\ast j^k/y^k$$
记${\sigma }_k(n)$为$n$的约数的$k$次幂之和，质因数分解$n$，$n=p_1^{a_1}{p}_2^{a_2}...p_t^{a_t}$
$$\begin{gathered} {\sigma }_k(n)=(1+p_1^{1k}+p_1^{2k}+\cdots +p_1^{a_{1}k}) \\ \ast (1+p_2^{1k}+p_2^{2k}+\cdots +p_2^{a_{2}k}) \\ \dots \\ \ast (1+p_t^{1k}+p_t^{2k}+\cdots +p_t^{a_{t}k}) \end{gathered}$$
考虑一个质因子$p$，它对${\sigma }_k(n)$的贡献为$1+p^k+p^{2k}+...+p^{ak}$
那么看${\sigma }_k(i\ast j)$，考虑它们的公共质因子$p$，在i中有a次，在j中有b次，造成的贡献就可以表示为$$\sum _{x=0}^a\sum _{y=0}^b[(p^x,p^y)==1]p^{xk}\ast p^{bk}/p^{yk}$$
根据乘法原理，有$$\begin{gathered} {\sigma }_k(i\ast j)=\sum _{x_1=0}^{a_1}\sum _{y_1=0}^{b_1}[(p^{x_1},p^{y_1})==1]p^{x_{1}k}\ast p^{b_{1}k}/p^{y_{1}k} \\ \ast \sum _{x_2=0}^{a_2}\sum _{y_2=0}^{b_2}[(p^{x_2},p^{y_2})==1]p^{x_{2}k}\ast p^{b_{2}k}/p^{y_{2}k} \\ \dots \\ \ast \sum _{x_t=0}^{a_t}\sum _{y_t=0}^{b_t}[(p^{x_t},p^{y_t})==1]p^{x_{t}k}\ast p^{b_{t}k}/p^{y_{t}k} \end{gathered}$$
只有当所有$(p^x,p^y)==1$都满足时才会对答案造成贡献，所以可变形为：
$${\sigma }_k(i\ast j)=\sum _{x|i}\sum _{y|j}[(x,y)==1]x^k\ast j^k/y^k$$

ok,正文开始

$$\begin{array}{rl}& \sum _{i=1}^n\sum _{j=1}^n{\sigma }_1(i\ast j)\\ & =\sum _{i=1}^n\sum _{j=1}^n\sum _{x|i}\sum _{y|j}[(x,y)==1]x\ast j/y\\ & =\sum _{k=1}^n\mu (k)\sum _{x=1}^{\frac{n}{k}}\lfloor \frac{n}{kx}\rfloor \sum _{y=1}^{\frac{n}{k}}\sum _{j=1}^{\frac{n}{ky}}xk\ast kyj/ky\\ & =\sum _{k=1}^n\mu (k)k\sum _{x=1}^{\frac{n}{k}}\lfloor \frac{n}{kx}\rfloor x\sum _{j=1}^{\frac{n}{k}}j\lfloor \frac{n}{kj}\rfloor \\ & =\sum _{k=1}^n\mu (k)k{\left(\sum _{x=1}^{\frac{n}{k}}\lfloor \frac{n}{kx}\rfloor x\right)}^2\\ & =\sum _{k=1}^n\mu (k)k{\left(\sum _{i=1}^{\frac{n}{k}}{\sigma }_1(i)\right)}^2\end{array}$$

线性筛$10^6$，$\mu (k)k$标准的杜教筛，${\sigma }_1(n)的前缀和可以\sqrt{n}$分块暴算
时间复杂度$O(n^{\frac{2}{3}})$

• 其实有种更简单清晰的推法，第一步先不要把x,y提到最前面
  $$\begin{array}{rl}& =\sum _{i=1}^n\sum _{j=1}^n\sum _{x|i}\sum _{y|j}[(x,y)==1]x\ast j/y\\ & =\sum _{k=1}^n\mu (k)\sum _{i=1}^{\frac{n}{k}}\sum _{j=1}^{\frac{n}{k}}\sum _{x|i}\sum _{y|i}xk\ast jk/yk\\ & =\sum _{k=1}^n\mu (k)k\sum _{i=1}^{\frac{n}{k}}\sum _{j=1}^{\frac{n}{k}}\sum _{x|i}x\sum _{y|j}\frac{j}{y}\\ & =\sum _{k=1}^n\mu (k)k\sum _{i=1}^{\frac{n}{k}}\sum _{j=1}^{\frac{n}{k}}{\sigma }_1(i)\ast {\sigma }_1(j)\\ & =\sum _{k=1}^n\mu (k)k{\left(\sum _{i=1}^{\frac{n}{k}}{\sigma }_1(i)\right)}^2\end{array}$$

#include<cstdio>
#include<map>
using namespace std;
const int N = 1000000, mod = 1e9+7, inv2 = 5e8+4;
int p[N/10],sm[N+5],sd[N+5],mp[N+5];
bool v[N+5];
map<int,int>F;
void Prime()
{
	sm[1]=sd[1]=1;int cnt=0;
	for(int i=2;i<=N;i++)
	{
		if(!v[i]) p[++cnt]=i,sm[i]=-i,sd[i]=mp[i]=i+1;
		for(int j=1,k;j<=cnt&&(k=p[j]*i)<=N;j++)
		{
			v[k]=1;
			if(i%p[j]==0) {sm[k]=0,mp[k]=mp[i]*p[j]+1,sd[k]=sd[i]/mp[i]*mp[k];break;}
			sm[k]=sm[i]*sm[p[j]];
			sd[k]=sd[i]*sd[p[j]];
			mp[k]=p[j]+1;
		}
	}
	for(int i=1;i<=N;i++) sm[i]=(sm[i]+sm[i-1])%mod,sd[i]=(sd[i]+sd[i-1])%mod;
}
int Summu(int n)
{
	if(n<=N) return sm[n];
	if(F.count(n)) return F[n];
	int ret=1;
	for(int i=2,j;i<=n;i=j+1)
	{
		j=n/(n/i);
		ret=(ret-1ll*(i+j)*(j-i+1)%mod*inv2%mod*Summu(n/i))%mod;
	}
	return F[n]=ret;
}
int calc(int n)
{
	if(n<=N) return sd[n];
	int ret=0;
	for(int i=1,j;i<=n;i=j+1)
	{
		j=n/(n/i);
		ret=(ret-1ll*(i+j)*(j-i+1)%mod*inv2%mod*(n/i))%mod;
	}
	return ret;
}
int main()
{
	Prime();
	int n,ans=0;
	scanf("%d",&n);
	for(int i=1,j,tmp,pre=0,tt;i<=n;i=j+1)
	{
		j=n/(n/i);
		ans=(ans+1ll*((tmp=Summu(j))-pre)*(tt=calc(n/i))%mod*tt)%mod;
		pre=tmp;
	}
	printf("%d",(ans+mod)%mod);
}