# 【BZOJ2693】jzptab （莫比乌斯反演）

## Output

$t$行，第$i$行为第$i$组询问的答案。

## HINT

$100000009$不是一个质数。

## 题解：

\begin{aligned}
ans&=\sum^n_{i=1}\sum^m_{j=1}\frac{ij}{gcd(i,j)}\\
&=\sum^n_d\sum^n_{i=1}\sum^m_{j=1}\frac{ij}{d}[gcd(i,j)=d]\\
&=\sum^n_d\sum^{\lfloor\frac{n}{d}\rfloor}_{i=1}\sum^{\lfloor\frac{m}{d}\rfloor}_{j=1}\frac{d^2ij}{d}[gcd(i,j)=1]\\
&=\sum^n_d\sum^{\lfloor\frac{n}{d}\rfloor}_{i=1}\sum^{\lfloor\frac{m}{d}\rfloor}_{j=1}dij\sum_{k\mid gcd(i,j)}\mu(k)&(\sum_{d\mid n}\mu(d)=[n=1])\\
&=\sum^n_k\mu(k)\sum^n_d\sum^{\lfloor\frac{n}{d}\rfloor}_{i=1}[k\mid i]\sum^{\lfloor\frac{m}{d}\rfloor}_{j=1}[k\mid j]dij\\
&=\sum^n_k\mu(k)\sum^n_d\sum^{\lfloor\frac{n}{kd}\rfloor}_{i=1}\sum^{\lfloor\frac{m}{kd}\rfloor}_{j=1}dk^2ij\\
&=\sum^n_k\mu(k)\sum^n_d\sum^{\lfloor\frac{n}{kd}\rfloor}_{i=1}i\sum^{\lfloor\frac{m}{kd}\rfloor}_{j=1}j\cdot dk^2\\
&设 T=kd\\
ans&=\sum^n_{T=1}\sum^{\lfloor\frac{n}{T}\rfloor}_{i=1}i\sum^{\lfloor\frac{m}{T}\rfloor}_{j=1}j\sum_{d\mid T}\mu(\frac{T}{d})\frac{T^2}{d}
\end{aligned}

$p\nmid i$，由于$i$和$p$互质而$g(x)$为积性函数，$g(x)=g(i\times p)=g(i)\times g(p)$

$p\mid i$，这个时候就有点不是很好搞了……

## CODE:

 1 #include<iostream>
2 #include<cstdio>
3 using namespace std;
4
5 #define mod 100000009LL
6 #define N 10000005
7 int t,n,m,cnt,ans;
8 long long pri[N],g[N],sum[N];
9 bool vis[N];
10
11 void init(){
12     sum[1]=g[1]=1;
13     for(int i=2;i<N;i++){
14         sum[i]=1LL*i*(i+1)/2%mod;
15         if(!vis[i]){
16             g[i]=(i-1LL*i*i%mod+mod)%mod;
17             pri[++cnt]=i;
18         }
19         for(int j=1;j<=cnt&&i*pri[j]<N;j++){
20             vis[i*pri[j]]=true;
21             if(i%pri[j])
22                 g[i*pri[j]]=g[i]*g[pri[j]]%mod;
23             else
24                 g[i*pri[j]]=g[i]*pri[j]%mod;
25         }
26     }
27     for(int i=1;i<N;i++)(g[i]+=g[i-1])%=mod;
28 }
29
30 int main(){
31     scanf("%d",&t);
32     init();
33     while(t--){
34         scanf("%d%d",&n,&m);
35         if(n>m)swap(n,m);
36         ans=0;
37         for(int i=1,pos=0;i<=n;i=pos+1){
38             pos=min(n/(n/i),m/(m/i));
39             ans+=sum[n/i]*sum[m/i]%mod*(g[pos]-g[i-1]+mod)%mod;
40             ans%=mod;
41         }
42         printf("%d\n",ans);
43     }
44 }

posted @ 2018-07-30 19:57  ezoiLZH  阅读(315)  评论(0编辑  收藏  举报