MyException - 我的异常网
当前位置:我的异常网» 数据结构与算法 » BZOJ3994 [SDOI2015]概数个数和

BZOJ3994 [SDOI2015]概数个数和

www.MyException.Cn  网友分享于:2013-08-22  浏览:0次
BZOJ3994 [SDOI2015]约数个数和

Description

 设d(x)为x的约数个数,给定N、M,求$\sum_{i=1}^n\sum_{j=1}^md(i, j)$

Input

输入文件包含多组测试数据。

第一行,一个整数T,表示测试数据的组数。
接下来的T行,每行两个整数N、M。

Output

 T行,每行一个整数,表示你所求的答案。

Sample Input

2
7 4
5 6

Sample Output

110
121

HINT

 1<=N, M<=50000

1<=T<=50000

Solution

计算$d(ij)$时,我们把$ij$的每个约数$d$映射到$(gcd(d, i), \frac{d}{gcd(d, i)})$,那么这两个数一定分别是$i, j$的因数,且$(a, b)$对应一个因数当且仅当$gcd(\frac ia, b) = 1$,所以

$$d(ij) = \sum_{x|i}\sum_{y|j} [gcd(\frac ix, y) = 1] = \sum_{x'|i}\sum_{y|j} [gcd(x', y) = 1]$$

于是(以下所有除法若不一定除尽,都指下取整)

$$\begin{aligned}
 \sum_{i=1}^N\sum_{j=1}^Md(ij)
&=\sum_{i=1}^N\sum_{j=1}^M\sum_{x|i}\sum_{y|j} [gcd(x, y) = 1]\\
&=\sum_{x=1}^N\sum_{y=1}^M [gcd(x, y) = 1]\left\lfloor\frac Nx\right\rfloor \left\lfloor\frac My\right\rfloor
\end{aligned}$$

我们令$f(d) = \sum_{x=1}^N\sum_{y=1}^M [gcd(x, y) = d]\left\lfloor\frac Nx\right\rfloor \left\lfloor\frac My\right\rfloor$,有

$$\begin{aligned}
\sum_{n|d}f(d)
&= \sum_{x=1}^N\sum_{y=1}^M [n|gcd(x, y)]\left\lfloor\frac Nx\right\rfloor \left\lfloor\frac My\right\rfloor\\
&= \sum_{i=1}^{\left\lfloor\frac Nn\right\rfloor}\sum_{j=1}^{\left\lfloor\frac Mn\right\rfloor} \left\lfloor\frac{\left\lfloor\frac Nn\right\rfloor}i \right\rfloor\left\lfloor\frac{\left\lfloor\frac Mn\right\rfloor}j\right\rfloor
\end{aligned}$$

如果我们令$t(n) = \sum_{i=1}^n \left\lfloor\frac ni\right\rfloor$,那上式就等于$t(\left\lfloor\frac Nn\right\rfloor)t(\left\lfloor\frac Mn\right\rfloor)$

于是$f(n) = \sum_{n|d} \mu(\frac dn)t(\left\lfloor\frac Nd\right\rfloor)t(\left\lfloor\frac Md\right\rfloor)$

$ans = f(1) = \sum_{n=1}^{min(N, M)}\mu(n)t(\left\lfloor\frac Nn\right\rfloor)t(\left\lfloor\frac Mn\right\rfloor)$

预处理$\mu(n)$的前缀和、$O(n\sqrt n)$预处理所有$t(n)$,查询时$O(\sqrt n)$即可。

代码:

#include <algorithm>
#include <cstdio>
typedef long long LL;
const int N = 50050;
LL t[N];
LL calcT(int n) {
  LL ans = 0;
  for (int i = 1, last; i <= n; i = last + 1) {
    last = n / (n / i);
    ans += n / i * (last - i + 1);
  }
  return ans;
}
bool mark[N];
int pr[N], pcnt = 0, mu[N], S[N];
void getMu() {
  mu[1] = 1;
  for (int i = 2; i < N; ++i) {
    if (!mark[i]) mu[pr[pcnt++] = i] = -1;
    for (int j = 0; j < pcnt && (LL)i * pr[j] < N; ++j) {
      mark[i * pr[j]] = 1;
      if (!(i % pr[j])) {
        mu[i * pr[j]] = 0;
        break;
      }
      mu[i * pr[j]] = -mu[i];
    }
  }
  for (int i = 1; i < N; ++i) S[i] = S[i - 1] + mu[i];
}
LL solve(int n, int m) {
  LL ans = 0;
  for (int i = 1, last; i <= n && i <= m; i = last + 1) {
    last = std::min(n / (n / i), m / (m / i));
    ans += t[n / i] * t[m / i] * (S[last] - S[i - 1]);
  }
  return ans;
}
int main() {
  for (int i = 1; i < N; ++i) t[i] = calcT(i);
  getMu();
  int T, n, m;
  scanf("%d", &T);
  while (T--) {
    scanf("%d%d", &n, &m);
    printf("%lld\n", solve(n, m));
  }
  return 0;
}

  

文章评论

旅行,写作,编程
旅行,写作,编程
2013年美国开发者薪资调查报告
2013年美国开发者薪资调查报告
程序员必看的十大电影
程序员必看的十大电影
科技史上最臭名昭著的13大罪犯
科技史上最臭名昭著的13大罪犯
如何区分一个程序员是“老手“还是“新手“?
如何区分一个程序员是“老手“还是“新手“?
程序员的一天:一寸光阴一寸金
程序员的一天:一寸光阴一寸金
那些争议最大的编程观点
那些争议最大的编程观点
做程序猿的老婆应该注意的一些事情
做程序猿的老婆应该注意的一些事情
鲜为人知的编程真相
鲜为人知的编程真相
程序员和编码员之间的区别
程序员和编码员之间的区别
我是如何打败拖延症的
我是如何打败拖延症的
那些性感的让人尖叫的程序员
那些性感的让人尖叫的程序员
为什么程序员都是夜猫子
为什么程序员都是夜猫子
10个帮程序员减压放松的网站
10个帮程序员减压放松的网站
5款最佳正则表达式编辑调试器
5款最佳正则表达式编辑调试器
 程序员的样子
程序员的样子
不懂技术不要对懂技术的人说这很容易实现
不懂技术不要对懂技术的人说这很容易实现
漫画:程序员的工作
漫画:程序员的工作
程序员都该阅读的书
程序员都该阅读的书
程序员最害怕的5件事 你中招了吗?
程序员最害怕的5件事 你中招了吗?
我跳槽是因为他们的显示器更大
我跳槽是因为他们的显示器更大
初级 vs 高级开发者 哪个性价比更高?
初级 vs 高级开发者 哪个性价比更高?
老美怎么看待阿里赴美上市
老美怎么看待阿里赴美上市
什么才是优秀的用户界面设计
什么才是优秀的用户界面设计
Google伦敦新总部 犹如星级庄园
Google伦敦新总部 犹如星级庄园
老程序员的下场
老程序员的下场
Java程序员必看电影
Java程序员必看电影
亲爱的项目经理,我恨你
亲爱的项目经理,我恨你
看13位CEO、创始人和高管如何提高工作效率
看13位CEO、创始人和高管如何提高工作效率
2013年中国软件开发者薪资调查报告
2013年中国软件开发者薪资调查报告
写给自己也写给你 自己到底该何去何从
写给自己也写给你 自己到底该何去何从
要嫁就嫁程序猿—钱多话少死的早
要嫁就嫁程序猿—钱多话少死的早
“肮脏的”IT工作排行榜
“肮脏的”IT工作排行榜
每天工作4小时的程序员
每天工作4小时的程序员
“懒”出效率是程序员的美德
“懒”出效率是程序员的美德
为啥Android手机总会越用越慢?
为啥Android手机总会越用越慢?
聊聊HTTPS和SSL/TLS协议
聊聊HTTPS和SSL/TLS协议
程序员应该关注的一些事儿
程序员应该关注的一些事儿
团队中“技术大拿”并非越多越好
团队中“技术大拿”并非越多越好
总结2014中国互联网十大段子
总结2014中国互联网十大段子
Web开发人员为什么越来越懒了?
Web开发人员为什么越来越懒了?
程序员周末都喜欢做什么?
程序员周末都喜欢做什么?
60个开发者不容错过的免费资源库
60个开发者不容错过的免费资源库
程序员眼里IE浏览器是什么样的
程序员眼里IE浏览器是什么样的
中美印日四国程序员比较
中美印日四国程序员比较
当下全球最炙手可热的八位少年创业者
当下全球最炙手可热的八位少年创业者
程序猿的崛起——Growth Hacker
程序猿的崛起——Growth Hacker
Java 与 .NET 的平台发展之争
Java 与 .NET 的平台发展之争
十大编程算法助程序员走上高手之路
十大编程算法助程序员走上高手之路
软件开发程序错误异常ExceptionCopyright © 2009-2015 MyException 版权所有