【YBT2023寒假Day15 C】缺口一样（数论）（莫队）（根号分治）

缺口一样

题目链接：YBT2023寒假Day15 C

题目大意

给你一个序列，多次询问，每次问你一个区间这里面所有非空点集的最大公约数之积，对质数取模。

思路

首先质因子之间是独立的，考虑对于每个质数分别计算再乘起来。
那对于一个质数，如果区间 $[l,r]$ 有恰好 $c_i$ 个是 $p^i$ 的倍数，那 $p$ 贡献的次数就相当于对于每个 $i$，选 $[l,r]$ 一个子集使得所有数 $p$ 的次数最小值恰好是 $i$ 的方案数乘倍数。
那要的是恰好，我们 $c_i$ 是倍数，也就是至少 $i$ 次，那答案就是：（因为是非空子集所以要减一）
$\sum\limits_{i\geqslant 1}i(2^{c_i}-1-(2^{c_{i+1}}-1))$
$=\sum\limits_{i\geqslant 1}(2^{c_i}-1)$

那这个 $p$ 的贡献就是 $p^{\sum\limits_{i\geqslant 1}(2^{c_i}-1)}$，根据欧拉定理上面的指数可以对 $mod-1$ 取模。
那一个区间的问题似乎解决了。

---

那多个区间要怎么多，那我们一是维护每个质数的 $c_i$，而是再维护对应给答案的贡献。
那维护数组这个不难想到莫队。
那每次就是枚举加入或删除那个数的因子 $p$，把 $p$ 的贡献改了。
那因子个数 $\log n$，那复杂度就是 $n\sqrt{n}\log n$，过不了。
瓶颈显然在每次插入的时候复杂度是 $O(\log n)$ 的。

---

考虑到一个事情是 $\geqslant \sqrt{R}$（$R$ 是值域）的因子只会有一个，而且题目的值域是跟 $n$ 同阶的。
考虑再来个根号分治，对于大因子我们直接像上面的方法暴力统计。
小的个数不超过 $\sqrt{R}$，可以直接预处理出每个前缀每个小质数 $p$ 每个 $k$ 这个前缀里有多少个数是 $p^k$ 的倍数。就不需要用根号分治来做，也是 $O(n\sqrt{n})$ 的。

那么题目最后就是 $O(n\sqrt{n})$ 解决了。

代码

#include<cmath>
#include<cstdio>
#include<vector>
#include<algorithm>
#define ll long long
#define mo 998244353

using namespace std;

const int N = 1e5 + 100;
const int K = 405;
int n, m, prime[N], a[N], sum[N][700];
int blo[N], bl[N], br[N], B, bg[N], np[N];
vector <pair <int, int> > id;
struct Quest {
	int l, r, id;
}q[N];
ll inv[N], now, ans[N];
vector <ll> tim[N], invs[N];

void Init() {
	inv[0] = inv[1] = 1; for (int i = 2; i < N; i++) inv[i] = inv[mo % i] * (mo - mo / i) % mo;
	for (int i = 2; i < N; i++) {
		if (!np[i]) {
			prime[++prime[0]] = i;
			if (i <= K) {
				int now = 1, cnt = 0;
				while (now * i < N) now *= i, id.push_back(make_pair(now, i));
			}
			tim[prime[0]].push_back(i);
			invs[prime[0]].push_back(inv[i]);
			np[i] = prime[0];
		}
		else np[i] = 0;
		for (int j = 1; j <= prime[0] && i * prime[j] < N; j++) {
			np[i * prime[j]] = 1; if (i % prime[j] == 0) break;
		}
	}
}

bool cmp(Quest x, Quest y) {
	if (blo[x.l] != blo[y.l]) return blo[x.l] < blo[y.l];
	return x.r < y.r;
}

void insert(int now) {
	for (int i = 1; i <= prime[0] && prime[i] * prime[i] <= now; i++)
		while (now % prime[i] == 0) {
			ll tmp = tim[i][tim[i].size() - 1];
			tim[i].push_back(tmp * tmp % mo);
			tmp = invs[i][invs[i].size() - 1];
			invs[i].push_back(tmp * tmp % mo);
			now /= prime[i];
		}
	if (now > 1) {
		ll tmp = tim[np[now]][tim[np[now]].size() - 1];
		tim[np[now]].push_back(tmp * tmp % mo);
		tmp = invs[np[now]][invs[np[now]].size() - 1];
		invs[np[now]].push_back(tmp * tmp % mo);
	}
}

int num[N];

void ins(int x) {
	if (bg[x] == 1) return ;
	now = now * tim[np[bg[x]]][num[np[bg[x]]]] % mo;
	num[np[bg[x]]]++;
} 

void dec(int x) {
	if (bg[x] == 1) return ;
	num[np[bg[x]]]--;
	now = now * invs[np[bg[x]]][num[np[bg[x]]]] % mo;
}

int main() {
//	freopen("C2.in", "r", stdin);
//	freopen("write.txt", "w", stdout);
	freopen("lhsx2023.in", "r", stdin);
	freopen("lhsx2023.out", "w", stdout);
	
	Init();
	
	scanf("%d %d", &n, &m); B = sqrt(n);
	for (int i = 1; i <= n; i++) {
		blo[i] = (i - 1) / B + 1;
		if (!bl[blo[i]]) bl[blo[i]] = i;
		br[blo[i]] = i;
	}
	
	for (int i = 1; i <= n; i++) {
		scanf("%d", &a[i]);
		int now = a[i];
		for (int j = 1; j <= prime[0] && prime[j] <= K; j++)
			while (now % prime[j] == 0) now /= prime[j];
		bg[i] = now;
		for (int j = 0; j < id.size(); j++)
			sum[i][j] = sum[i - 1][j] + (a[i] % id[j].first == 0);
		insert(a[i]);
	}
	
	for (int i = 1; i <= m; i++) {
		scanf("%d %d", &q[i].l, &q[i].r); q[i].id = i;
	}
	sort(q + 1, q + m + 1, cmp);
	
	int l = 1, r = 0; now = 1;
	for (int i = 1; i <= m; i++)  {
		while (l > q[i].l) l--, ins(l);
		while (r < q[i].r) r++, ins(r);
		while (l < q[i].l) dec(l), l++;
		while (r > q[i].r) dec(r), r--;
		ans[q[i].id] = now;
		for (int j = 0; j < id.size(); j++) {
			int num = sum[r][j] - sum[l - 1][j];
			(ans[q[i].id] *= tim[np[id[j].second]][num] * inv[id[j].second] % mo) %= mo;
		}
	}
	
	for (int i = 1; i <= m; i++) printf("%lld\n", ans[i]);
	
	return 0;
}

__EOF__

本文作者：あおいSakura
本文链接：https://www.cnblogs.com/Sakura-TJH/p/YBT2023Day15_C.html
关于博主：评论和私信会在第一时间回复。或者直接私信我。
版权声明：本博客所有文章除特别声明外，均采用 BY-NC-SA 许可协议。转载请注明出处！
声援博主：如果您觉得文章对您有帮助，可以点击文章右下角【推荐】一下。您的鼓励是博主的最大动力！