这个题让我们求的就是这个柿子

\[\sum_{i=1}^n\sum_{j=1}^m[i\perp j][j\perp k]\]

非常简单啊，\(i\perp j\)保证了\(\frac{i}{j}\)是一个最简分数，\(j\perp k\)保证\(\frac{1}{j}\)是一个无限不循环小数

于是我们要求的就是上面那个柿子啦

首先先看一个简单的问题，求\(1\)到\(n\)中和\(k\)互质的数的个数

我们设

\[F(d)=\sum_{i=1}^n[d|(i,k)]=[d|k]\times \left \lfloor \frac{n}{d} \right \rfloor\]

反演可得，我们要求的东西就是

\[\sum_{d|k}\mu(d)\left \lfloor \frac{n}{d} \right \rfloor\]

我们可以在\(\sigma(k)\)的时间内求出这个东西

现在我们可以来求我们的柿子啦

交换一下求和符号

\[\sum_{j=1}^m[j\perp k]\sum_{i=1}^n[i\perp j]\]

于是后面可以写成

\[\sum_{j=1}^m[j\perp k]\sum_{d|j}\mu(d)\left \lfloor \frac{n}{d} \right \rfloor\]

到前面来枚举\(d\)

\[\sum_{d=1}^n\left \lfloor \frac{n}{d} \right \rfloor\mu(d)\sum_{d|j,j<=m}[j\perp k] \]

\[=\sum_{d=1}^n\left \lfloor \frac{n}{d} \right \rfloor\mu(d)\sum_{t=1}^{\left \lfloor \frac{m}{d} \right \rfloor}[td\perp k]\]

我们发现\(td\)和\(k\)互质完全可以表示成\(t\perp k\)且\(d\perp k\)，于是我们又能把柿子写成这个样子

\[\sum_{d=1}^n\left \lfloor \frac{n}{d} \right \rfloor S(\left \lfloor \frac{m}{d} \right \rfloor)\mu(d)[d\perp k]\]

其中\(S(n)=\sum_{i=1}^n[i\perp k]\)

我们发现我们现在只需要求出\(\mu(d)[d\perp k]\)的前缀和就可以整除分块啦

数据范围不小啊，考虑杜教筛

。。。之后就考虑不出来了，之后成爷随便一化就秒掉了

我们设\(f(i)=\mu(i)[i\perp k]\)，\(g(i)=[i\perp k]\)

我们考虑一下\(f\times g\)

\[(f\times g)(n)=\sum_{d|n}[d\perp k][\frac{n}{d}\perp k]\mu(d)\]

我们发现只有当\(n\perp k\)的时候才能保证\(d\perp k\)和\(\frac{n}{d}\perp k\)同时满足，于是我们可以直接把\([d\perp k][\frac{n}{d}\perp k]\)变成\([n\perp k]\)提到外面来，变成

\[(f\times g)(n)=[n\perp k]\sum_{d|n}\mu(d)=[n=1]\]

至于\(g(i)\)的前缀和，发现那不就是\(S(i)\)吗

于是我们现在就可以杜教筛啦，复杂度是\(O(\sigma(k)n^{\frac{2}{3}})\)

不过有一点需要注意，就是当同时对\(n,m\)整除分块的时候不能使用类似于min_25筛的那一个\(trick\)，只能使用\(map\)或\(hash\)来记忆化

代码

#include
    
     #include
     
      #include
      
       #include
       
        #include
        
         #include
         
          #define re register#define LL long long#define max(a,b) ((a)>(b)?(a):(b))#define min(a,b) ((a)<(b)?(a):(b))using namespace std::tr1;const int maxn=2000005;int f[maxn],p[maxn>>1],mu[maxn],pre[maxn];int d[2005],n,m,k,H,M,Sqr,tot;unordered_map
          
            vis,ma,F;inline int calc(int n) { if(F[n]) return F[n]; int ans=0; for(re int i=1;i<=tot;i++) ans+=mu[d[i]]*(n/d[i]); return F[n]=ans;}inline int gcd(int a,int b) {return !b?a:gcd(b,a%b);}inline int solve(int n) { if(n<=M) return pre[n]; if(vis[n]) return ma[n]; int ans=1; for(re int l=2,r;l<=n;l=r+1) { r=n/(n/l); ans-=solve(n/l)*(calc(r)-calc(l-1)); } vis[n]=1,ma[n]=ans; return ans;}int main() { scanf("%d%d%d",&n,&m,&k);M=min(max(n,m),2000000); f[1]=1,mu[1]=1; for(re int i=2;i<=M;i++) { if(!f[i]) p[++p[0]]=i,mu[i]=-1; for(re int j=1;j<=p[0]&&p[j]*i<=M;j++) { f[p[j]*i]=1; if(i%p[j]==0) break; mu[p[j]*i]=-1*mu[i]; } } for(re int i=1;i<=k;i++) { if(k%i||!mu[i]) continue; d[++tot]=i; } for(re int i=1;i<=M;i++) if(gcd(i,k)!=1) pre[i]=pre[i-1]; else pre[i]=pre[i-1]+mu[i]; LL ans=0; H=min(n,m); for(re int l=1,r;l<=H;l=r+1) { r=min(n/(n/l),m/(m/l)); ans+=1ll*(n/l)*calc(m/l)*(solve(r)-solve(l-1)); } printf("%lld\n",ans); return 0;}