题目链接

题目描述

有一张N*m的数表，其第i行第j列（1 < =i < =n，1 < =j < =m）的数值为能同时整除i和j的所有自然数之和。给定a，计算数表中不大于a的数之和。

题解

假设没有a的限制
每次是求：
$\sum^n_{i=1}\sum^m_{j=1}\sum_{d|gcd(i,j)}d\;\;\;$ 枚举d
$=\sum_{d=1}^{min(n,m)}d*\dfrac n d*\dfrac m d$

是不是化的太快了，好像没什么用。加上a的限制就不会做了。
看一下原来每个位置的数是gcd(i,j)的约数和。$设F(x)表示x的约数和$
那么要求的是:
$\sum^n_{i=1}\sum^m_{j=1}F(gcd(i,j))$
$=\sum_d^{min(n,m)}F(d)\sum_{i=1}^{\frac nd}\sum_{j=1}^{\frac md}[gcd(i,j)=1]$
$=\sum_d^{min(n,m)}F(d)\sum_{p=1}^{min(\frac nd,\frac md)}\mu(p)\dfrac n {pd}*\dfrac m {pd}$
$=\sum_{T=1}^{min(n,m)}\dfrac nT*\dfrac mT\sum_{d|T}\mu(\frac Td)F(d)$

P.S:你有没有发现$\sum_{d|T}\mu(\frac Td)F(d)=d$
这下有a的限制也好办了，我们不就是要快速搞一个后面那一坨的和嘛，有a的限制的话，把a排序，把F(d)按顺序插到一个树状数组里就可以了。

$F(x)可以线性筛。$

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstdlib>
#include<cstring>
using namespace std;
const int N=1e5+10;
typedef long long ll;
const ll mod=(1ll<<31);
inline int read()
{
    int x=0;char ch=getchar();int t=1;
    for(;ch>'9'||ch<'0';ch=getchar()) if(ch=='-') t=-1;
    for(;ch>='0'&&ch<='9';ch=getchar()) x=(x<<1)+(x<<3)+(ch-48);
    return x*t;
}
inline ll readl()
{
    ll x=0;char ch=getchar();ll t=1;
    for(;ch>'9'||ch<'0';ch=getchar()) if(ch=='-') t=-1;
    for(;ch>='0'&&ch<='9';ch=getchar()) x=(x<<1)+(x<<3)+(ch-48);
    return x*t;
}
int pri[N];int cnt=0;
int mu[N];
int c[N];//某数的质因子中当前最小的次数的幂
int g[N];//某数的质因子中当前最小的次数的数的幂的和
int sum[N];//某数的约数个数和
bool vis[N];
int id[N];
//若n=p1^a1 p2^a2 p3^a3...pn^an
//那么因为质因子是从小到大来枚举的,那么可以保证每次往数中添加的是最小的质因子
inline void prepare(int MAXN)
{
    mu[1]=1;vis[1]=1;id[1]=1;sum[1]=1;
    for(register int i=2;i<=MAXN;++i)
    {
        id[i]=i;
        if(!vis[i]) {pri[++cnt]=i;mu[i]=-1;c[i]=i;g[i]=1+i;sum[i]=i+1;}
        for(register int j=1;j<=cnt&&(1ll*i*pri[j]<=MAXN);++j)
        {
            register int x=i*pri[j];
            vis[x]=1;
            if(i%pri[j]==0){
                mu[x]=0;c[x]=c[i]*pri[j];
                g[x]=g[i]+c[x];
                sum[x]=(sum[i]/g[i])*g[x];
                break;
            }
            mu[x]=-mu[i];c[x]=pri[j];
            g[x]=pri[j]+1;
            sum[x]=sum[i]*g[x];
        }
    }
}
struct que{
    int n;int m;ll a;int id;
    inline bool operator <(que b)const{
        return a<b.a;
    }
}Q[20020];
int ans[20020];
inline bool cmp(int a,int b){return sum[a]<sum[b];}//按sum排序
int tr[N];int mn;
#define lowbit(a) ((a)&(-a))
inline void insert(int p,int x){while(p<=mn) {tr[p]+=x;p+=lowbit(p);}}
inline int Query(int p){int res=0;while(p){res+=tr[p];p-=lowbit(p);}return res;}
inline int solve(int n,int m)
{
    register int lst=0;register int l,r;
    register int res=0;
    for(l=1;l<=n;l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        register int now=Query(r);
        res+=(n/l)*(m/l)*(now-lst);
        lst=now;
    }
    return res;
}
int main()
{
    int T=read();register int n,m;register ll a;
    for(register int i=1;i<=T;++i) {
        n=read(),m=read(),a=readl();if(n>m) swap(n,m);
        Q[i]=(que){n,m,a,i};mn=max(mn,n);
    }
    prepare(mn);
    sort(Q+1,Q+1+T);sort(id+1,id+1+mn,cmp);
    register int h=1;
    for(register int i=1;i<=T;i++){
        for(;h<=mn&&sum[id[h]]<=Q[i].a;h++) {
            for(register int k=id[h];k<=mn;k+=id[h])
                if(mu[k/id[h]]) insert(k,sum[id[h]]*mu[k/id[h]]);
        }
        ans[Q[i].id]=solve(Q[i].n,Q[i].m);
    }
    for(register int i=1;i<=T;i++)
    {
        if(ans[i]<0) ans[i]+=2147483647,ans[i]++;
        printf("%d\n",ans[i]);
    }
}