题目描述

定义 \(d(n)\) 为 \(n\) 的正因数的个数,比如 \(d(2) = 2, d(6) = 4\)。

令 $ S_1(n) = \sum_{i=1}^n d(i) $

给定 \(n\),求 \(S_1(n)\)。

输入格式

第一行包含一个正整数 \(T\) (\(T \leq 10^5\)),表示数据组数。

接下来的 \(T\) 行,每行包含一个正整数 \(n\) (\(n < 2^{63}\))。

输出格式

对于每个 \(n\),输出一行一个整数,表示 \(S_1(n)\) 的值。

题解

显然 $ S_1(n) = \sum_{i=1}^n \left \lfloor \frac ni \right \rfloor $
整除分块可以做到\(O \left( \sqrt n \right)\),这是要死的复杂度。
但实际上,这玩意是有个优化套路的:
显然有$ \sum_{i=1}^n \left \lfloor \frac ni \right \rfloor = 2 \sum_{i=1}^{\lfloor \sqrt n \rfloor} \left \lfloor \frac ni \right \rfloor - \left \lfloor \sqrt n \right \rfloor^2 $
那么只需求 \(\sum\limits_{i=1}^{\lfloor \sqrt n \rfloor} \left \lfloor \dfrac ni \right \rfloor\) 的值
它显然等于 \(f(x) = \dfrac n x\)、直线 \(y = \sqrt n\)、\(x\) 轴和 \(y\) 轴所围成的区域 \(R\) 中整点 (格点) 的个数 (含边界,\(x\) 轴除外)。我们把它变成了一个区间数点问题。

但实际上这个函数是个凸函数,因此区域 \(R\) 的补集 \(R'\) 为一个凸集,我们很显然的想到用凸包去拟合,然后Pick定理随便搞一搞。

怎么找凸包呢?我们用SBT大力求斜率,然后弄个单调栈去搞。
记 \(B = \left \lfloor \sqrt n \right \rfloor\),我们从点 \(P_0 \left( \left \lfloor \dfrac nB \right \rfloor, B+1 \right)\) 开始,一步一步寻找凸包上所有的点。
(步骤一)我们在栈中加入两个分数 (为斜率的绝对值) \(\dfrac 01\) 和 \(\dfrac 11\),代表向量 \((1, 0)\) 和 \((1, 1)\),由于 \(f' \left( \sqrt n \right) = -1\),因此这部分所有斜率都在 \(-1 \sim 0\) 之间。 此外,这个栈需要满足自顶向上是单调递增的。
(步骤二)接下来,我们取出栈顶向量,将 \(P_0\) 持续与这个向量 (关于 \(x\) 轴的对称,下略) 相加,直到这个点 \(P_k\) 不在区域 \(R'\) 中 (即在区域 \(R\) 中)。由于这些分数是在 Stern-Brocot 树中产生的,因此一定是既约分数 (即分子与分母互素)。
因此每加一步,我们可以计算出这个横条的面积:设向量为 \((u, v)\),上一个点 (\(P_{k-1}\)) 的横坐标为 \(x\),则面积为 \(x v + \dfrac 12 (v + 1) (u - 1)\)。
(步骤三)接着我们要对栈进行调整。由于函数的斜率的绝对值单调递减,因此栈中的分数也需要单调递减。
故我们需要把值过大的分数弹出栈外,直到栈顶和它下面的元素与 \(P_k\) 相加后,前者在 \(R'\) 外,后者在 \(R'\) 内。把这两个向量记作 \(l\) 和 \(r\)。
(步骤四)然后我们在SBT上二分(不是说大力求嘛).
记 \(l = \dfrac {y_l} {x_l}, r = \dfrac {y_r} {x_r}\) (\(l > r\)),则 \(m = \dfrac {y_m} {x_m} = \dfrac {y_l + y_r} {x_l + x_r}\)。令 \(M = P_k + m\) (即 \(P_k\) 按照向量 \(m\) 平移后的点),如果 \(M\) 在 \(R'\) 中,则将 \(r\) 压入栈后令 \(r \gets m\),继续二分;否则,分以下两种情况讨论(二轮判断):
一:如果 \(\left| f' \left( x_k + x_m \right) \right| \leq r\) (其中 \(x_k\) 为 \(P_k\) 的横坐标)——由 \(f'(x) = - \dfrac n {x^2}\) 可得该条件等价于 \(n x_r \leq (x_k + x_m)^2 y_r\) ——则容易得到如果再迭代下去的话,所有的 \(P_k + m\) 都不会落在 \(R'\) 内,也就能说已经二分完毕了,因此只需保留当前的栈重新回到步骤 2 进入下一轮迭代。
二:如果 \(\left| f' \left( x_k + x_m \right) \right| > r\),则接下来的向量还有可能落入 \(R'\) 中,因此令 \(l \gets m\) 后继续二分。

这个做法的复杂度上界不超过 \(O \left( n^{1/3} \log n \right)\) ,因为斜率只有不超过\(O \left( n^{1/3} \right)\)个。(至于怎么证明文文也不会)。
后期函数的斜率非常小,都是 \(\dfrac 1k\) 的形式,会退化为 \(O(y)\),因此可以在 \(y \leq \sqrt[3]n\) 的部分中直接暴力计算

Code

#include <bits/stdc++.h>
const int N=1000005;
typedef long long ll;
struct qwq {ll x, y;qwq (ll x0 = 0, ll y0 = 0) : x(x0), y(y0) {}inline qwq operator + (const qwq &B) const {return qwq(x + B.x, y + B.y);}
};
ll n;
int top;
qwq stk[N];
inline void push(qwq x) { stk[++top]=x;}
inline void putint(__int128 x) {static char buf[36];if (!x) {putchar(48); return;} int i = 0;for (; x; buf[++i] = x % 10 | 48, x /= 10);for (; i; --i) putchar(buf[i]);
}inline bool check(ll x, ll y) {return n < x * y;}inline bool judge(ll x, qwq v) {return (__int128)n * v.x <= (__int128)x * x * v.y;}__int128 S1() {int i, crn = cbrt(n);ll srn = sqrt(n), x = n / srn, y = srn + 1;__int128 ret = 0;qwq L, R, M;push(qwq(1, 0)); push(qwq(1, 1));while(true) {for (L = stk[top--]; check(x + L.x, y - L.y); x += L.x, y -= L.y)ret += x * L.y + (L.y + 1) * (L.x - 1) / 2;if (y <= crn) break;for (R = stk[top]; !check(x + R.x, y - R.y); R = stk[--top]) L = R;for (; M = L + R, 1; )if (check(x + M.x, y - M.y)) push(R = M);else {if (judge(x + M.x, R)) break;L = M;}}for (i = 1; i < y; ++i) ret += n / i;return ret*2-srn*srn;
}
int main() {int T;for (scanf("%d", &T); T; --T) {scanf("%lld", &n); putint(S1());puts("");}return 0;
}

转载于:https://www.cnblogs.com/Syameimaru/p/10197934.html

【SP26073】DIVCNT1 - Counting Divisors 题解相关推荐

  1. SPOJ 20713 DIVCNT2 - Counting Divisors (square)

    DIVCNT2 - Counting Divisors (square) #sub-linear #dirichlet-generating-function Let \sigma_0(n)σ​0​​ ...

  2. SP34096 DIVCNTK - Counting Divisors (general)(Min_25筛)

    题面 洛谷 \(\sigma_0(i)\) 表示\(i\) 的约数个数 求\(S_k(n)=\sum_{i=1}^n\sigma_0(i^k)\mod 2^{64}\) 多测,\(T\le10^4,n ...

  3. Counting Divisors HDU - 6069

    设n=p_1^{c_1}p_2^{c_2}...p_m^{c_m}n=p​1​c​1​​​​p​2​c​2​​​​...p​m​c​m​​​​,则d(n^k)=(kc_1+1)(kc_2+1)...( ...

  4. [Spoj]Counting Divisors (cube)

    来自FallDream的博客,未经允许,请勿转载,谢谢. 设d(x)表示x的约数个数,求$\sum_{i=1}^{n}d(i^{3})$ There are 5 Input files. - Inpu ...

  5. 【线性筛】【质因数分解】【约数个数定理】hdu6069 Counting Divisors

    d(x)表示x的约数个数,让你求(l,r<=10^12,r-l<=10^6,k<=10^7) #include<cstdio> using namespace std; ...

  6. 2017 Multi-University Training Contest - Team 4:1003. Counting Divisors(积性函数)

    公式: 上面的a|b表示a是b的约数,ai是每个质因数p分解后的项数 枚举每个x(l<=x<=r) 对于每个数暴力分解小于sqrt(r)的质数就好了, 因为大于sqrt(r)的质数最多只有 ...

  7. 【数学1】基础数学问题 - 题单 - 洛谷

    这里写目录标题 [[数学1]基础数学问题 - 题单 - 洛谷](https://www.luogu.com.cn/training/117) [P1143 进制转换](https://www.luog ...

  8. ICPC程序设计题解书籍系列之九:罗勇军《算法竞赛入门到进阶》

    罗书<算法竞赛入门到进阶>题目一览 第1章 算法竞赛概述 HDU1000 HDU1089-HDU1096 A+B for Input-Output Practice (I)-(VIII)( ...

  9. 浅谈积性函数求前缀和

    转载至https://blog.csdn.net/skywalkert/article/details/50500009 前置技能 积性函数的定义 若f(n)f(n)的定义域为正整数域,值域为复数,即 ...

最新文章

  1. boost之asio同步io使用实例
  2. Spring Data Jpa 实体类自动创建数据库表失败解决
  3. 我ABAP开发生涯中搜集的一些有意思的数据库表
  4. teamviewer无法建立连接原因未知_3389远程无法连接的5种原因分析
  5. 回文判断--链表实现
  6. phpcms能做什么呢?有什么作用呢?
  7. HYSBZ2761 不重复数字【序列处理】(BZOJ2761)
  8. 学python有前途吗-三十岁了,从零开始学python还有前途吗?
  9. Url 助手类(Url Helper)
  10. 飞盘比赛(入门oj Problem 5961)
  11. SysTick系统滴答定时器
  12. 计算机英语pork,[语音]各种肉的英文
  13. android9模拟刘海屏,刘海屏  |  Android 开源项目  |  Android Open Source Project
  14. google、bing ,baidu 等搜索引擎 查询参数
  15. C语言 关系运算符与关系表达式
  16. k宝无法连接计算机,农行K宝无法识别
  17. 设计模式 策略模式 以角色游戏为背景
  18. 【PYTHON】正则匹配的贪婪模式和懒惰模式
  19. 单曲循环 翻译_Loop, repeat and replay “单曲循环” 用哪个词?
  20. EGL, GLX/WGL/AGL 和GL之间的关系

热门文章

  1. elementUI的DatePicker+DateTimePicker组件的自定义日期禁用
  2. CSS布局奇淫巧计之-强大的负边距
  3. Python-100 | 练习题 01 列表推导式
  4. 类的构造器-init和new
  5. 南京市儿童医院用医保身份(医保通道)网上预约挂号以及取号、付费看病流程...
  6. Android获取最新发送短信的基本信息,没有之一
  7. maven常见问题问答 收藏
  8. EM Algorithm
  9. 对lIKE语句的优化
  10. Winform登录窗体登录仿asp.net验证成功后进入主界面