loj#6229. 这是一道简单的数学题 (??反演+杜教筛)

题目链接

题意：给定(nle 10^9)，求：(F(n)=sum_{i=1}^nsum_{j=1}^ifrac{mathrm{lcm}(i,j)}{mathrm{gcd}(i,j)})，对1e9+7取模

推式子：

(F(n)=sum_{i=1}^nsum_{j=1}^ifrac{mathrm{lcm}(i,j)}{mathrm{gcd}(i,j)})

(=sum_{i=1}^nsum_{j=1}^ifrac{ij}{gcd^2(i,j)})

(=sum_{p=1}^nfrac1{p^2}sum_{i=1}^nsum_{j=1}^iij[gcd(i,j)=p])

(=sum_{p=1}^nsum_{i=1}^{n/p}sum_{j=1}^{i/p}ij[gcd(i,j)=1])

(=sum_{p=1}^nsum_{i=1}^{n/p}isum_{j=1}^{i/p}j[gcd(i,j)=1])

根据一个经典式子：(sum_{i=1}^ni[gcd(i,n)=1]=frac {[n=1]+nvarphi(n)}{2})

(=sum_{p=1}^nsum_{i=1}^{n/p}ifrac{[i=1]+ivarphi(i)}{2})

(=frac {n+sum_{i=1}^nlfloorfrac n i floor i^2varphi(i)}2)

现在我们考虑求(sum_{i=1}^nlfloorfrac n i floor i^2varphi(i))，暴力线性筛是(O(n))的，显然会T，考虑用杜教筛优化

设(f(i)=i^2varphi(i),S(n)=sum_{i=1}^nf(i))

原式=(sum_{i=1}^nlfloorfrac n i floor f(i))

由于(lfloorfrac ni floor)的取值只有(O(sqrt n))种，所以可以计算每一块比下一块多多少，然后计算下前缀和

这块可能不太好理解，例如n=11时

ans=11f(1)+5f(2)+3f(3)+2(f(4)+f(5))+1(f(6)+...+f(11))，替换成S则有

ans=(11-5)S(1)+(5-3)S(2)+(3-2)S(3)+(2-1)S(5)+S(11)

直接上杜教筛就行了，因为这些取值在计算S(11)时候都要用到，并且每个n/x只有一种，开个数组记忆化即可

代码里枚举的顺序不太一样，是从S大到S小一样，两块的差就相当于某一块的大小

#include <cstdio>
using namespace std;

const int p = 1000000007, fuck = 1000000;

bool vis[fuck + 233];
int prime[fuck], tot;
int phi[fuck + 233];

int qpow(int x, int y)
{
	int res = 1;
	while (y > 0)
	{
		if (y & 1) res = res * (long long)x % p;
		x = x * (long long)x % p;
		y >>= 1;
	}
	return res;
}

int inv2 = qpow(2, p - 2), inv6 = qpow(6, p - 2);

int s1(int x) { x %= p; return x* (long long)(x + 1) % p * inv2 % p; }
int s2(int x) { x %= p; return x * (long long)(x + 1) % p * (x * 2 + 1) % p * inv6 % p; }
int s3(int x) { return qpow(s1(x), 2); }

bool count[4000010];
int ans[400010], n;

int chuans(int x)
{
	if (x <= fuck) { return phi[x]; }
	if (count[n / x]) return ans[n / x];
	count[n / x] = 1;
	int res = s3(x);
	for (int i = 2, j; i <= x; i = j + 1)
	{
		j = x / (x / i);
		res = ((res - (s2(j) - s2(i - 1)) * (long long)chuans(x / i) % p) % p + p) % p;
	}
	return ans[n / x] = res;
}

int main()
{
	phi[1] = 1;
	for (int i = 2; i <= fuck; i++)
	{
		if (vis[i] == false) prime[++tot] = i, phi[i] = i - 1;
		for (int j = 1; j <= tot && i * prime[j] <= fuck; j++)
		{
			vis[i * prime[j]] = true;
			if (i % prime[j] == 0) { phi[i * prime[j]] = phi[i] * prime[j]; break; }
			else phi[i * prime[j]] = phi[i] * (prime[j] - 1);
		}
		phi[i] = (phi[i - 1] + phi[i] * (long long)i % p * i % p) % p;
	}
	scanf("%d", &n);
	int sum = 0;
	for (int i = 1, j; i <= n; i = j + 1)
	{
		j = n / (n / i);
		sum = (sum + (j - i + 1) * (long long)chuans(n / i) % p) % p;
	}
	sum = (sum + n) % p * (long long)inv2 % p;
	printf("%d
", sum);
	return 0;
}

Min_25筛表示他需要重复计算根N次，就GG了

Min_25筛可以做，不过设S时候要用递推法，而不是递归计算，这里懒得写了