题解 SP20174 【DIVCNT3 - Counting Divisors (cube)】
ywy_c_asm
2018-12-01 15:31:15
首先这题是一道洲阁筛板子题，就当做是一篇洲阁筛讲解吧~~不要跟我提啥min25~~，另外部分内容参考自[__debug大神的博客](http://debug18.com/posts/calculate-the-sum-of-multiplicative-function/)在此%%%OrzstO%%%一下。
洲阁筛和杜教筛一样，是用来在低于$O(n)$的时间内求一类积性函数的前缀和的。然后这些函数$F(x)$必须在当$x$是质数的时候是容易求的，emm……或者应该这样说，当$x$是质数的时候，$F(x)$是关于$x$的一个常数次的多项式。洲阁筛的关键就是利用了在这些质数位置的函数值，然后再利用积性函数的性质钦定这些质数是$mindiv$或者$maxdiv$从而把答案递推出来。
好了，这个题让我们求这个：
$\sum_{i=1}^nσ_0(i^3)$
然后我们不妨以$\sqrt n$为界限把这$n$个$F(x)=σ_0(x^3)$的求和给他分开，就是这$n$个数里的每个数都最多只会有一个大于$\sqrt n$的质因子吧，如果真的存在的话把这个大质数除去剩下的部分显然是不超过$\sqrt n$的，我们不妨在$[1,\sqrt n]$内枚举这个部分，那么由于$F(x)=σ_0(x^3)$是积性函数，当$x\in[1,\sqrt n]$或者$x$存在大于$\sqrt n$的质因子，这些$F(x)$的和可以表示为：
$\sum_{i=1}^nF(i)(1+\sum_{j=\sqrt n+1}^{\lfloor\frac{n}{i}\rfloor}[j\text{为质数}]F(j))……(1)$
然后剩下的数既大于$\sqrt n$，又没有超过$\sqrt n$的质因子，我们把他们的$F(x)$单独求个和：
$\sum_{i=\sqrt n+1}^n[maxdiv(i)<=\sqrt n]F(i)……(2)$
对于①式因为$i<=\sqrt n$我们可以把他们的$F$预处理出来，然后对于每个$i\in[1,\sqrt n]$我们需要求出$[\sqrt n+1,\lfloor\frac n i \rfloor]$内所有质数的$F(x)$之和，在最开头我们说了$F(x)$需要在当$x$是质数的时候是一个关于$x$的常数次多项式，实际上我们发现$F(x)=σ_0(x^3)$当$x$是质数的时候恒等于4……于是我们的问题其实就是求这段区间内的质数个数。
首先这个区间内的质数它必须是大于$\sqrt n$的，这是一个很厉害的条件，那么我们就可以把这个问题等价转化为在$[1,\lfloor\frac n i \rfloor]$内有多少数和$[1,\sqrt n]$内的质数互质，显然其他合数是达不到这个条件的。
于是就可以有这样一个递推的思路（我发现这种把问题递推出来的思路我自己肯定想不出来……），就是设$g(i,j)$为在$[1,j]$内与$P_{1…i}$（$P$是质数数组）的质数都互质的数的个数。我们肯定要考虑从$g(i-1,j)->g(i,j)$的转移，即新加了一个质数$P_i$，显然随着$i$的不断增大满足条件的数会不断变少，这之前$g(i-1,j)$是仅与前$i-1$个质数互质的数，这其中势必会有$P_i$的倍数，显然$P_i$是这些数的$mindiv$，那么我们要把这些数减掉（也就是“筛”），已经确定了他们的$mindiv$是$P_i$，那么剩下的部分显然就在$[1,\lfloor\frac j {P_i} \rfloor]$内，也不能出现小于$P_i$的质因子（等于$P_i$也没事），那么就是$g(i-1,\lfloor\frac j {P_i} \rfloor)$，于是$g$的递推式就是：
$g(i,j)=g(i-1,j)-g(i-1,\lfloor\frac j {P_i}\rfloor)$
注意$g(0,j)$就是没有互质限制的数，显然就是$j$了，这是边界问题。
然后我们来考虑他的优化。先明确一下我们的目标是对每个$i\in [1,\sqrt n]$都求一个$g(m,\lfloor\frac n i\rfloor)$（我们令$P_m$为$<=\sqrt n$的最大的质数，注意是$<=$而不是$<$我被这个坑了好久……），然后我们用递归的那种方式去想，发现$j$一开始就是$n$除一个整数，然后在递归的过程中$j$将会不断的除一些整数，根据等式$\lfloor\frac{\lfloor\frac n {i}\rfloor}{j}\rfloor=\lfloor\frac n {ij}\rfloor$可以发现$j$实际上是$n$不断除一个整数除出来的，我们知道$n$除$[1,n]$内的整数只有$O(\sqrt n)$级别个不同的结果，再多了也就是0，所以我们断言$j$只有$O(\sqrt n)$个，我们可以先把这些$j$求出来用哈希表离散化一下再去递推。
然而这样做复杂度仍然是暴力，但是我们可以接着用一开头就用过的根号思想，就是在$[1,j]$内的数最多只会有一个超过$\sqrt j$的质因子，那么当$P_i^2>j$的时候我们不是要把$P_i$的倍数筛掉吗？显然最多只有一个数会被筛掉，而且$i$越大$P_i$也是越大的，那么当我们在递推的时候发现$P_i^2>j$的时候就在以后不再枚举这个$j$了，而是记录一下最后访问他的位置$i_{last}$与$g(i_{last},j)$，当我们以后要访问一个$g(i,j)$的时候就用$g(i_{last},j)$减去$P[i_{last}+1,i]$内有多少个$<=j$的质数（这个比较坑因为这些质数不一定都在$[1,j]$内，可能就没有这个质数不需要筛掉），这个乱搞一下还是很好求的。
然后递推出这个$g$数组的复杂度据说是$O(\frac{n^{\frac 3 4}}{logn})$的，下面除一个$logn$大概是因为$i$是枚举的质数密度会除个$logn$……算了不口胡了反正我不会证……感性理解一下就好
然后我们就成功地解决了①式。这个②式我们可以采用和上面一样的递推的做法，显然我们可以把这个式子前缀和相减一下，那么就是求一段区间内$maxdiv$不超过$\sqrt n$的所有数的$F(x)$之和，于是我们照葫芦画瓢的设$f(i,j)$为$[1,j]$中$maxdiv$不超过$P_i$的数的$\sum F(x)$，边界就是$f(0,j)=1$显然一个质因子都不能有只有1能够满足。这个的转移方式跟上面的$g$正好相反，我们考虑$f(i-1,j)->f(i,j)$的转移时是考虑在当前质数集合里添加一个$P_i$，我们把$maxdiv=P_i$的所有数都添加进来，但是我们在这里就不能让去掉$P_i$之后剩下的那部分带有$P_i$，因为那个也要从$f(i-1,?)$转移过来，所以我们还要枚举一下$P_i$在这个数里的次数，这一点跟$g$还是有所不同的。于是再根据积性函数的性质（所以求积性函数的和还是相对较容易的……）$f$的递推式就长成这样：
$f(i,j)=f(i-1,j)+\sum_{q=1}^{P_i^q<=j}F(P_i^q)f(i-1,\lfloor\frac{j}{P_i^q}\rfloor)$
然后这个式子应该也可以用和$g$一样的优化？当时我在学洲阁筛的时候说$f$需要把$i$的意义改为后$i$个质数然后就可以用上面那个优化了，现在突然发现好像直接套用$g$的方法好像也是可以的……不过并没有试过qwq……把$i$倒过来做就是让$f(i,j)$变成不超过$\sqrt n$的后$i$个质数，我们在递推的时候只枚举$P_i^2<=j$的$j$，当他第一次被枚举到就把没把他枚举到的质数添加进去（在$[1,j]$中只出现一次），然后复杂度也是$O(\frac{n^{\frac 3 4}}{logn})$的。
然后这题就做完了……
然后你在理论上分析完这题可做，然后一写代码发现越写越恶心……因为全都是一堆细节问题……建议在写代码之前在一张纸上把整个过程一点不漏的写出来。然后下面是一些实现上的问题。
1.我们对每次询问需要的$f$和$g$的$i$都是$m$就是不超过$\sqrt n$的最大质数的下标，那么我们可以先把表打出来，并且$i$显然是可以用滚动数组滚掉的。
2.我们要把$j$用哈希表离散化……不要引入带$log$的二分或$map$……对了不要忘记清空哈希表。
3.在求一段区间内的质数有多少个小于$x$的时候，可以不用二分，因为这都是$\sqrt n$级别的，所以我们在最开始线性筛预处理质数的时候可以对每个数记一下不超过他的最大的质数的下标是啥就可以直接$O(1)$算了。
4.$n$比较小的时候可以直接用前缀和，这样会快一点。
另外不要高估了这题的数据……虽然有$1e4$的询问但是这个复杂度还是可以过的……数据里不会全是极端数据。
上代码~
```cpp
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define int unsigned long long
#define p 19260817
using namespace std;
typedef unsigned long long ll;
namespace ywy{
typedef struct _n{
ll data;int dp;int nxt;
}
node;node memchi[2000001];
int gn=1;int heads[19260818];
inline void insert(ll data,int dp){
memchi[gn].data=data;
memchi[gn].dp=dp;
memchi[gn].nxt=heads[data%p];
heads[data%p]=gn;gn++;
}
inline int get(ll data){
for(register int i=heads[data%p];i;i=memchi[i].nxt){
if(memchi[i].data==data)return(memchi[i].dp);
}return(0);
}
inline ll get(){
ll n=0;char c;while((c=getchar())||23333){
if(c>='0'&&c<='9')break;if(c=='-')goto s;
}n=c-'0';while((c=getchar())||23333){
if(c>='0'&&c<='9')n=n*10+c-'0';else return(n);
}s:while((c=getchar())||23333){
if(c>='0'&&c<='9')n=n*10-c+'0';else return(n);
}return(0);
}
int nearby[400001];
void iprint(ll num){
if(num>=10)iprint(num/10);
putchar(num%10+'0');
}
int prime[400001];
unsigned char bv[400001];
int ptr=1;
int mindiv[400001];
ll sig[400001],sums[400001];
inline void prim(int n){
for(register int i=2;i<=n;i++){
if(!bv[i]){
nearby[i]=ptr;
prime[ptr]=i;ptr++;
mindiv[i]=i;
}for(register int j=1;j<ptr;j++){
int cjr=i*prime[j];
if(cjr>n)break;
mindiv[cjr]=prime[j];
bv[cjr]=1;
if(i%prime[j]==0)break;
}
}for(register int i=1;i<=n;i++){
if(!nearby[i])nearby[i]=nearby[i-1];
ll ans=1;int tmp=i;
while(tmp!=1){
int d=mindiv[tmp];int k=1;
while(mindiv[tmp]==d)k+=3,tmp/=d;ans*=k;
}
sig[i]=ans;
sums[i]=sums[i-1]+ans;
}
}
ll f[2][4000001],g[2][4000001];
inline int erfen(ll n){
int ans=0,l=1,r=ptr-1;
while(l<=r){
int mid=(l+r)>>1;
if((ll)prime[mid]*(ll)prime[mid]<=n)ans=mid,l=mid+1;
else r=mid-1;
}return(ans);
}
ll ints[4000001];
ll n;ll sqr;int sqrp;int pos=0;
int i0[1000001];ll rg[1000001];
inline int how(int l,int r,ll x){
if(l>r)return(0);
if(prime[r]<=x)return(r-l+1);
if(prime[l]>x)return(0);
return(nearby[x]-l+1);
}
inline void getg(){
for(register int i=0;i<pos;i++)g[0][i]=ints[i],i0[i]=0,rg[i]=0;
for(register int i=1;i<=sqrp;i++){
for(register int j=pos-1;~j&&!i0[j];j--){
ll ans=g[(i-1)&1][j];
if(prime[i]>ints[j]/prime[i]){
i0[j]=i;ans-=(ints[j]>=prime[i]);rg[j]=ans;
}else{
int dst=get(ints[j]/prime[i]);
if(i0[dst])ans-=(rg[dst]-how(i0[dst]+1,i-1,ints[j]/prime[i]));
else ans-=g[(i-1)&1][dst];
}g[i&1][j]=ans;
}
}
}
int fi0[4000001];
inline void getf(){
for(register int i=0;i<pos;i++)fi0[i]=0,f[0][i]=f[1][i]=0;
for(register int i=1;i<=sqrp;i++){
for(register int j=pos-1;~j&&prime[sqrp-i+1]<=ints[j]/prime[sqrp-i+1];j--){
ll ans=f[(i-1)&1][j];
if(!fi0[j])fi0[j]=i,ans=1ll+4ll*(ll)how(sqrp-i+2,sqrp,ints[j]);
ll ji=prime[sqrp-i+1],F=4;
while(ji<=ints[j]){
int dst=get(ints[j]/ji);
if(!fi0[dst]){
ans+=F*(1ll+4ll*(ll)how(sqrp-i+2,sqrp,ints[j]/ji));
}
else{
ans+=F*f[(i-1)&1][dst];
}
ji*=prime[sqrp-i+1];F+=3;
}
f[i&1][j]=ans;
}
}
}
void ywymain(){
prim(400000);int t=get();
while(t){
t--;n=get();
if(n<=400000){
iprint(sums[n]);putchar('\n')；continue;
}
sqr=sqrt(n),sqrp=erfen(n);pos=0;
for(register ll i=1;i<=n;){
ints[pos]=n/i;pos++;i=n/(n/i)+1;
}
for(register int i=1;i<=sqr;){
ints[pos]=sqr/i;pos++;i=sqr/(sqr/i)+1;
}
sort(ints,ints+pos);
pos=unique(ints,ints+pos)-ints;
for(register int i=0;i<pos;i++){
insert(ints[i],i);
}
getg();getf();
ll a=f[sqrp&1][get(n)],b=f[sqrp&1][get(sqr)];
if(!fi0[get(n)])a=1ll+4ll*how(1,sqrp,n);
if(!fi0[get(sqr)])b=1ll+4ll*how(1,sqrp,sqr);
ll ans=a-b;
for(register int i=1;i<=sqr;i++){
int ywy=get(n/i);
ll cjr=g[sqrp&1][ywy];
if(i0[ywy])cjr=rg[ywy]-how(i0[ywy]+1,sqrp,n/i);
ans+=((cjr-1ll)*4ll+1ll)*sig[i];
}
iprint(ans);putchar('\n');
for(register int i=0;i<pos;i++){
heads[ints[i]%p]=0;
}gn=1;
}
}
}
signed main(){
ywy::ywymain();return(0);
}
```