题解:P1221 最多因子数

· · 题解

求一个数的因数的个数

由算数基本定理 知,对于任意一个正整数 n,若其有 m 个质因子,则 n 可表示为 \prod_{i=1}^{m}p_i^{k_i}

考虑 n 的每一个质因子 p_i 其指数的取值范围是 [0,k_i],因此满足是该质因子的自然数次幂的 n 的因数的个数就是 k_i+1

m 个质因子相当于把类似的事情分成了 m 步,因此由分步乘法原理,n 的因数的个数为 \prod_{i=1}^{m}(k_i+1)

注:为了方便表述,下面用 r 表示原题中的 u.

暴力求解 时间复杂度 O((r-l)\sqrt{r})

首先容易想到暴力,枚举区间 [l,r] 内的每一个数 n,对 n 做质因子分解,套用上述公式计算每个 n 的因子数,取最大值即可,代码实现如下。

#include<iostream>
using namespace std;
int l,r,ans,cnt;
inline void brute_force(){
    int i,n,p,res,mul;
    for(i=l;i<=r;i++){
        n=i,res=1;
        for(p=2;p*p<=n;p++){
            mul=1;
            while(!(n%p))
                n/=p,mul++;
            res*=mul;
        }if(n!=1) res<<=1;
        if(res>cnt) ans=i,cnt=res;
    }
}int main(){
    cin.tie(0),cout.tie(0);
    ios::sync_with_stdio(false);
    cin>>l>>r,brute_force();
    cout<<"Between "<<l<<" and "<<r<<", "<<ans;
    cout<<" has a maximum of "<<cnt<<" divisors.";
    return 0;
}

正解

上述的解法效率显然十分低下,因此考虑从最本质的定义入手,来求解因数最多的数。

还是考虑一个正整数 n,显然其可表示为 \prod_{i=1}^{m}p_i^{k_i}

那么我们考虑枚举质数 p 与该质数对应的指数 k 来找到其对应的数 n

这样子的话我们只需要先用欧拉筛筛出前 \sqrt{r} 个质数就可以找出区间 [1,r] 内因数最多的数(注意这里的表述是 [1,r] 而不是 [l,r],具体原因后面会进行解释)。

因此筛出质数后,即可使用 dfs 进行搜索,每一层枚举当前质因子的指数即可,部分代码实现如下。

void dfs(int x,int now,int mul){
    int i,m,&p=prm[x];
    if(l<=now&&now<=r){
        if(mul>cnt||mul==cnt&&now<ans)
            ans=now,cnt=mul;
    }if(x==n) return;
    for(i=0,m=1;now*m<=r;i++,m*=p)
        dfs(x+1,now*m,mul*(i+1));
}

然而发现运行效率还是很低,这个时候就要使用搜索算法中的优化技巧——剪枝。

剪枝1

枚举到当前的数 now 时,需要判断区间 [l,r] 内是否存在可以被 now 整除的数。如果不存在,那么就没有继续搜索的必要了。

具体的,我们考虑这样一件事,如果区间 [l,r] 内的一个数 n 可以被 now 整除,那么一定有 \lfloor\frac{n-1}{now}\rfloor\neq\lfloor\frac{n}{now}\rfloor

又因为 l-1\le n-1n<=r 所以有 \lfloor\frac{l-1}{now}\rfloor\neq\lfloor\frac{r}{now}\rfloor

由此我们得出了第一个剪枝方案:若 \lfloor\frac{l-1}{now}\rfloor=\lfloor\frac{r}{now}\rfloor,则直接返回。

剪枝2

假如我们枚举到当前的数 now 算出来 now 的因子数为 mul,又知道当前算出来的因子数最多的数的因子数为 cnt,则可考虑最优性剪枝。

显然,为了不让我们枚举的数 now 超过 r 我们最多还可以对 now 乘上 \lfloor\frac{r}{now}\rfloor

又因为我们当前枚举到的质因子为 p,考虑这个因数的指数 k,最大也只能是 \lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor

由于后面质因子都要比当前的 p 更大,因此 now 最多还能乘上 \lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor 个质数,这样的话当前的 mul 如果继续搜下去的话最大也不会超过 mul\times2^{\lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor},如果这个数都比 cnt 小的话那么就可以直接剪枝。

剪枝3

考虑我们枚举到了质因子 p,发现这个时候 p\times now>r 那么就说明了当前质因子 p 的指数 k 就只能为 0 了(因为如果指数 k>0 的话枚举到下一层的 now 一定会比 r 大)。

而又因为我们是按照从小到大的顺序枚举每一个质因子 p 的,所以如果当前这一层的指数 k 都只能为 0 的话,那么后面的指数也一定只能为 0,因此没有必要再往下枚举,直接返回。

而在具体实现上有个小技巧,由前文知我们在每一层枚举的指数 k 最大为 \lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor,因此直接判断 \lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor 是否为 0 ,如果是则说明指数 k 最大为 0 ,即 p\times now>r,直接返回。

小优化

由上文知我们当前枚举的质因子 p 的指数 k 最大只能为 \lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor 又发现取较小的数作为质因子显然对搜出因数更多的数更有利。

因此我们可以考虑将指数 k\lfloor\log_p{\lfloor\frac{r}{now}\rfloor}\rfloor 向下枚举到 0,这样就能更快地搜出较大的答案,优化后面的剪枝。

部分代码实现如下

inline int lg(int x,int y){
    return log(y)/log(x)+eps;
}inline int qpow(int x,int y){
    int ret=1;
    while(y){
        if(y&1) ret*=x;
        x*=x,y>>=1;
    }return ret;
}void dfs(int x,int now,int mul){
    if(l<=now&&now<=r){
        if(mul>cnt||mul==cnt&&now<ans)
            ans=now,cnt=mul;
    }if(x==n||(l-1)/now==r/now) return;
    int m,&p=prm[x],i=lg(p,r/now);
    if(!i||mul<<i<cnt) return;
    for(m=qpow(p,i);~i;i--,m/=p)
        dfs(x+1,now*m,mul*(i+1));
}

上文中就提到了我们讨论的是处理区间 [1,r] 的方案,为什么不是 [l,r] 呢?

因为在特殊的输入数据下(比如 lr 都是一个大质数)可能就挂掉了,因为一个大质数无法通分解成小的质因子来将其搜索到。

由此,我们还需要结合上文给出的暴力求解方案。

我们考虑先筛出区间 [1,10^7] 以内所有的质数,用这些质数来做 dfs,假如有没有搜索到的数,那么这个数也一定含有一个大于 10^7 的质因子。

不妨假设这个数最大的质因子 p 刚好大 10^7 一点,而 n 的数量级达到 10^9,则 \frac{n}{10^7} <100

所以最坏的情况下(即 n 无法搜到且因数最多时),n 会有 24 个因数。因为在 100 以内因数最多的数为 60,有 12 个因子,再算上 n 的那个最大的质因子,所以 n 最多只有 2\times12=24 个因数。

综上,在 dfs 搜索结束后,如果答案 cnt<=24 则意味着可能会有没搜索到的更优解,在这种情况下 lr 的相差也不会非常大,使用暴力算法再跑一遍即可。

这样,我们就既保证了算法的正确性,也在最大程度上优化了算法的复杂度,完整代码如下。

#include<cmath>
#include<bitset>
#include<iostream>
using namespace std;
constexpr int N=1e7;
constexpr double eps=1e-14;
int prm[N];
bitset<N>vis;
int i,j,l,r,n,ans,cnt;
inline int lg(int x,int y){
    return log(y)/log(x)+eps;
}inline int qpow(int x,int y){
    int ret=1;
    while(y){
        if(y&1) ret*=x;
        x*=x,y>>=1;
    }return ret;
}void dfs(int x,int now,int mul){
    if(l<=now&&now<=r){
        if(mul>cnt||mul==cnt&&now<ans)
            ans=now,cnt=mul;
    }if(x==n||(l-1)/now==r/now) return;
    int m,&p=prm[x],i=lg(p,r/now);
    if(!i||mul<<i<cnt) return;
    for(m=qpow(p,i);~i;i--,m/=p)
        dfs(x+1,now*m,mul*(i+1));
}inline void brute_force(){
    int i,j,n,p,res,mul;
    for(i=r;i>=l;i--){
        n=i,res=1;
        for(j=0;p*(p=prm[j])<=n;j++){
            mul=1;
            while(!(n%p))
                n/=p,mul++;
            res*=mul;
        }if(n!=1) res<<=1;
        if(res>=cnt) ans=i,cnt=res;
    }
}int main(){
    cin.tie(0),cout.tie(0);
    ios::sync_with_stdio(false);
    for(i=2;i<N;i++){
        if(!vis[i]) prm[n++]=i;
        for(j=0;j<n&&i*prm[j]<N;j++){
            vis.set(i*prm[j]);
            if(!(i%prm[j])) break;
        }
    }cin>>l>>r,dfs(0,1,1);
    if(cnt<25) brute_force();
    cout<<"Between "<<l<<" and "<<r<<", "<<ans;
    cout<<" has a maximum of "<<cnt<<" divisors.";
    return 0;
}