hdu6390 GuGuFishtion 题解（Mobius反演+欧拉函数+除法分块）

题目：hdu6390.
题目大意：给定$n,m,p$，$p$为素数，求：
$$\sum _{i=1}^n\sum _{j=1}^m\frac{\varphi (ij)}{\varphi (i)\varphi (j)}(modp)$$

$1\le n,m\le min(10^6,p),1\le p\le 10^9+7$，最多三组数据.

先考虑如何化简$\frac{\varphi (ij)}{\varphi (i)\varphi (j)}$，显然除了$i,j$的共有质因子外其它因子都会被约分掉，所以我们先取出共有的部分$g=gcd(i,j)$.

考虑把共有部分唯一分解$g={\prod }_{l=1}^k{p}_l^{c_l}$，显然最后的答案就会变成：
$$\frac{\varphi (ij)}{\varphi (i)\varphi (j)}=\frac{g}{g{\prod }_{l=1}^k\frac{p_l-1}{p_l}}=\frac{gcd(i,j)}{\varphi (gcd(i,j))}$$

这个式子看起来比原来优美多了，所以开始推导：
$$\begin{gathered} \sum _{i=1}^n\sum _{j=1}^m\frac{\varphi (ij)}{\varphi (i)\varphi (j)} \\ =\sum _{i=1}^n\sum _{j=1}^m\frac{gcd(i,j)}{\varphi (gcd(i,j))} \\ =\sum _{d=1}^{min(n,m)}\frac{d}{\varphi (d)}\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }[gcd(i,j)=1] \\ =\sum _{d=1}^{min(n,m)}\frac{d}{\varphi (d)}\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\sum _{t|i\cup t|j}\mu (t) \\ =\sum _{d=1}^{min(n,m)}\frac{d}{\varphi (d)}\sum _{t=1}^{\lfloor \frac{min(n,m)}{d}\rfloor }\mu (t)\lfloor \frac{n}{dt}\rfloor \lfloor \frac{m}{dt}\rfloor \\ =\sum _{c=1}^{min(n,m)}\lfloor \frac{n}{c}\rfloor \lfloor \frac{m}{c}\rfloor \sum _{d|c}\frac{d}{\varphi (d)}\mu (\frac{c}{d}) \end{gathered}$$

我们设$f(n)={\sum }_{d|n}\frac{d}{\varphi (d)}\mu (\frac{n}{d})$，发现这是个积性函数，可以线性筛：
1.$i=1$，则$f(i)=1$.
2.$i$为素数，则$f(i)=\frac{1}{\varphi (1)}\mu (i)+\frac{i}{\varphi (i)}\mu (1)=\frac{1}{i-1}$.
3.$gcd(i,p)=1$，则$f(ip)=f(i)f(p)$.
4.当$p$为素数且$p|i$时，先考虑$i=p^c$的情况：
$$\begin{gathered} f(ip)=f(p^{c+1})=\sum _{d|p^{c+1}}\frac{d}{\varphi (d)}\mu (\frac{p^{c+1}}{d}) \\ =-\frac{p^c}{\varphi (p^c)}+\frac{p^{c+1}}{\varphi (p^{c+1})} \\ =-\frac{p^c}{(p-1)p^{c-1}}+\frac{p^{c+1}}{(p-1)p^c} \\ =0 \end{gathered}$$

然后我们用除法分块处理一下，就可以做到$O(Tn)$解决这个问题了.

代码如下：

#include<bits/stdc++.h>
  using namespace std;

#define Abigail inline void
typedef long long LL;

const int N=1000000;

int n,m;
LL inv[N+9],mod;

void pre_inv(int n,LL mod){
  inv[1]=1;
  for (int i=2;i<=n;++i)
    inv[i]=(mod-mod/i)*inv[mod%i]%mod;
}

int tp,pr[N+9],b[N+9];
LL f[N+9];

void sieve(int n){
  f[1]=1;tp=0;
  for (int i=2;i<=n;++i) b[i]=1;
  for (int i=2;i<=n;++i){
  	if (b[i]) pr[++tp]=i,f[i]=inv[i-1];
  	for (int j=1;j<=tp&&i*pr[j]<=n;++j){
  	  b[i*pr[j]]=0;
  	  if (i%pr[j]==0){f[i*pr[j]]=0;break;}
  	  f[i*pr[j]]=f[i]*f[pr[j]]%mod;
  	}
  }
}

LL Query(int n,int m){
  LL ans=0;
  for (int l=1,r;l<=n;l=r+1){
  	r=min(n/(n/l),m/(m/l));
  	ans=(ans+(f[r]-f[l-1]+mod)*(n/l)%mod*(m/l))%mod;
  }
  return ans;
}

Abigail into(){
  scanf("%d%d%lld",&n,&m,&mod); 
}

Abigail work(){
  if (n>m) swap(n,m);
  pre_inv(n,mod);
  sieve(n);
  for (int i=1;i<=n;++i){ 
    f[i]+=f[i-1];
    if (f[i]>=mod) f[i]-=mod;
  } 
}

Abigail outo(){
  printf("%lld\n",Query(n,m));
}

int main(){
  int T;
  scanf("%d",&T);
  while (T--){
  	into();
    work();
  	outo();
  }
  return 0;
}