P3312 [SDOI2014]数表 解题笔记

原文创建时间：2021-03-07 13:31:05

这道题是需要离线询问的 数论 + 数据结构 题目，比较有意思。

这道题我本来建议改成黑题的，但是 nzhtl1477 不认为这题算得上黑题。嘤嘤嘤。

题目链接

P3312

题意简述

给出 $n,m,a$，求
$$\displaystyle \sum_{i=1}^n\sum_{j=1}^m\sigma(\gcd(i,j))[\sigma(\gcd(i,j)) \leq a]$$

多组询问，每次询问 $n,m,a$ 均改变。

解

不妨令 $n \leq m$。

先推式子：

$\displaystyle \ \ \ \ \sum_{i=1}^n\sum_{j=1}^m\sigma(\gcd(i,j))[\sigma(\gcd(i,j)) \leq a]$

$\displaystyle = \sum_{d=1}^n\sigma(d)[\sigma(d) \leq a]\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[i \perp j]$

$\displaystyle = \sum_{d=1}^n\sigma(d)[\sigma(d) \leq a]\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}\sum_{k|\gcd(i,j)}\mu(k)$

$\displaystyle = \sum_{d=1}^n\sigma(d)[\sigma(d) \leq a]\sum_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(k)\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[i|k]\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[j|k]$

$\displaystyle = \sum_{d=1}^n\sigma(d)[\sigma(d) \leq a]\sum_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(k)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor$

$\displaystyle = \sum_{t=1}^n\sum_{d|t}\sigma(d)[\sigma(d) \leq a]\mu\left(\frac{t}{d}\right)\left\lfloor\frac{n}{t}\right\rfloor\left\lfloor\frac{m}{t}\right\rfloor$

$\displaystyle = \sum_{t=1}^n\left\lfloor\frac{n}{t}\right\rfloor\left\lfloor\frac{m}{t}\right\rfloor\sum_{d|t}[\sigma(d) \leq a]\sigma(d)\mu\left(\frac{t}{d}\right)$

推导到这一步，如果没有 $\sigma(d) \leq a$ 的限制，就可以直接喜闻乐见地 $\mathrm{O} (n + \sqrt{n})$ 计算了。

考虑如何处理 $a$ 的限制。

数论分块计算答案时，每一块中都使用的是这一块中函数的前缀和。

那么，对于有序的 $a$ 限制，每次扩大范围时，可以将新增的部分加到前缀和中。

依照此思路，我们离线询问，按 $a$ 递增排序，每次将新增的部分加入前缀和，再计算回答当次询问。易见，前缀和宜用树状数组维护。

怎样计算新增部分呢？

记 $f = \mu \ast \sigma$，则对于每一个 $\sigma(d) \leq a$ 的 $d$，对所有 $dk \leq n$ 的 $f(dk)$ 加上 $\sigma(d)\mu(k)$。

这一部分的时间复杂度是 $\mathrm{O} (\log n \cdot n \ln n) = \mathrm{O} (n \log^2 n)$，一个 $\log$ 来自树状数组，另一个是枚举倍数（调和级数复杂度）。

不计离散询问和排序，算法总时间复杂度为 $\mathrm{O} (n \log^2 n + T\sqrt n \log n)$，每次回答时的 $\log$ 也来自树状数组。

代码

#include<bits/stdc++.h>
#define REG register

using namespace std;

namespace ACAH
{
	const int maxt = 2e4 + 7;
	const int maxn = 1e5 + 7;
	
	template<typename T> inline void qread(T &x)
	{
		x = 0; REG int f = 0; REG char ch = getchar();
		while(!isdigit(ch)) f |= (ch == '-'), ch = getchar();
		while(isdigit(ch)) x = (x << 1) + (x << 3) + (ch ^ 48), ch = getchar();
		x = f ? -x : x;
	}
	
	int T, bd;
	int ans[maxt];
	
	struct qust{
		int n, m, a, ser;
		const bool operator < (const qust &th) const
		{return a == th.a ? ser < th.ser : a < th.a;}
	} t[maxt];
	
	struct dpair{
		int v;
		int s;
		const bool operator < (const dpair &th) const
		{return s == th.s ? v < th.v : s < th.s;}
	} d[maxn];
	int pt;
	
	int pr[maxn], pc;
	bool ip[maxn];
	int mu[maxn];
	int sgm[maxn], pw[maxn];
	
	void linear_sieve()
	{
		mu[1] = 1, sgm[1] = 1, pw[1] = 1;
		for(REG int i = 2; i <= bd; i++) {
			if(!ip[i]) pr[++pc] = i, mu[i] = -1, sgm[i] = i + 1, pw[i] = i * i;
			for(REG int j = 1; j <= pc && 1ll * i * pr[j] <= bd; j++) {
				REG int c = i * pr[j]; ip[c] = true;
				if(i % pr[j]) mu[c] = -mu[i], sgm[c] = sgm[i] * (pr[j] + 1), pw[c] = pr[j] * pr[j];
				else {pw[c] = pw[i] * pr[j], sgm[c] = 1ll * sgm[i] * (pw[c] - 1) / (pw[i] - 1); break;}
			}
		}
		for(REG int i = 1; i <= bd; i++) d[i].v = i, d[i].s = sgm[i];
	}
	
	int c[maxn];
	
	inline int lowbit(int p) {return p & (-p);}
	inline void upd(int p, int x) {while(p <= bd) c[p] += x, p += lowbit( p );}
	inline int qry(int p) {REG int res = 0; while(p > 0) res += c[p], p -= lowbit( p ); return res;}
	inline int qry(int l, int r) {return qry(r) - qry(l - 1);}
	
	int work()
	{
		qread(T);
		
		for(REG int i = 1; i <= T; i++)
			qread(t[i].n), qread(t[i].m), qread(t[i].a), t[i].ser = i,
			bd = max(bd, min(t[i].n, t[i].m));
		
		linear_sieve();
		sort(t + 1, t + T + 1); sort(d + 1, d + bd + 1);
		
		for(REG int i = 1; i <= T; i++) {
			
			while(pt < bd && d[pt + 1].s <= t[i].a) {
				REG int cur = d[++pt].v; REG int curd = d[pt].s;
				for(REG int k = 1; 1ll * cur * k <= bd; k++)
					upd(cur * k, curd * mu[k]);
			}
			
			REG int N = t[i].n, M = t[i].m, ser = t[i].ser; if(N > M) swap(N, M);
			for(REG int l = 1, r; l <= N; l = r + 1) {
				r = min(N / (N / l), M / (M / l));
				ans[ser] += (N / l) * (M / l) * qry(l, r);
			}
		}
		
		for(REG int i = 1; i <= T; i++) printf("%lld\n", ans[i] & 2147483647);
		return 0;
	}
}

int main() {return ACAH::work();}