【ybtoj高效进阶 21177】小小网格(杜教筛)(数论分块)(莫比乌斯反演)

小小网格

题目链接:ybtoj高效进阶 21177

题目大意

给你求 ∑i=1~n∑j=1~mφ(gcd(i,j))。

思路

重新写好看点:
i=1nj=1mφ(gcd(i,j))
dφ(d)i=1nj=1m[gcd(i,j)=d]
dφ(d)i=1ndj=1md[gcd(i,j)=1]
dφ(d)i=1ndj=1mdp|gcd(i,j)μ(p)
dφ(d)pi=1,i|pndj=1,j|pmdμ(p)
dφ(d)pμ(p)ndpmdp

然后设 T=pd,然后有:
dφ(d)pμ(p)nTmT
T=1min(n,m)nTmT(d|Tφ(d)μ(Td))

然后你会发现右边就是一个狄利克雷卷积。
那我们设 g=φμ(积性函数,因为是两个积性函数相乘),然后式子就是:
T=1min(n,m)nTmTg(d)

然后你发现左边两个可以数论分块,然后你考虑右边的能不能求前缀和。
看到是积性函数,考虑用杜教筛,设 S(n)=i=1ng(i)

然后你就要找一个合适的 f,不难看出可以用 II=d
然后 gd=φμII=(φI)(μI)=idϵ=id

然后就有 d(1)S(n)=i=1nid(i)i=2nd(i)S(ni)=n(n+1)2i=2nd(i)S(ni)

然后 d(i) 的前缀和可以数论分块,小于 n23 的预处理前缀和,大于的直接根号求。

然后小小卡个常即可。

代码

#include<map> #include<cstdio> #define ll long long #define mo 1000000007 using namespace std; const int Maxn = 2000001; int n, m; ll d[Maxn + 1], g[Maxn + 1], ans; int low[Maxn + 1], phi[Maxn + 1], prime[Maxn + 1]; map <int, ll> ans_g; void init() {//杜教筛的预处理 g[1] = 1; phi[1] = 1; low[1] = 1; for (int i = 2; i <= Maxn; i++) { if (!low[i]) phi[i] = i - 1, low[i] = i, prime[++prime[0]] = i, g[i] = i - 2; for (int j = 1; j <= prime[0] && i * prime[j] <= Maxn; j++) { if (i % prime[j]) low[i * prime[j]] = prime[j], phi[i * prime[j]] = phi[i] * (prime[j] - 1), g[i * prime[j]] = g[i] * g[prime[j]]; else { low[i * prime[j]] = low[i] * prime[j], phi[i * prime[j]] = phi[i] * prime[j]; if (low[i] != i) g[i * prime[j]] = g[i / low[i]] * g[low[i] * prime[j]]; else g[i * prime[j]] = phi[i * prime[j]] + phi[i] * (-1); break; } } } for (int i = 1; i <= Maxn; i++) for (int j = 1; i * j <= Maxn; j++) d[i * j]++; for (int i = 1; i <= Maxn; i++) d[i] += d[i - 1], g[i] += g[i - 1]; } ll sum_d(ll n) {//数论分块 if (n <= Maxn) return d[n]; ll re = 0; for (ll l = 1, r; l <= n; l = r + 1) { r = n / (n / l); re += (r - l + 1) * (n / l); } return re; } ll sum_g(ll n) {//杜教筛 if (n <= Maxn) return g[n]; if (ans_g[n]) return ans_g[n]; ll re = n * (n + 1) / 2; for (ll l = 2, r; l <= n; l = r + 1) { r = n / (n / l); re -= (sum_d(r) - sum_d(l - 1)) * sum_g(n / l); } return ans_g[n] = re; } int main() { // freopen("mesh.in", "r", stdin); // freopen("mesh.out", "w", stdout); scanf("%d %d", &n, &m); init(); for (int l = 1, r; l <= n && l <= m; l = r + 1) {//数论分块 r = min(n / (n / l), m / (m / l)); ans = (ans + 1ll * (n / l) * (m / l) * (sum_g(r) - sum_g(l - 1))) % mo; } printf("%lld", ans); return 0; }

__EOF__

本文作者あおいSakura
本文链接https://www.cnblogs.com/Sakura-TJH/p/YBT_GXJJ_21177.html
关于博主:评论和私信会在第一时间回复。或者直接私信我。
版权声明:本博客所有文章除特别声明外,均采用 BY-NC-SA 许可协议。转载请注明出处!
声援博主:如果您觉得文章对您有帮助,可以点击文章右下角推荐一下。您的鼓励是博主的最大动力!
posted @   あおいSakura  阅读(31)  评论(0编辑  收藏  举报
编辑推荐:
· go语言实现终端里的倒计时
· 如何编写易于单元测试的代码
· 10年+ .NET Coder 心语,封装的思维:从隐藏、稳定开始理解其本质意义
· .NET Core 中如何实现缓存的预热?
· 从 HTTP 原因短语缺失研究 HTTP/2 和 HTTP/3 的设计差异
阅读排行:
· 分享一个免费、快速、无限量使用的满血 DeepSeek R1 模型,支持深度思考和联网搜索!
· 基于 Docker 搭建 FRP 内网穿透开源项目(很简单哒)
· ollama系列01:轻松3步本地部署deepseek,普通电脑可用
· 25岁的心里话
· 按钮权限的设计及实现
点击右上角即可分享
微信分享提示