思路:
莫比乌斯反演+整除分块+树状数组维护前缀和
代码:
#pragma GCC optimize(2)#pragma GCC optimize(3)#pragma GCC optimize(4)#include using namespace std;#define y1 y11#define fi first#define se second#define pi acos(-1.0)#define LL long long#define u unsigned int//#define mp make_pair#define pb push_back#define ls rt<<1, l, m#define rs rt<<1|1, m+1, r#define ULL unsigned LL#define pll pair #define pli pair #define pii pair #define piii pair #define pdd pair #define mem(a, b) memset(a, b, sizeof(a))#define debug(x) cerr << #x << " = " << x << "\n";#define fio ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);//headconst int N = 1e5 + 10;const LL MOD = 1LL<<31;int T, n, m;int prime[N], mu[N], g[N], gg[N], cnt;bool not_p[N];u ans[N];struct A { int id, g; bool operator < (const A & rhs) const { return g < rhs.g; }};struct B { int n, m, a, id; bool operator < (const B & rhs) const { return a < rhs.a; }}q[N];vector vc;void seive() { mu[1] = 1; g[1] = 1; for (int i = 2; i < N; ++i) { if(!not_p[i]) prime[++cnt] = i, mu[i] = -1, g[i] = 1+i, gg[i] = 1+i; for (int j = 1; j <= cnt && i*prime[j] < N; ++j) { not_p[i*prime[j]] = true; if(i%prime[j] == 0) { mu[i*prime[j]] = 0; g[i*prime[j]] = g[i]/gg[i]; gg[i*prime[j]] = gg[i]*prime[j]+1; g[i*prime[j]] *= gg[i*prime[j]]; break; } mu[i*prime[j]] = -mu[i]; g[i*prime[j]] = g[i]*(1+prime[j]); gg[i*prime[j]] = 1+prime[j]; } } for (int i = 1; i < N; ++i) vc.pb(A{i, g[i]}); sort(vc.begin(), vc.end());}struct BIT { u bit[N]; inline void add(int x, u a) { while(x < N) bit[x] = bit[x]+a, x += x&-x; } inline LL sum(int x) { LL res = 0; while(x) res = res+bit[x], x -= x&-x; return res; }}B;int Q;int main() { seive(); scanf("%d", &Q); for (int i = 1; i <= Q; ++i) { scanf("%d %d %d", &q[i].n, &q[i].m, &q[i].a); q[i].id = i; } sort(q+1, q+Q+1); int now = 0; for (int i = 1; i <= Q; ++i) { while(now < vc.size() && vc[now].g <= q[i].a) { for (int j = vc[now].id; j < N; j += vc[now].id) { B.add(j, vc[now].g*1u*mu[j/vc[now].id]); } ++now; } int n = q[i].n, m = q[i].m; int up = min(n, m); for (int l = 1, r; l <= up; l = r+1) { r = min(n/(n/l), m/(m/l)); u tmp = (n/l)*1u*(m/l)*1u*(B.sum(r)-B.sum(l-1)); ans[q[i].id] += tmp; } } for (int i = 1; i <= Q; ++i) printf("%lld\n", ans[i]%MOD); return 0;}