题意:给定\(n\le 10^9\),求:\(F(n)=\sum_{i=1}^n\sum_{j=1}^i\frac{\mathrm{lcm}(i,j)}{\mathrm{gcd}(i,j)}\),对1e9+7取模
推式子:
\(F(n)=\sum_{i=1}^n\sum_{j=1}^i\frac{\mathrm{lcm}(i,j)}{\mathrm{gcd}(i,j)}\)
\(=\sum_{i=1}^n\sum_{j=1}^i\frac{ij}{\gcd^2(i,j)}\)
\(=\sum_{p=1}^n\frac1{p^2}\sum_{i=1}^n\sum_{j=1}^iij[\gcd(i,j)=p]\)
\(=\sum_{p=1}^n\sum_{i=1}^{n/p}\sum_{j=1}^{i/p}ij[\gcd(i,j)=1]\)
\(=\sum_{p=1}^n\sum_{i=1}^{n/p}i\sum_{j=1}^{i/p}j[\gcd(i,j)=1]\)
根据一个经典式子:\(\sum_{i=1}^ni[\gcd(i,n)=1]=\frac {[n=1]+n\varphi(n)}{2}\)
\(=\sum_{p=1}^n\sum_{i=1}^{n/p}i\frac{[i=1]+i\varphi(i)}{2}\)
\(=\frac {n+\sum_{i=1}^n\lfloor\frac n i\rfloor i^2\varphi(i)}2\)
现在我们考虑求\(\sum_{i=1}^n\lfloor\frac n i\rfloor i^2\varphi(i)\),暴力线性筛是\(O(n)\)的,显然会T,考虑用杜教筛优化
设\(f(i)=i^2\varphi(i),S(n)=\sum_{i=1}^nf(i)\)
原式=\(\sum_{i=1}^n\lfloor\frac n i\rfloor f(i)\)
由于\(\lfloor\frac ni\rfloor\)的取值只有\(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\n", sum);
return 0;
}
Min_25筛表示他需要重复计算根N次,就GG了
Min_25筛可以做,不过设S时候要用递推法,而不是递归计算,这里懒得写了