P3312 [SDOI2014]数表

题意描述：

洛谷

有一个 \(n\times m\) 的表格，第 \(i\) 行第 \(j\) 列的权值为 \(\sigma(\gcd(i,j))\) , 有 \(Q\) 组询问，每次问你一个 \(n\times m\) 的表格中，不大于 \(a\) 的数字之和。

数据范围： \(n,m\leq 10^5, Q\leq 2e4\)

前置芝士：

1.线性筛约数个数

设 \(n = p_1^{k_1}p_2^{k_2}...p_n^{k_n}\) ， 那么 \(\sigma(n) = (1+p_1^1+p_1^2+..+p_1^{k_1})(1+p_2^1+p_2^2+...p_2^{k2}) ...\)

考虑维护一下两个数组：

• \(sd(i)\) 表示 \(i\) 的约数个数和
• \(d(i)\) 表示 \(i\) 的最小质因子的 \((1+p_1^1+p_1^2+..p_1^{k_1})\) （实际上是个等比数列求和的形式）

还是考虑分三种情况来考虑：

• \(i\) 为质数，显然 \(sd(i) = i+ 1\) , \(d(i) = i+1\)

• 当 \(i \% prime[j] \neq 0\) 的时候, \(prime[j]\) 是 \(i\times prime[j]\) 的最小质因子，那么

  \(sd(i\times prime[j]) = sd(i) \times sd(prime[j])\)

  \(d(i\times prime[j]) = prime[j] + 1\)

• 当 \(i\% \ prime[j] = 0\) 的时候，\(i\times prime[j]\) 的最小质因子次数为 \(i\) 的最小质因子次数加1。

  \(sd(i\times prime[j]) = {sd(i)\over d(i)} \times(d(i)\times prime[j]) + 1\)

  \(d(i\times prime[j]) = d(i) \times prime[j] + 1\)

然后在线性筛的时候分情况讨论即可。

solution：

我们要求的其实是：

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

还是先枚举一下 \(gcd\) :

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

先把 \(\gcd(i,j) == d\) 反演掉：

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

\(\displaystyle\sum_{d=1}^{n}\sigma(d)[\sigma(d) \leq a] \sum_{i=1}^{n\over d}\sum_{j=1}^{m\over d} \sum_{p|i,p|j}\mu(p)\)

\(\displaystyle\sum_{d=1}^{n}\sigma(d) [\sigma(d) \leq a] \sum_{p=1}^{n\over d} \mu(p) {n\over {dp}} {m\over dp}\)

\(dp\) 化 \(Q\) 可得：

\(\displaystyle\sum_{Q=1}^{n} {n\over Q}{m\over Q} \sum_{d\mid Q} \sigma(d)[\sigma(d)\leq a] \mu({n\over Q})\)

对于每组询问的话，我们很难直接做，考虑离线下来，按 \(a\) 从小到大排一下序。

前面的那一坨可以用数论分块来解决，设 \(\displaystyle f(n) = \sum_{d\mid n} \sigma(d) \mu({n\over d})\) 。

因为我们把询问按 \(a\) 从小到大排序，枚举询问时，\(a\) 的增大会使得一些 \(\sigma(d) \leq a\) ,这时候我们需要把他的贡献算上去。

对于 \(\sigma(d)\) 不难发现它会对他的倍数 \(i\times d\) 的贡献为 \(\sigma(d) \mu(i)\) 。

调和级数来枚举 \(d\) 的倍数，在拿个树状数组维护一下 $f(n) $ 前缀和即可。

然后前面的整除分块，后面的就是 \(f(r)-f(l-1)\) ，在树状数组上查询即可。

#include<iostream>
#include<cstdio>
#include<algorithm>
using namespace std;
#define int long long
const int p = (1<<31);
const int N = 1e5+10;
int cntq,maxn,tot,l = 1;
int prime[N],mu[N],tr[N],omg[N],d[N],ans[N];
bool check[N];
struct QAQ
{
    int w,id;
    QAQ(){}
    QAQ(int a,int b){id = a, w = b;}
}e[N];
struct node
{
    int n,m,a,id;
}q[N];
inline int read()
{
    int s = 0,w = 1; char ch = getchar();
    while(ch < '0' || ch > '9'){if(ch == '-') w = -1; ch = getchar();}
    while(ch >= '0' && ch <= '9'){s = s * 10 + ch - '0'; ch = getchar();}
    return s * w;
}
int lowbit(int x){return x & -x;}
void chenge(int x,int val)
{
    for(; x <= N-5; x += lowbit(x)) tr[x] += val;
}
int query(int x)
{
    if(x == 0) return 0;
    int res = 0;
    for(; x; x -= lowbit(x)) res += tr[x];
    return res;
}
bool cmp(QAQ a,QAQ b)
{
    return a.w < b.w;
}
bool comp(node a,node b)
{
    return a.a < b.a;
}
void Add(int x)
{
    for(int i = 1; i * e[l].id <= maxn; i++) chenge(i*e[l].id,omg[e[l].id] * mu[i]);
}
void YYCH()
{
    mu[1] = 1; omg[1] = 1;
    for(int i = 2; i <= maxn; i++)
    {
        if(!check[i])
        {
            prime[++tot] = i;
            mu[i] = -1;
            omg[i] = i+1;
            d[i] = i+1;
        }
        for(int j = 1; j <= tot && i * prime[j] <= maxn; j++)
        {
            check[i * prime[j]] = 1;
            if(i % prime[j] == 0)
            {
                mu[i * prime[j]] = 0;
                d[i * prime[j]] = (d[i] * prime[j]) + 1;
                omg[i * prime[j]] = (omg[i] / d[i]) * (d[i] * prime[j] + 1);
                break;
            }
            else
            {
                mu[i * prime[j]] = -mu[i];
                d[i * prime[j]] = prime[j] + 1;
                omg[i * prime[j]] = omg[i] * omg[prime[j]];
            }
        } 
    }
    for(int i = 1; i <= maxn; i++) e[i] = QAQ(i,omg[i]);
    sort(e+1,e+maxn+1,cmp);
}
int slove(int n,int m,int id)
{
    int res = 0;
    while(l <= maxn && e[l].w <= q[id].a) Add(l), l++;
    for(int l = 1, r; l <= n; l = r+1)
    {
        r = min(n/(n/l),m/(m/l));
        res += (n/l) * (m/l) * (query(r)-query(l-1));
    }
    return res;
}
signed main()
{
    cntq = read();
    for(int i = 1; i <= cntq; i++)
    {
        q[i].n = read();
        q[i].m = read();
        q[i].a = read();
        q[i].id = i;
        if(q[i].n > q[i].m) swap(q[i].n,q[i].m);
        maxn = max(maxn,q[i].n);
    }
    YYCH();
    sort(q+1,q+cntq+1,comp);
    for(int i = 1; i <= cntq; i++) ans[q[i].id] = slove(q[i].n,q[i].m,i);
    for(int i = 1; i <= cntq; i++) printf("%lld\n",ans[i] % p);
    return 0;
}