[bzoj2693&bzoj2194]Crash的数字表格&jzptab

题目大意

$$ans=\sum_{i=1}^n\sum_{j=1}^m[i,j]$$
求出ans并模一个数mo（两道题mo不同但都是常数，其中一个是质数另外一个不是）
两题的区别在于询问是否多组。
n,m<=10^7

单组询问

首先假设$n<=m$

$$ans=\sum_{i=1}^n\sum_{j=1}^m[i,j]$$
$$ans=\sum_{i=1}^n\sum_{j=1}^m\frac{ij}{(i,j)}$$
令 $d=(i,j)$
$$ans=\sum_{d=1}^nf[d]*d^{-1}$$
其中
$$f[d]=\sum_{i=1}^n\sum_{j=1}^mi*j[(i,j)=d]$$
$$f[d]=\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}d^2ij[(i,j)=1]$$
这条式子可以让我们分块每个块的f值和可以轻易统计。
我们显然可以莫比乌斯反演然后变成
$$f[d]=\frac{1}{4}d^2\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i^2*sum(\lfloor\frac{n}{d*i}\rfloor,\lfloor\frac{m}{d*i}\rfloor)*\mu[i]$$
所以
$$ans=\sum_{d=1}^n\frac{1}{4}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i^2*sum(\lfloor\frac{n}{d*i}\rfloor,\lfloor\frac{m}{d*i}\rfloor)*\mu[i]$$
其中 $sum(x,y)=x*(x+1)*y*(y+1)$
然后显然分块做可以o(n)解决单次询问。

多组询问

设$T=d*i$
那么

$$ans=\sum_{T=1}^nsum(\lfloor\frac{n}{T}\rfloor,\lfloor\frac{m}{T}\rfloor)*\sum_{d|T}\frac{1}{4}d*(\frac{T}{d})^2*\mu[\frac{T}{d}]$$
即
$$ans=\frac{1}{4}*\sum_{T=1}^nsum(\lfloor\frac{n}{T}\rfloor,\lfloor\frac{m}{T}\rfloor)*\sum_{d|T}T*d*\mu[d]$$
于是我们发现后面的东西只与T，有关，那么我们可以设那玩意是a[T]，然后预处理a数组，接着每次询问只需要根号复杂度。
a数组的预处理用线性筛。
具体的，对于一个比T的最小质因数小的数p，如果p|T，那么a[T*p]=a[T]*p。
否则a[T*p]=a[T]*(1-p)*p。

#include<cstdio>
#include<algorithm>
#define fo(i,a,b) for(i=a;i<=b;i++)
using namespace std;
typedef long long ll;
const int mo=100000009,maxn=10000000+10;
int pri[maxn],a[maxn],sum[maxn];
bool bz[maxn];
int i,j,k,l,t,n,m,ans,top,ca,inv;
int main(){
    fo(i,2,maxn-10){
        if (!bz[i]) pri[++top]=i,a[i]=1-i;
        fo(j,1,top){
            if ((ll)i*pri[j]>maxn-10) break;
            bz[i*pri[j]]=1;
            if (i%pri[j]==0){
                a[i*pri[j]]=a[i];
                break;
            }
            a[i*pri[j]]=(ll)(1-pri[j])*a[i]%mo;
        }
    }
    a[1]=1;
    fo(i,1,maxn-10) a[i]=(ll)a[i]*i%mo,a[i]=(ll)(a[i]+mo)%mo;
    fo(i,1,maxn-10) sum[i]=(ll)(sum[i-1]+a[i])%mo;
    inv=75000007;
    scanf("%d",&ca);
    while (ca--){
        ans=0;
        scanf("%d%d",&n,&m);
        if (n>m) swap(n,m);
        i=1;
        while (i<=n){
            j=min(n/(n/i),m/(m/i));
            ans=(ll)(ans+(ll)(n/i)*(n/i+1)%mo*(m/i)%mo*(m/i+1)%mo*(sum[j]-sum[i-1])%mo)%mo;
            ans=(ll)(ans+mo)%mo;
            i=j+1;
        }
        printf("%d\n",(ll)ans*inv%mo);
    }
}