[Luogu 4176] Lucas的数论

BZOJ传送门

Description

去年的Lucas非常喜欢数论题，但是一年以后的Lucas却不那么喜欢了。

在整理以前的试题时，发现了这样一道题目“求$\sum f(i)$,其中$1\le i\le N$”，其中$f(i)$表示$i$的约数个数。他现在长大了，题目也变难了。

求如下表达式的值：

其中$f(ij)$表示$ij$的约数个数。

他发现答案有点大，只需要输出模$1000000007$的值。

Input

第一行一个整数$n$。

Output

一行一个整数$ans$，表示答案模$1000000007$的值。

Sample Input

2

Sample Output

8

HINT

对于100%的数据$n\le 10^9$。

解题分析

首先有个妙妙的结论：
$${\sigma }_0(n\ast m)=\sum _{i|n}\sum _{j|m}[gcd(i,j)=1]$$
设$n=p_1^{a_1}\ast p_2^{a_2}\ast p_3^{a_3}\ast ......\ast p_k^{a_k}$, $m=p_1^{b_1}\ast p_2^{b_2}\ast p_3^{b_3}\ast ......\ast p_k^{b_k}$，那么如果$gcd(i,j)=1$， 对于一个质数$p_1$来说， 要么是$i$里含有$p_1$， 共有$a_1$种情况， 要么是$j$里含有$p_1$，共有$b_1$种情况， 要么是两个都没有， 共一种情况， 所以加起来一共是$a_1+b_1+1$种， 刚好符合${\sigma }_0$的性质。

大力推式子：
$$\begin{gathered} \sum _{i=1}^n\sum _{j=1}^n{\sigma }_0(i,j) \\ =\sum _{i=1}^n\sum _{j=1}^n\sum _{p|i}\sum _{d|j}[gcd(i,j)=1] \\ =\sum _{q=1}^n\mu (q)\sum _{q|p}^n\sum _{q|d}^n\lfloor \frac{n}{d}\rfloor \lfloor \frac{n}{p}\rfloor \\ =\sum _{q=1}^n\mu (q)\sum _{i=1}^{\lfloor \frac{n}{q}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{q}\rfloor }\lfloor \frac{n}{qi}\rfloor \lfloor \frac{n}{qj}\rfloor \end{gathered}$$
设$g(x)={\sum }_{i=1}^x\lfloor \frac{x}{i}\rfloor$，然后发现$g(x)={\sigma }_0(x)$ ，就可以把上面这玩意化成：
$$\sum _{q=1}^n\mu (q){\sigma }_0^2(\lfloor \frac{n}{q}\rfloor )$$
${\sum }_{i=1}^n\mu (i)$我们已经可以用杜教筛筛出来了， 现在考虑后面那个玩意。由于我们知道$I(n)\ast I(n)={\sigma }_0(n)$， 那么莫比乌斯反演一下即可得到${\sigma }_0(n)\ast \mu (n)=I(n)$。 这样就可以科学筛出${\sigma }_0$前缀和了。

代码如下：

#include <cstdio>
#include <cmath>
#include <cstring>
#include <cctype>
#include <cstdlib>
#include <algorithm>
#include <tr1/unordered_map>
#define R register
#define IN inline
#define W while
#define gc getchar()
#define ll long long
#define MOD 1000000007ll
#define MX 5000050
std::tr1::unordered_map <int, int> mp1;
std::tr1::unordered_map <int, ll> mp2;
int pcnt, n;
int pri[MX], miu[MX], sig0[MX], sum1[MX], cnt[MX];
bool npr[MX];
ll sum2[MX];
void get()
{
	miu[1] = 1, sig0[1] = 1;
	R int i, j, tar;
	for (i = 2; i <= 5e6; ++i)
	{
		if (!npr[i]) miu[i] = -1, sig0[i] = 2, cnt[i] = 1, pri[++pcnt] = i;
		for (j = 1; j <= pcnt; ++j)
		{
			tar = i * pri[j];
			if (tar > 5e6) break;
			npr[tar] = true;
			if (!(i % pri[j]))
			{
				miu[tar] = 0;
				sig0[tar] = sig0[i] / (cnt[i] + 1) * (cnt[i] + 2);
				cnt[tar] = cnt[i] + 1;
				break;
			}
			miu[tar] = -miu[i];
			sig0[tar] = 2 * sig0[i];
			cnt[tar] = 1;
		}
	}
	for (R int i = 1; i <= 5e6; ++i)
	sum1[i] = sum1[i - 1] + miu[i], sum2[i] = (sum2[i - 1] + sig0[i]) % MOD;
}
ll solve1(R int pos)//for miu
{
	if (pos <= 5e6) return sum1[pos];
	if (mp1.count(pos)) return mp1[pos];
	ll ret = 1;
	for (R int lef = 2, rig; lef <= pos; lef = rig + 1)
	{
		rig = pos / (pos / lef);
		ret -= (rig - lef + 1) * solve1(pos / lef) % MOD;
		ret %= MOD;
	}
	return mp1[pos] = ret;
}
ll solve2(R int pos)
{
	if (pos <= 5e6) return sum2[pos];
	if (mp2.count(pos)) return mp2[pos];
	ll ret = pos;
	for (R int lef = 2, rig; lef <= pos; lef = rig + 1)
	{
		rig = pos / (pos / lef);
		ret -= (solve1(rig) - solve1(lef - 1)) * solve2(pos / lef) % MOD;
		ret %= MOD;
	}
	return mp2[pos] = ret;
}
int main(void)
{
	ll ans = 0;
	scanf("%d", &n); get();
	for (R int lef = 1, rig; lef <= n; lef = rig + 1)
	{
		rig = n / (n / lef);
		ans = (ans + (solve1(rig) - solve1(lef - 1)) * solve2(n / lef) % MOD * solve2(n / lef) % MOD) % MOD;
	}
	printf("%lld", (ans + MOD) % MOD);
}