QOJ.ac

QOJ

IDProblemSubmitterResultTimeMemoryLanguageFile sizeSubmit timeJudge time
#684869#8327. 积性函数求和 $10^{13}$ 方便 FFT 版zydyWA 225ms3900kbC++179.5kb2024-10-28 16:18:022024-10-28 16:18:02

Judging History

你现在查看的是最新测评结果

  • [2024-10-28 16:18:02]
  • 评测
  • 测评结果:WA
  • 用时:225ms
  • 内存:3900kb
  • [2024-10-28 16:18:02]
  • 提交

answer

#include <iostream>
#include <vector>
#include <cmath>
#include <functional>
#include <algorithm>
using namespace std;
using u32 = unsigned int;
using u64 = unsigned long long;
using i64 = long long;
constexpr u32 mod = 1000000007;
inline constexpr u32 norm(const u32 x) { return x < mod ? x : x - mod; }
struct m32 {
	u32 x;
	m32() { }
	constexpr m32(const u32 _x) : x(_x) { }
};
inline constexpr m32 operator + (const m32 x1, const m32 x2) { return norm(x1.x + x2.x); }
inline constexpr m32 operator - (const m32 x1, const m32 x2) { return norm(x1.x + mod - x2.x); }
inline constexpr m32 operator - (const m32 x) { return x.x ? mod - x.x : 0; }
inline constexpr m32 operator * (const m32 x1, const m32 x2) { return static_cast<u64>(x1.x) * x2.x % mod; }
inline m32& operator += (m32& x1, const m32 x2) { return x1 = x1 + x2; }
inline m32& operator -= (m32& x1, const m32 x2) { return x1 = x1 - x2; }
inline m32& operator *= (m32& x1, const m32 x2) { return x1 = x1 * x2; }
inline bool operator == (const m32 x1, const m32 x2) { return x1.x == x2.x; }
inline bool operator != (const m32 x1, const m32 x2) { return x1.x != x2.x; }

struct block {
	static int v;
	vector<m32> sv;
	vector<m32> lv;
	block() : sv(v + 1, 0), lv(v + 1, 0) {}
};
int block::v;

inline void add(m32& x, const u32 a, const u32 b) { x = (x.x + 1ULL * a * b) % mod; }
inline void add(m32& x, const m32 a, const m32 b) { add(x, a.x, b.x); }
inline void sub(m32& x, const m32 a, const m32 b) { add(x, a.x, mod - b.x); }

block solve(const i64 N, const m32 A, const m32 B) {
	const int v = sqrt(N + 0.5);
	const int n_4 = sqrt(v + 0.5);
	const int n_8 = sqrt(n_4 + 0.5);
	block::v = v;

	vector<int> primes;
	vector<int> pi(v + 1);
	vector<bool> is_prime(v + 1);
	primes.push_back(1);
	is_prime[2] = true;
	for (int i = 3; i <= v; i += 2) is_prime[i] = true;
	for (int i = 3; i * i <= v; i += 2)
		for (int j = i * i; is_prime[i] && j <= v; j += (i << 1))
			is_prime[j] = false;
	for (int i = 2; i <= v; ++i) {
		pi[i] = pi[i - 1] + is_prime[i];
		if (is_prime[i]) primes.push_back(i);
	}

	vector<m32> sup;
	sup.resize(primes.size());
	u32 rec[4];
	rec[1] = 1;
	for (int i = 2; i <= 3; ++i)
		rec[i] = (i64)(mod - mod / i) * rec[mod % i] % mod;
	m32 inv3 = m32(rec[3]);

	const auto divide = [](i64 n, i64 d) -> i64 {return double(n) / d; };

	auto calc_medium = [&](const function<m32(u64)>& fp) {
		sup.clear();
		sup[0] = m32(0);
		for (int i = 1; i <= pi[v]; ++i) sup[i] = fp(primes[i]);

		vector<m32> lq(v + 1, 0);
		for (int i = 1; i <= pi[v]; ++i) lq[primes[i]] += sup[i];
		for (int i = 1; i <= v; ++i) lq[i] += lq[i - 1];

		block f;
		i64 r = pi[v];
		for (int i = pi[n_4] + 1; i <= pi[v]; ++i) {
			const i64 m = divide(N, primes[i]);
			if (i * primes[i] > m) break;
			while (r * primes[r] > m) r--;
			for (int j = i; j <= r; ++j) f.lv[divide(m, primes[j])] += sup[i] * sup[j];
		}
		for (int i = v - 1; i; --i) f.lv[i] += f.lv[i + 1];
		r = pi[v];
		for (int i = pi[n_4] + 1; i <= pi[v]; ++i) {
			const i64 m = divide(N, primes[i]);
			while (r * primes[r] > m) r--;
			const int j = max(primes[r], primes[i] - 1), t1 = divide(m, j), t0 = divide(v, primes[i]);
			for (int k = 1; k <= t0; ++k)
				add(f.lv[k], sup[i], lq[v] - lq[j]);
			for (int k = t0 + 1; k <= t1; ++k)
				add(f.lv[k], sup[i], lq[divide(m, k)] - lq[j]);
		}
		for (int k = 1; k <= n_4; ++k) {
			int t = v / k;
			i64 m = N / k;
			m32 ans = m32(0);
			for (int i = pi[n_4] + 1; i <= pi[t]; ++i) ans += sup[i] * f.lv[primes[i] * k];
			for (int i = pi[n_4] + 1; i <= pi[t]; ++i) {
				i64 q = (i64)primes[i] * primes[i];
				if (q * n_4 > m) break;
				ans += sup[i] * sup[i] * (lq[divide(m, q)] - lq[n_4]);
			}
			t = cbrt(m + 0.5);
			for (int i = pi[n_4] + 1; i <= pi[t]; ++i) ans += sup[i] * sup[i] * sup[i];
			f.lv[k] += ans * inv3;
		}
		for (int i = 1; i <= v; ++i) f.lv[i] += lq[v] - lq[n_4] + 1;
		for (int i = 1; i <= n_4; ++i) f.sv[i] = 1;
		for (int i = n_4 + 1; i <= v; ++i) f.sv[i] = lq[i] - lq[n_4] + 1;
		for (int i = 0; i <= v; ++i) lq[i] = 0;

		vector<m32> sq(v + 1, 0);
		int mm = v * max(log(N) / 10, 1.);
		int K = min(N, (i64)(mm * pow(N, 0.125))), B = N / K;
		m32 sum_s = 0;
		const auto add_s = [&](int x, m32 cnt) -> void {
			sum_s += cnt;
			while (x <= v) sq[x] += cnt, x += x & -x;
			};
		const auto add_l = [&](int x, m32 cnt) -> void {
			x = v + 1 - x;
			while (x <= v) lq[x] += cnt, x += x & -x;
			};
		function <void(int, int, m32)> dfs = [&](int n, int id, m32 fn) {
			if (n <= v) add_s(n, fn);
			else add_l(divide(N, n), fn);
			for (int i = id; i <= pi[v]; ++i) {
				i64 q = (i64)n * primes[i];
				if (q > K) break;
				dfs(q, i, fn * sup[i]);
			}
			};
		auto query_s = [&](int x) -> m32 {
			m32 ans = f.sv[x];
			while (x) ans += sq[x], x ^= x & -x;
			return ans;
			};
		auto query_l = [&](int x) -> m32 {
			x = v + 1 - x;
			m32 ans = sum_s;
			while (x) ans += lq[x], x ^= x & -x;
			return ans;
			};
		int K2, B2;
		for (int id = pi[n_4]; id > pi[n_8]; --id) {
			const int p = primes[id];
			const u64 m = N / p;
			dfs(p, id, sup[id]);
			const int t0 = B / p, t1 = min(B, v / p);
			for (int i = B; i > t1; --i) add(f.lv[i], sup[id], query_s(divide(m, i)));
			for (int i = t1; i > t0; --i) add(f.lv[i], sup[id], f.lv[i * p] + query_l(i * p));
			for (int i = t0; i; --i) add(f.lv[i], sup[id], f.lv[i * p]);
			K2 = mm * sqrt(p);
			B2 = N / K2;
			for (int i = B2; i > B; --i) f.lv[i] += query_l(i);
			K = K2, B = B2;
		}
		for (int i = B + 1; i <= v; ++i) f.lv[i] += query_l(i);
		for (int i = 1; i <= v; ++i)
			if (i & (i - 1))
				sq[i] += sq[i & (i - 1)];
		for (int i = 1; i <= v; ++i) f.sv[i] += sq[i];

		return f;
		};

	auto attach_small = [&](block&& f, const function<m32(u64)>& fp) {
		for (int id = pi[n_8]; id; --id) {
			const int p = primes[id], t = v / p;
			const i64 m = N / p;
			for (int j = 1, i = p; j <= t; ++j) {
				const m32 c1 = sup[id] * f.sv[j];
				for (int e = min(v + 1, i + p); i < e; ++i) f.sv[i] += c1;
			}
			for (int i = v; i > t; --i) add(f.lv[i], sup[id], f.sv[divide(m, i)]);
			for (int i = t; i >= 1; --i) add(f.lv[i], sup[id], f.lv[i * p]);
		}
		return move(f);
		};

	auto calc_large = [&](const function<m32(u64)>& fp, const function<m32(u64)>& sum_fp) {
		block f = attach_small(calc_medium(fp), fp);
		block res;
		for (int i = v; i >= 1; --i) {
			m32 ans = sum_fp(N / i) - f.lv[i];
			for (int j = 2; i * j <= v; ++j) sub(ans, fp(j), res.lv[i * j]);
			res.lv[i] = ans;
		}
		return res;
		};

	auto mult_large = [&](block&& f, const block& l) {
		for (int i = 1; i <= v; ++i) {
			for (int j = 1; i * j <= v; ++j)
				if (f.sv[j] != f.sv[j - 1])
					add(f.lv[i], f.sv[j] - f.sv[j - 1], l.lv[i * j]);
		}
		return move(f);
		};

	auto mult_powerful = [&](block&& f, const function<m32(u32, u32)>& fpp) {
		block h;
		function< void(u64, int, m32)> dfs = [&](u64 n, int beg, m32 coeff) -> void {
			if (n <= v) h.sv[n] += coeff;
			else h.lv[divide(N, n)] += coeff;
			u64 t = divide(N, n);
			for (int i = beg; i <= pi[v]; ++i) {
				const int p = primes[i];
				u64 q = 1ULL * p * p;
				if (q > t) break;
				for (int e = 2; q <= t; q *= p, ++e)
					dfs(n * q, i + 1, coeff * (fpp(p, e) - fpp(p, 1) * fpp(p, e - 1)));
			}
			};
		dfs(1, 1, 1);

		block res;
		for (int i = 1; i <= v; ++i)
			if (h.sv[i].x) {
				const i64 m = divide(N, i);
				const int t0 = sqrt(m + 0.5);
				for (int k = 1; k * i <= v; ++k)
					add(res.lv[k], h.sv[i], f.sv[v] - f.sv[t0]);
				for (int k = v / i + 1; k <= t0; ++k)
					add(res.lv[k], h.sv[i], f.sv[divide(m, k)] - f.sv[t0]);
			}
		for (int i = 1; i < v; ++i) res.lv[i] -= res.lv[i + 1];

		m32 sum_s = f.sv[v];
		f.sv[v] -= f.sv[v - 1];
		for (int i = v - 1; i; --i) h.lv[i] += h.lv[i + 1], f.sv[i] -= f.sv[i - 1];
		for (int j = 1; j <= v; ++j) {
			const int t0 = v / j;
			for (int t = 1; t < t0; ++t)
				add(res.lv[t], f.sv[j], h.lv[j * t] - h.lv[j * (t + 1)]);
			add(res.lv[t0], f.sv[j], h.lv[j * t0]);
		}
		for (int i = 1; i <= v; ++i)
			if (h.sv[i].x) {
				const i64 m = divide(N, i);
				const int t = sqrt(m + 0.5), t0 = v / i;

				for (int t = 1; t < t0; ++t) add(res.lv[t], h.sv[i], f.lv[t * i] - f.lv[i * (t + 1)]);
				add(res.lv[t0], h.sv[i], f.lv[t0 * i] - sum_s);
				for (int j = t0 + 1; j <= t; ++j) add(res.lv[divide(m, j)], h.sv[i], f.sv[j]);

				for (int j = 1; j <= t0; ++j) add(res.sv[i * j], h.sv[i], f.sv[j]);
			}

		for (int i = 1; i <= v; ++i) res.sv[i] += res.sv[i - 1];
		res.lv[v] += res.sv[v];
		for (int i = v - 1; i; --i) res.lv[i] += res.lv[i + 1];
		return res;
		};
	block l0 = calc_large([&](u64 n) { return 1; }, [&](u64 n) { return m32(n % mod); });
	block l1 = calc_large([&](u64 n) { return m32(n % mod); }, [&](u64 n) { return n %= mod, m32(n * (n + 1) / 2 % mod); });
	for (int i = 1; i <= v; ++i) l1.lv[i] = A * l0.lv[i] + B * l1.lv[i], l1.sv[i] = 0;

	auto fp = [&](u32 p) { return A + B * p; };
	auto fpp = [&](u32 p, u32 e) { return A * e + B * p; };
	return mult_powerful(attach_small(mult_large(calc_medium(fp), l1), fp), fpp);
}
signed main() {
	i64 T, n;
	m32 A, B;
	cin >> T;
	while (T--) {
		cin >> n >> A.x >> B.x;
		block f = solve(n, A, B);
		vector<m32> res(f.sv.begin() + 1, f.sv.end());
		res.insert(res.end(), f.lv.begin() + 1, f.lv.end());
		sort(res.begin(), res.end(), [](const m32 a, const m32 b) { return a.x < b.x; });
		res.erase(unique(res.begin(), res.end()), res.end());
		u32 ans = 0;
		for (auto x : res) ans ^= x.x;
		cout << ans << endl;
	}
	return 0;
}

Details

Tip: Click on the bar to expand more detailed information

Test #1:

score: 0
Wrong Answer
time: 225ms
memory: 3900kb

input:

10000
988 56395756 60780067
7923 293552717 438195956
4847 24236686 75287211
6694 74889751 64994726
3720 385482711 188638093
6021 2928896 248853035
6808 310612405 330739724
4062 15289930 175596707
9583 56394035 335888448
9798 151396947 371306315
4365 216662501 351771943
1359 165179730 80942360
1436 3...

output:

327673201
479926096
663026977
111144426
319473635
714222341
258160171
554751614
388475468
727146564
636145309
38152978
93106134
830427145
786158021
246674358
869531160
688788719
14789346
231357574
320951046
800315648
671005040
699086312
918860003
831589515
647662768
983625594
84728540
297605582
1307...

result:

wrong answer 1st numbers differ - expected: '6702293', found: '327673201'