Description

有一张$n\times m$的数表,其第$i$行第$j$列($1 \le i \le n, 1 \le j \le m$)的数值为能同时整除$i$和$j$的所有自然数之和。给定$a$,计算数表中不大于$a$的数之和。

Input

输入包含多组数据。 输入的第一行一个整数$Q$表示测试点内的数据组数,接下来$Q$行,每行三个整数$n,m,a$描述一组数据。

$Q \le 2 \times 10^4, a \le 10^9, n, m \le 10^5$

Output

对每组数据,输出一行一个整数,表示答案$\mod 2^{31}$的值。

Solution

好题!

首先忽略$a$的限制,设$F(i)$为$i$的约数之和,那么答案就是

令$k = id$,原式化为

由于$F(i)$可以直接线性求出,暴力算出$G(k)$前缀和即可,那么每次询问只需要$O(\sqrt{min(n, m)})$。

接下来考虑$a$的限制,不妨离线处理询问,将询问按$a$排序。我们每次将在接下来的询问中需要计入答案的$F(i)$计入贡献,即维护$G(k)$的前缀和,用一个树状数组即可。

总的复杂度为$O(Q \sqrt n log n + n \ln n \log n)$

Source

#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cmath>

#define For(i, j, k) for(int i = j; i <= k; i++)
#define Forr(i, j, k) for(int i = j; i >= k; i--)

using namespace std;

const int N = 100010, M = 20010;

typedef unsigned long long LL;

int mu[N], prm[N], ind[N], n, tot;
LL sumf[N], F[N];

bool cmp(int x, int y){
return F[x] < F[y];
}

void init(){
n = N - 5;
mu[1] = F[1] = 1;
For(i, 2, n){
if(!prm[i]) prm[++tot] = i, mu[i] = -1, F[i] = i + 1, sumf[i] = 1LL * i * i;
for(int j = 1; j <= tot && 1LL * i * prm[j] <= n; j++){
int x = i * prm[j];
if(i % prm[j]){
F[x] = F[i] * (1 + prm[j]), sumf[x] = F[i] * prm[j] * prm[j];
mu[x] = -mu[i], prm[x] = true;
}
else{
F[x] = F[i] + sumf[i], sumf[x] = sumf[i] * prm[j];
mu[x] = 0, prm[x] = true;
break;
}
}
}
For(i, 1, n) ind[i] = i;
sort(ind + 1, ind + n + 1, cmp);
}

LL c[N];

inline int lowbit(int x){
return x & (-x);
}

void add(int x, LL v){
while(x <= n){
c[x] += v;
x += lowbit(x);
}
}

LL sum(int x){
LL ret = 0;
while(x){
ret += c[x];
x -= lowbit(x);
}
return ret;
}

struct Query{
int n, m, a, id;

bool operator < (const Query& A) const{
return a < A.a;
}
}A[M];

LL Ans[M];

int main(){
init();
int q;
scanf("%d", &q);
For(i, 1, q) scanf("%d%d%d", &A[i].n, &A[i].m, &A[i].a), A[i].id = i;
sort(A + 1, A + q + 1);
int p = 1;
For(i, 1, q){
while(p <= n && F[ind[p]] <= A[i].a){
int x = ind[p++];
for(int j = x; j <= n; j += x) add(j, F[x] * mu[j / x]);
}
int x = A[i].n, y = A[i].m;
for(int j = min(x, y), nxt; j; j = nxt){
nxt = max(x / (x / j + 1), y / (y / j + 1));
Ans[A[i].id] += 1LL * (x / j) * (y / j) * (sum(j) - sum(nxt));
}
}
For(i, 1, q) printf("%llu\n", Ans[i] & ((1u << 31) - 1));
return 0;
}