洛谷 P1829 [国家集训队]Crash的数字表格 / JZPTAB

最新推荐文章于 2024-05-02 13:48:13 发布

原创最新推荐文章于 2024-05-02 13:48:13 发布 · 210 阅读

1 ·

CC 4.0 BY-SA版权

文章标签：

#莫比乌斯反演

数论专栏收录该内容

12 篇文章

订阅专栏

本文介绍了一种使用数论分块技巧解决求和∑i=1n∑j=1mlcm(i,j)(mod20101009)的问题，通过将原问题转换为求gcd(i,j)=1时的i*j之和，再利用莫比乌斯函数进行化简，最终实现了高效求解。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题意

求

$20101009)\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\text{lcm}(i,j)(\bmod 20101009)$

思路

容易想到原式等价于

$∑i=1n∑j=1mi∗jgcd⁡(i,j)\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\frac{i* j}{\gcd(i,j)}$

枚举 $i, j$ 的最大公约数 $d$ ，显然 $gcd⁡(id,jd)=1\gcd(\frac id,\frac jd)=1$ ，即 $id\frac id$ 和 $jd\frac jd$ 互质

$∑i=1n∑j=1m∑d∣i,d∣j,gcd⁡(id,jd)=1i∗jd\sum\limits_{i=1}^{n}\sum\limits_{j=1}^m\sum\limits_{d|i,d|j,\gcd(\frac id,\frac jd)=1}\frac{i*j}d$

变换求和顺序

$∑d=1nd∑i=1⌊nd⌋∑j=1⌊md⌋[gcd⁡(i,j)=1]i∗j\sum\limits_{d=1}^{n}d\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}[\gcd(i,j)=1]i*j$

记 $sum(n,m)=∑i=1n∑j=1m[gcd⁡(i,j)=1]i∗jsum(n,m)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j)=1]i*j$

对其进行化简，用 $ε(gcd⁡(i,j))\varepsilon(\gcd(i,j))$ 替换 $[gcd⁡(i,j)=1][\gcd(i,j)=1]$

$∑i=1n∑j=1m∑d∣gcd⁡(i,j)μ(d)∗i∗j\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum\limits_{d|\gcd(i,j)}\mu(d)*i*j$

转化为首先枚举约数

$∑d=1min⁡(n,m)∑d∣in∑d∣jmμ(d)∗i∗j\sum\limits_{d=1}^{\min(n,m)}\sum\limits_{d|i}^{n}\sum\limits_{d|j}^{m}\mu(d)*i*j$

设 $i = i^{'} * d, j = j^{'} * d$ ，则可以进一步转化

$∑d=1min⁡(n,m)μ(d)∗d2∗∑i=1⌊nd⌋∑j=1⌊md⌋i∗j\sum\limits_{d=1}^{\min(n,m)}\mu(d)*d^2*\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}i*j$

前半段可以处理前缀和，后半段可以 $O (1)$ 求，设

$Q(n,m)=∑i=1n∑j=1mi∗j=n∗(n+1)2∗m∗(m+1)2Q(n,m)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}i*j=\frac{n*(n+1)}{2}*\frac{m*(m+1)}{2}$

显然可以 $O (1)$ 求解

到现在

$sum(n,m)=∑d=1min⁡(n,m)μ(d)∗d2∗Q(⌊nd⌋,⌊md⌋)sum(n,m)=\sum\limits_{d=1}^{\min(n,m)}\mu(d)*d^2*Q(\lfloor\frac nd \rfloor,\lfloor\frac md\rfloor)$

可以用数论分块求解

回带到原式中

$∑d=1min⁡(n,m)d∗sum(⌊nd⌋,⌊md⌋)\sum\limits_{d=1}^{\min(n, m)}d*sum(\lfloor\frac nd \rfloor,\lfloor\frac md\rfloor)$

又可以数论分块求解了

然后就做完啦

代码

/*
Author:loceaner
*/
#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define int long long
using namespace std;

const int A = 1e7 + 11;
const int B = 1e6 + 11;
const int mod = 20101009;
const int inf = 0x3f3f3f3f;

inline int read() {
	char c = getchar(); int x = 0, f = 1;
	for( ; !isdigit(c); c = getchar()) if(c == '-') f = -1;
	for( ; isdigit(c); c = getchar()) x = x * 10 + (c ^ 48);
	return x * f;
}

bool vis[A];
int n, m, mu[A], p[B], sum[A], cnt;

void getmu() {
    mu[1] = 1;
    int k = min(n, m);
    for (int i = 2; i <= k; i++) {
        if (!vis[i]) p[++cnt] = i, mu[i] = -1;
        for (int j = 1; j <= cnt && i * p[j] <= k; ++j) {
            vis[i * p[j]] = 1;
            if (i % p[j] == 0) break;
            mu[i * p[j]] = -mu[i];
        }
    }
    for (int i = 1; i <= k; i++) sum[i] = (sum[i - 1] + i * i % mod * mu[i]) % mod;
}

int Sum(int x, int y) { 
	return (x * (x + 1) / 2 % mod) * (y * (y + 1) / 2 % mod) % mod; 
}

int solve2(int x, int y) {
    int res = 0;
    for (int i = 1, j; i <= min(x, y); i = j + 1) {
        j = min(x / (x / i), y / (y / i));
        res = (res + 1LL * (sum[j] - sum[i - 1] + mod) * Sum(x / i, y / i) % mod) % mod;
    }
    return res;
}

int solve(int x, int y) {
    int res = 0;
    for (int i = 1, j; i <= min(x, y); i = j + 1) {
        j = min(x / (x / i), y / (y / i));
        res = (res + 1LL * (j - i + 1) * (i + j) / 2 % mod * solve2(x / i, y / i) % mod) % mod;
    }
    return res;
}

signed main() {
    n = read(), m = read();
    getmu();
    cout << solve(n, m) << '\n';
}