原题链接

我天哪 大大的庆祝一下:

数论黑题 \(T1\) 达成!

激动地不行

记住套路:乱推 \(\gcd\),欧拉筛模板,然后乱换元,乱换式子,完了整除分块,欧拉筛和前缀和就解决了!

\[\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p\gcd(i\cdot j,i\cdot k,j\cdot k)\times \gcd(i,j,k)\times \left(\frac{\gcd(i,j)}{\gcd(i,k)\times \gcd(j,k)}+\frac{\gcd(i,k)}{\gcd(i,j)\times \gcd(j,k)}+\frac{\gcd(j,k)}{\gcd(i,j)\times \gcd(i,k)}\right)
\]

\[= \sum_{i=1}^n \sum_{j=1}^m \sum_{k=1}^p \frac{\gcd(i,j) \times \gcd(i,k) \times \gcd(j,k)}{\gcd(i,j,k)} \times \gcd(i,j,k) \times \left(\frac{\gcd(i,j)}{\gcd(i,k)\times \gcd(j,k)}+\frac{\gcd(i,k)}{\gcd(i,j)\times \gcd(j,k)}+\frac{\gcd(j,k)}{\gcd(i,j)\times \gcd(i,k)}\right)
\]

\[= \sum_{i=1}^n \sum_{j=1}^m \sum_{k=1}^p \gcd(i,j)^2 + \gcd(i,k)^2 + \gcd(j,k)^2
\]

\[= p \times \sum_{i=1}^n \sum_{j=1}^m \gcd(i,j)^2 + m \times \sum_{i=1}^n \sum_{k=1}^p \gcd(i,k)^2 + n \times \sum_{j=1}^m \sum_{k=1}^p \gcd(j,k)^2
\]

下面我们只需考虑 \(\sum_{i=1}^n \sum_{j=1}^m \gcd(i,j)^2\)的值即可。

\[\sum_{i=1}^n \sum_{j=1}^m \gcd(i,j)^2
\]

\[= \sum_{d=1}^{\min(n,m)} d^2 \sum_{i=1}^n \sum_{j=1}^m [\gcd(i,j) == d]
\]

\[= \sum_{d=1}^{\min(n,m)} d^2 \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{d} \rfloor} [gcd(i,j) == 1]
\]

这里我们要知道 \(\sum_{i=1}^a \sum_{j=1}^b [\gcd(i,j) == 1]\).

然后开始莫比乌斯反演。

令:

\[f_x = \sum_{i=1}^a \sum_{j=1}^b [\gcd(i,j) == x]
\]

\[F_x = \sum_{i=1}^a \sum_{j=1}^b [x | \gcd(i,j)]
\]

显然有:

\[F_x = \sum_{x|d} f_d
\]

\[f_1 = \sum_{d=1} F_d \times \mu_d
\]

\[= \sum_{d=1}^{\min(n,m)} \mu_d \times \lfloor \frac{n}{d} \rfloor \times \lfloor \frac{m}{d} \rfloor
\]

回到原式:

\[= \sum_{d=1}^{\min(n,m)} d^2 \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)} d^2 \sum_{g=1}^{\min(\lfloor \frac{n}{d} \rfloor,\lfloor \frac{m}{d} \rfloor)} \mu_d \lfloor \frac{n}{g \times d} \rfloor \lfloor \frac{m}{g \times d} \rfloor
\]

设 \(T = g \times d\) ,换 \(d|T\) 枚举。

\[= \sum_{T=1}^{\min(n,m)} \lfloor \frac{n}{T} \rfloor \times \lfloor \frac{m}{T} \rfloor \times \sum_{d|T} d^2 \mu{\lfloor \frac{T}{d} \rfloor}
\]

对于前面的,我们整除分块。

后面的,是个卷积,两个还都是积性函数,用欧拉筛一下。

(欧拉筛能解决所有的积性函数。如果不能,就再筛一遍。)

提醒一句:我们是要对 \(n\) , \(m\) ,\(p\) 求三次,不要只求一次啊。

突然感觉黑题也不怎难,会了套路就好

#pragma GCC optimize(2)
#include<bits/stdc++.h>
using namespace std; typedef long long ll;
const ll MOD=1e9+7;
const int N=2e7+1; inline ll read(){char ch=getchar();int f=1;while(ch<'0' || ch>'9') {if(ch=='-') f=-f; ch=getchar();}
ll x=0;while(ch>='0' && ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar();return x*f;} ll T,n,m,p,cnt=0; bool h[N];
ll prime[N],sum[N],f[N]; inline ll solve(ll n,ll m) {
ll ans=0;
for(ll i=1,t;i<=min(n,m);i=t+1) {
t=min(n/(n/i),m/(m/i));
ans=(ans+(n/i)*(m/i)%MOD*((sum[t]-sum[i-1]+MOD)%MOD)%MOD)%MOD;
} return ans;
} //整除分块 inline void Euler() {
f[1]=1;
for(ll i=2;i<N;i++) {
if(!f[i]) prime[++cnt]=i,f[i]=(i*i%MOD-1+MOD)%MOD;
for(ll j=1;j<=cnt && i*prime[j]<N;j++) {
f[i*prime[j]]=1;
if(i%prime[j]==0) {
f[i*prime[j]]=f[i]*prime[j]%MOD*prime[j]%MOD;
break;
} else f[i*prime[j]]=f[i]*f[prime[j]]%MOD;
}
}
for(ll i=1;i<N;i++) sum[i]=(sum[i-1]+f[i])%MOD;
} //欧拉筛,做前缀和 int main(){
Euler(); T=read(); while(T--) {
n=read(),m=read(),p=read();
printf("%lld\n",(p*solve(n,m)%MOD+m*solve(n,p)%MOD+n*solve(m,p)%MOD)%MOD);
} //别忘了3次轮换
return 0;
}

最新文章

  1. Oracle错误:动态执行表不可访问,本会话自动统计被禁止,关闭自动统计之后的问题
  2. [SAP ABAP开发技术总结]日期函数
  3. ASP.NET MVC 中使用 AjaxFileUpload 插件时,上传图片后不能显示(预览)
  4. 四则运算小程序测试--c++--软件工程课
  5. Tesseract-OCR牛刀小试:模拟请求时的验证码识别
  6. Android发送通知栏通知
  7. stm32之Systick(系统时钟)
  8. SilkTest Q&amp;A 2
  9. 关于Python 解包,你需要知道的一切
  10. 一起学python-语法
  11. winform无需安装pdf阅读器打开pdf文件
  12. Flume架构以及应用介绍
  13. [问题解决]基于注解配置dubbo遇到ConnectionLoss for /dubbo/xxx问题解决
  14. web技术栈中不可或缺的Linux技术
  15. win10下安装scala流程及问题
  16. BZOJ 4407 于神之怒加强版 (莫比乌斯反演 + 分块)
  17. 33个与众不同的Web表单设计
  18. AllowOverride None
  19. JavaScript - Scope of variables
  20. System V 共享内存 和 系列函数

热门文章

  1. 醉酒驾驶VS睡眠不足,哪个更危险
  2. Android Studio NDK编程初探
  3. 无线个人区域网WPAN
  4. Cenots 7 通过Yum 安装Node.js 报错问题
  5. 关于Newtonsoft.Json引用报错
  6. PAT-进制转换
  7. sql--自链接(推荐人)
  8. 网络编程概念 和OSI七层结构简介
  9. Layabox UI 笔记
  10. 将url转化成file文件