[SPOJ] DIVCNT2 - Counting Divisors (square) (平方的约数个数前缀和容斥卡常)

# 题目
vjudge URL:[Counting Divisors (square) ](https://vjudge.net/problem/SPOJ-DIVCNT2)
Let $\sigma_0(n)$ be the number of positive divisors of $n$.

For example, $\sigma_0(1) = 1$, $\sigma_0(2) = 2$ and $\sigma_0(6) = 4$.

Let $$S_2(n) = \sum _{i=1}^n \sigma_0(i^2).$$

Given $N$, find $S_2(N)$.

###### Input
First line contains $T$ ($1 \le T \le 10000$), the number of test cases.

Each of the next $T$ lines contains a single integer $N$. ($1 \le N \le 10^{12}$)

###### Output
For each number $N$, output a single line containing $S_2(N)$.

###### Example
- Input
>5
1
2
3
10
100
- Output
> 1
4
7
48
1194
- Explanation for Input
$~~~~S_2(3) = \sigma_0(1^2) + \sigma_0(2^2) + \sigma_0(3^2) = 1 + 3 + 3 = 7$

###### Information
>There are 6 Input files.
>- Input #1: $1 \le N \le 10000$, TL = 1s.
>- Input #2: $1 \le T \le 800,\ 1 \le N \le 10^{8}$, TL = 20s.
>- Input #3: $1 \le T \le 200,\ 1 \le N \le 10^{9}$, TL = 20s.
>- Input #4: $1 \le T \le 40,\ 1 \le N \le 10^{10}$, TL = 20s.
>- Input #5: $1 \le T \le 10,\ 1 \le N \le 10^{11}$, TL = 20s.
>- Input #6: $T = 1,\ 1 \le N \le 10^{12}$, TL = 20s.

***My C++ solution runs in 5.3 sec. (total time)*** //呵呵

***Source Limit is 6 KB.****

#### 题目分析
$\sigma_0(n)$表示$n$的约数个数,即求
$$\large \sum_{i=1}^n
\sigma_0(i^2)$$
$1<=n<=10^{12}$

##### 前言
- mdzz，写了1h的Latex公式没保存。。因为机房电脑的烂CPU，本地测这道题的极限数据时崩溃了c
- 这道题是真的卡常，最后点线性筛必须筛到很大，5e7能过，1e7、2e7都不行(可能因为我是大常数选手吧，悄悄打上卡常的FLAG)
- 本地评测炸电脑，真的无语。。
- 刚开始模了1e9+7 WA了好久。。。
- 由于是求函数前缀和，又和杜教筛的题一起做的，就放在杜教筛/莫比乌斯反演里吧
##### 正文
设$\large n=\prod_{i=1}^kp_i^{a_i}$
则$$\large \sigma_0(n^2)=\prod_{i=1}^k(2a_i+1)\\=\sum_{S\in\{1,2,...,k\}}2^{|S|}\cdot\prod_{i\in S}a_i$$
将$\prod_{i\in S}a_i$看作是只由$\color{blue}S内的下标所对应的素数$(乘起来)构成的$n$的约数的**个数**,因为每一个约数都对答案造成了对应的$2^{|S|}$的贡献,那么
$$\large \sigma_0(n^2)=\sum_{d|n}2^{\omega(d)}$$其中$\omega(d)$表示d的质因子的个数(想想)
与此同时，我们又发现$2^{\omega(d)}$实质上是$d$的质因子选或不选的方案数,也就是$n$的**无平方因子的约数**的个数，则
$$\large \sigma_0(n^2)=\sum_{d|n}\sum_{k|d}|\mu(k)|$$因为根据$\mu$函数的定义，只有无平方因子数的函数值才为$1$或$-1$，加上绝对值就相当于统计了个数(有的题解也写的是$\mu(k)^2$,个人认为第一眼看到这个平方会懵一会)
$$\large \sum_{i=1}^n\sigma_0(i^2)=\sum_{i=1}^n\sum_{d|i}\sum_{k|d}|\mu(k)|\\=\sum_{k=1}^n|\mu(k)|\sum_{k|d}\sum_{d|i}1\\=\sum_{k=1}^n|\mu(k)|\sum_{k|d}\lfloor\frac nd\rfloor\\=\sum_{k=1}^n|\mu(k)|\sum_{d=1}^{\lfloor\frac nk\rfloor}\lfloor\frac n{dk}\rfloor\\=\sum_{k=1}^n|\mu(k)|\sum_{d=1}^{\lfloor\frac nk\rfloor}\lfloor\frac {\lfloor\frac nk\rfloor}d\rfloor\\=\sum_{k=1}^n|\mu(k)|\sum_{d=1}^{\lfloor\frac nk\rfloor}\lfloor\frac {\lfloor\frac nk\rfloor}d\rfloor$$
- 先看第二个$\large\sum$，对于某一个$\large{\lfloor\frac nk\rfloor}$的取值，把它记作$N$，就以$N$的范围做整除分块优化，$\large\Theta(\sqrt N)$的时间复杂度，那么外层还有一个求和，于是在外面也套一层整除分块优化，预处理出前$\large n^{\frac 23}$后时间复杂度为$\large\Theta(n^{\frac23})$
- 此处预处理为线性筛，考虑变换，$\large\sum_{i=1}^n\large{\lfloor\frac ni\rfloor}$实际可看作枚举$i$后看$n$以内有多少个数能被$i$整除，这不就是$\large\sum_{i=1}^n\sigma_0(i)$吗？
于是我们只需要筛出约数个数在累加就行了，线性筛时存一下当前数的最小质因子的次数就可以愉快的线性筛了
- 由于在外面一层套上了整除分块优化，则需要求出$\large |\mu(k)|$的前缀和，也就是$n$以内的无平方因子数
- 这里处理无平方因子数时用容斥原理，有
$$\large\sum_{i=1}^n|\mu(i)|=\sum_{i=1}^{\sqrt n}\mu(i)\cdot\lfloor\frac n{i^2}\rfloor$$想想$\mu$函数的定义，这个容斥还是比较好理解的
$\large \Theta(\sqrt n)$可处理出来

综上，各种操作之后把时间复杂度降到了$\large\Theta(n^{\frac 23})$

等等，真的降到了吗？？！看看降到$\large\Theta(n^{\frac 23})$的条件？
- 预处理出前$\large n^{\frac 23}$

然而$1<=n<=10^{12}$(掀桌)
所以只能尽可能的接近，实测$5e7$能过，$2e7$都会TLE

##### AC code

#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
using namespace std;
typedef long long LL;
const int MAXN = 5e7 + 1;//!!!
int Prime[MAXN], mu[MAXN], d[MAXN], Min_a[MAXN], Cnt;
bool IsnotPrime[MAXN];
LL sum_d[MAXN], sum_mu[MAXN];
void init(int n)//线性筛，Min_a[i]存的是i最小质因子的次数
{
    d[1] = mu[1] = 1;
    for(int i = 2; i <= n; ++i)
    {
    
        if(!IsnotPrime[i])
            Prime[++Cnt] = i, mu[i] = -1, d[i] = 2, Min_a[i] = 1;
        for(int j = 1, v; j <= Cnt && i * Prime[j] <= n; ++j)
        {
            v = i * Prime[j];
            IsnotPrime[v] = 1; Min_a[v] = 1;
            if(i % Prime[j] == 0)
            {
                Min_a[v] += Min_a[i];
                mu[v] = 0;
                d[v] = d[i] / Min_a[v] * (Min_a[v] + 1);
                break;
            }
            mu[v] = -mu[i];
            d[v] = d[i]<<1;
        }
    }
    for(int i = 1; i <= n; ++i)
        sum_d[i] = sum_d[i-1] + d[i],
        sum_mu[i] = sum_mu[i-1] + mu[i]*mu[i];
}
inline LL Sum_mu(LL n)//莫比乌斯函数的绝对值的前缀和/[1,n]无平方因子数个数
{
    if(n < MAXN) return sum_mu[n];
    LL ret = 0;
    for(LL i = 1; i*i <= n; ++i)
        ret += mu[i] * (n/(i*i));
    return ret;
}
inline LL Sum_d(LL n) //约数个数前缀和
{
    if(n < MAXN) return sum_d[n];
    LL ret = 0;
    for(LL i = 1, j; i <= n; i=j+1)
    {
        j = n/(n/i);
        ret += (n/i) * (j-i+1);
    }
    return ret;
}

inline LL solve(LL n)
{
    LL ret = 0;
    for(LL i = 1, j; i <= n; i=j+1)
    {
        j = n/(n/i);
        ret += (Sum_mu(j)-Sum_mu(i-1)) * Sum_d(n/i);
    }
    return ret;
}

int main ()
{
    LL T, n;
    scanf("%lld", &T);
    init(T > 800 ? 10000 : MAXN-1); //优化
    while(T--)
    {
        scanf("%lld", &n);
        printf("%lld\n", solve(n));
    }
}

View Code

参见 [传送门:大佬博客](https://www.cnblogs.com/qzqzgfy/p/5600088.html)
...
...
...
...
...

***再次吐槽数学公式的难打***

posted @ 2019-12-14 14:52 _Ark 阅读(315) 评论(1) 编辑收藏举报

会员力量，点亮园子希望

刷新页面返回顶部

[SPOJ] DIVCNT2 - Counting Divisors (square) (平方的约数个数前缀和 容斥 卡常)

公告

[SPOJ] DIVCNT2 - Counting Divisors (square) (平方的约数个数前缀和容斥卡常)