P2522 HAOI2011 Problem b [莫比乌斯反演,数论分块]
P2522 HAOI2011
题意
对于给出的n个询问,每次求有多少个数对(x,y)(x,y)(x,y),满足a≤x≤ba≤x≤ba≤x≤b,c≤y≤dc≤y≤dc≤y≤d,且gcd(x,y)=kgcd(x,y) = kgcd(x,y)=k,gcd(x,y)gcd(x,y)gcd(x,y)函数为xxx和yyy的最大公约数.
题解
即求式子∑x=ab∑y=cd[gcd(x,y)=k]\sum_{x=a}^b\sum_{y=c}^d[gcd(x,y)=k]∑x=ab∑y=cd[gcd(x,y)=k].
记f(n,m)=∑x=1n∑y=1m[gcd(x,y)=k]f(n,m)=\sum_{x=1}^n\sum_{y=1}^m[gcd(x,y)=k]f(n,m)=∑x=1n∑y=1m[gcd(x,y)=k]
根据二维前缀和公式,可以将式子转换成:
∑x=ab∑y=cd[gcd(x,y)=k]=f(b,d)+f(a−1,c−1)−f(a−1,d)−f(b,c−1)\sum_{x=a}^b\sum_{y=c}^d[gcd(x,y)=k]=f(b,d)+f(a-1,c-1)-f(a-1,d)-f(b,c-1)∑x=ab∑y=cd[gcd(x,y)=k]=f(b,d)+f(a−1,c−1)−f(a−1,d)−f(b,c−1)
因此我们只要能得到f(n,m)f(n,m)f(n,m)的计算方法即可.
求f(n,m)f(n,m)f(n,m)的套路非常明显:莫比比乌斯反演
由于∑k∣d∑x=1n∑y=1m[gcd(x,y)=d]=⌊nk⌋⌊mk⌋\sum_{k|d}\sum_{x=1}^n\sum_{y=1}^m[gcd(x,y)=d] = \lfloor \frac{n}{k} \rfloor \lfloor \frac{m}{k} \rfloor∑k∣d∑x=1n∑y=1m[gcd(x,y)=d]=⌊kn⌋⌊km⌋.
反演得到f(n,m)=∑k∣dμ(dk)⌊nd⌋⌊md⌋f(n,m)=\sum_{k|d}\mu(\frac{d}{k})\lfloor \frac{n}{d} \rfloor \lfloor \frac{m}{d} \rfloorf(n,m)=∑k∣dμ(kd)⌊dn⌋⌊dm⌋
令t=d/kt = d/kt=d/k.
f(n,m)=∑t=1n/dμ(t)⌊nkt⌋⌊mkt⌋f(n,m)=\sum_{t=1}^{n/d}\mu(t)\lfloor \frac{n}{kt} \rfloor \lfloor \frac{m}{kt} \rfloorf(n,m)=∑t=1n/dμ(t)⌊ktn⌋⌊ktm⌋
对上式子进行分块计算,可以将时间复杂度从O(n)O(n)O(n)降至O(n)O(\sqrt{n})O(n).
总结
对于形如f(x)=∑x∣dμ(dx)g(⌊nd⌋)f(x)=\sum_{x|d}\mu(\frac{d}{x})g(\lfloor \frac{n}{d} \rfloor)f(x)=∑x∣dμ(xd)g(⌊dn⌋)这样的式子,我们都可以用t=d/xt=d/xt=d/x代换后数论分块进行加速.
f(x)=∑t=1n/xμ(t)g(⌊nxt⌋)f(x)=\sum_{t=1}^{n/x}\mu(t)g(\lfloor \frac{n}{xt} \rfloor)f(x)=∑t=1n/xμ(t)g(⌊xtn⌋)
时间复杂度从O(n)O(n)O(n)降至O(n)O(\sqrt{n})O(n).
代码
#include <iostream>
#include <algorithm>
#include <cstring>
#define pr(x) std::cout << #x << ':' << x << std::endl
#define rep(i,a,b) for(int i = a;i <= b;++i)const int N = 50000;int prime[N+10],zhi[N+10],mu[N+10],pcnt;void sieve() {zhi[1] = mu[1] = 1;for(int i = 2;i <= N;++i) {if(!zhi[i]) {mu[i] = -1;prime[pcnt++] = i;}for(int j = 0;j < pcnt && prime[j]*i <= N;++j) {zhi[i*prime[j]] = 1;if(i % prime[j] == 0) {mu[i*prime[j]] = 0;break;}else mu[i*prime[j]] = -mu[i];}}for(int i = 1;i <= N;++i) mu[i] += mu[i-1];
}
int a,b,c,d,k,T;
int calc(int n,int m) {int ans = 0;int lim = std::min(n/k,m/k);for(int i = 1,nx1,nx2,nxt;i <= lim;i=nxt+1) {nx1 = n/(n/i);nx2 = m/(m/i);nxt = nx1>nx2?nx2:nx1;ans += (mu[nxt]-mu[i-1])*(n/i/k)*(m/i/k);}return ans;
}
int main() {std::ios::sync_with_stdio(false);sieve();std::cin >> T;while(T--) {std::cin >> a >> b >> c >> d >> k;int ans = calc(b,d)+calc(a-1,c-1)-calc(a-1,d)-calc(b,c-1);std::cout << ans << std::endl;}return 0;
}
P2522 HAOI2011 Problem b [莫比乌斯反演,数论分块]相关推荐
- BZOJ 1101 Luogu P3455 POI 2007 Zap (莫比乌斯反演+数论分块)
BZOJ 1101 Luogu P3455 POI 2007 Zap (莫比乌斯反演+数论分块) 手动博客搬家: 本文发表于20171216 13:34:20, 原地址https://blog.csd ...
- BZOJ2301: [HAOI2011]Problem b[莫比乌斯反演 容斥原理]【学习笔记】
2301: [HAOI2011]Problem b Time Limit: 50 Sec Memory Limit: 256 MB Submit: 4032 Solved: 1817 [Submi ...
- [POI2007]ZAP-Queries [HAOI2011]Problem b 莫比乌斯反演
1,[POI2007]ZAP-Queries ---题面--- 题解: 首先列出式子:$$ans = \sum_{i = 1}^{n}\sum_{j = 1}^{m}[gcd(i, j) == d]$ ...
- HDU 6833 莫比乌斯反演 + 数论分块
给定下列式子: ∑a1=1n∑a2=1n...∑ax=1n(∏j=1xajk)f(gcd(a1,a2...,ax))⋅gcd(a1,a2...,ax)\sum_{a_1=1}^n\sum_{a_2 ...
- BZOJ2301:[HAOI2011]Problem b(莫比乌斯反演,容斥)
Description 对于给出的n个询问,每次求有多少个数对(x,y),满足a≤x≤b,c≤y≤d,且gcd(x,y) = k,gcd(x,y)函数为x和y的最大公约数. Input 第一行一个整数 ...
- BZOJ2301: [HAOI2011]Problem b(莫比乌斯反演)
Description 对于给出的n个询问,每次求有多少个数对(x,y),满足a≤x≤b,c≤y≤d,且gcd(x,y) = k,gcd(x,y)函数为x和y的最大公约数. Input 第一行一个整数 ...
- P2522 [HAOI2011]Problem b (莫比乌斯反演)
P2522 [HAOI2011]Problem b 题意: n个询问,在 a ≤ x ≤ b , c ≤ y ≤ d a\le x\le b, c \le y \le d a≤x≤b,c≤y≤d范围内 ...
- 欧拉函数+狄利克雷卷积+莫比乌斯函数+莫比乌斯反演+整除分块+杜教筛
Powered by:NEFU AB-IN 文章目录 欧拉函数 狄利克雷卷积 莫比乌斯函数 莫比乌斯反演 P3455 [POI2007]ZAP-Queries 整除分块 P2522 [HAOI2011 ...
- 小A的数学题(莫比乌斯反演数论)
小A的数学题(莫比乌斯反演&数论) 1.容斥 原始化简为 ∑ d = 1 n d 2 ∑ i = 1 n ∑ j = 1 m [ g c d ( i , j ) = d ] \sum\limi ...
最新文章
- 阿里 双11 同款,流量防卫兵 Sentinel go 源码解读
- ftm模块linux驱动,飞思卡尔k系列_ftm模块详解.doc
- 内存泄漏以及常见的解决方法
- loadrunner ftp linux,loadrunner ftp协议使用详细步骤
- linux 超好用的命令行工具
- 轻量级小型网站导航条
- 修改网课播放器倍速播放
- 【MySQL】轻松学习 唯一索引
- 【C语言】快速排序psort函数的用法及自主实现psort函数
- Rayman的绝顶之路——Leetcode每日一题打卡15
- MATLAB生成一段音乐《小星星》
- 作业一:我与软件工程相约
- 明德扬XILINX-K7-325T/410T核心板数据手册
- 【免费制作电子杂志】云展网教程 | 上传双页的PDF文档进行切割使杂志单页显示?
- python学习笔记之模块爬虫实战-爬取网站图片(十四)
- DBCP连接池常用参数详解
- webrtc不支持怎么调用设想吐_关于webRTC中video的使用实践
- 远程抄表-LoRa水电气无线抄表解决方案
- 高中学历父亲自学基因编辑,看五六百篇论文,自制药用级化合物救治罕见病儿子!...
- 用Python采集热门城市景点数据并简单制作数据可视化图
热门文章
- 怎么删除mysql的所有文件内容_mysql删除全部数据库
- java不朽神迹,不朽的神迹 Eternal Legacy HD v1.0.8
- java static new_java静态类new的对象是否能被回收?
- 7-52 两个有序链表序列的交集 (20 分)(思路加详解尾插法)come Boby!
- 文件的创建与读取 文件的数据添加
- 网络编程-网络分层的意义
- android创建类的包名称,如何知道/配置Xamarin Android生成的程序包名...
- 安装 java decompiler_Eclipse离线安装Java Decompiler插件(反编译)
- 2019-03-06-算法-进化(最接近的三数之和)
- CF 1638 E. Colorful Operations set 区间平推