题解 CF235E 【Number Challenge】

· · 题解

CF235E Number Challenge解题报告：

更好的阅读体验

题意

分析

第一眼就看出来好像见过这道题，然后想起来它是P4619 [SDOI2018]旧试题的弱化版到底是CCF抄CF的还是CF抄CCF的呢：

我们先证明一个毒瘤式子d(x\cdot y\cdot y)=\sum_{i\mid x}\sum_{j\mid y}\sum_{k\mid z}[\gcd(i,j,k)==1]

当x=\prod_{i=1}^s p_i^{k_i}因为d(x)=\prod_{i=1}^s(k_i+1)为x的因数个数，因此考虑将x\cdot y\cdot z的每个因数对应一组互质的因数

对于每个质数p\mid x\cdot y\cdot z，且x中p的幂次为a，y中p的幂次为b，z中p的幂次为c，则x\cdot y\cdot z中p的幂次为a+b+c。

根据d的积性与\gcd(p^{a+b+c},\frac{x\cdot y\cdot z}{p^{a+b+c}})==1，得：d(x\cdot y\cdot z)=d(\frac{x\cdot y\cdot z}{p^{a+b+c}})\cdot d(p^{a+b+c})=d(\frac{x\cdot y\cdot z}{p^{a+b+c}})\cdot(a+b+c+1)，即p对左式的贡献为a+b+c+1

令i\mid x,j\mid y,k\mid z且\gcd(i,j,k)=1，则i中p的幂次u不超过a，j中p的幂次v不超过b，k中p的幂次w不超过c，且u，v，w中起码有两个为0

即每个p可以对右式同样产生a+b+c+1的贡献，因此对于每个x\cdot y\cdot z的质因子p，其会对左右两式均产生a+b+c+1即p的指数+1的贡献

因此我们可以把原式化为：\sum_{i=1}\sum_{j=1}^B\sum_{k=1}^C\sum_{x\mid i}\sum_{y\mid j}\sum_{z\mid k}[\gcd(x,y,z)==1]

=\sum_{i=1}\sum_{j=1}^B\sum_{k=1}^C\sum_{x\mid i}\sum_{y\mid j}\sum_{z\mid k}[\gcd(x,y)==1][\gcd(y,z)==1][\gcd(x,z)==1]

改变枚举顺序可得：\sum_{x=1}^A\lfloor\frac{A}{x}\rfloor\sum_{y=1}^B\lfloor\frac{B}{y}\rfloor\sum_{z=1}^C\lfloor\frac{C}{z}\rfloor[\gcd(x,y)==1][\gcd(y,z)==1][\gcd(x,z)==1]

使用莫比乌斯反演得原式可得：\sum_{x=1}^A\lfloor\frac{A}{x}\rfloor\sum_{y=1}^B\lfloor\frac{B}{y}\rfloor\sum_{z=1}^C\lfloor\frac{C}{z}\rfloor\sum_{u\mid\gcd(x,y)}\mu(u)\sum_{v\mid\gcd(y,z)}\mu(v)\sum_{w\mid\gcd(x,z)}\mu(w)

继续交换枚举顺序得：\sum_{u=1}^{min(A,B)}\mu(u)\sum_{v=1}^{min(B,C)}\mu(v)\sum_{w=1}^{min(A,C)}\mu(w)\sum_{i=1}^{\lfloor\frac{A}{lcm(u,v)}\rfloor}\lfloor\frac{A}{lcm(u,v)\cdot i}\rfloor\sum_{j=1}^{\lfloor\frac{B}{lcm(v,w)}\rfloor}\lfloor\frac{B}{lcm(v,w)\cdot j}\rfloor\sum_{k=1}^{\lfloor\frac{C}{lcm(u,w)}\rfloor}\lfloor\frac{C}{lcm(u,w)\cdot k}\rfloor

由引理\lfloor\frac{a}{b\cdot c}\rfloor=\lfloor\frac{\lfloor\frac{a}{b}\rfloor}{c}\rfloor得：\sum_{i=1}^{\lfloor\frac{x}{y}\rfloor}\lfloor\frac{x}{y\cdot i}\rfloor=\sum_{i=1}^{\lfloor\frac{x}{y}\rfloor}\lfloor\frac{\lfloor\frac{x}{y}\rfloor}{i}\rfloor，原式化为：\sum_{u=1}^{min(A,B)}\mu(u)\sum_{v=1}^{min(B,C)}\mu(v)\sum_{w=1}^{min(A,C)}\mu(w)\sum_{i=1}^{\lfloor\frac{A}{lcm(u,v)}\rfloor}\lfloor\frac{\lfloor\frac{A}{lcm(u,v)}\rfloor}{i}\rfloor\sum_{j=1}^{\lfloor\frac{B}{lcm(v,w)}\rfloor}\lfloor\frac{\lfloor\frac{B}{lcm(v,w)}\rfloor}{j}\rfloor\sum_{k=1}^{\lfloor\frac{C}{lcm(u,w)}\rfloor}\lfloor\frac{\lfloor\frac{C}{lcm(u,w)}\rfloor}{k}\rfloor

设f(x)=\sum_{i=1}^x\lfloor\frac{x}{i}\rfloor，考虑如何预处理f(x)：对于每一个i\leqslant x，x都会得到i的贡献\lfloor\frac{x}{i}\rfloor，即在[1,x]的范围内i的倍数个数。

因此f(x)=\sigma(x)=\sum_{d\mid x}d即约数和函数，这个函数通过可以枚举每个数的倍数在O(n\log n)的时间内预处理。

原式可化为：\sum_{u=1}^{min(A,B)}\sum_{v=1}^{min(B,C)}\sum_{w=1}^{min(A,C)}\mu(u)\cdot\mu(v)\cdot\mu(w)\cdot\sigma(\frac{A}{lcm(u,v)})\cdot\sigma(\frac{B}{lcm(v,w)})\cdot\sigma(\frac{C}{lcm(u,w)})

但是最后的式子还是O(n^3)的，就比暴力快一点点。

首先，我们可以把所有\mu(x)=0的x删去，这种x有不少，因此可以快很多。

其次，因为\sigma(0)=0，因此所有的lcm(u,v)>A,lcm(v,w)>B,lcm(u,w)>C可以剪掉。（其实不剪掉也没有关系，对答案没有影响，因此我们将A,B,C统一改为max(A,B,C)）

为了达到更低的复杂度，我们可以使用一种毒瘤算法：三元环计数。感性理解一下，枚举三个数u,v,w，且满足lcm(u,v)\leqslant A,lcm(v,w)\leqslant B,lcm(u,w)\leqslant C就相当于在一张只要任意两个点满足其最小公倍数小于等于max(A,B,C)边会有一条边的图中枚举三元环。

三元环计数算法可以在O(m\sqrt{m})的时间内求出所有三元环，具体见nekko的博客：不常用的黑科技——「三元环」

时间复杂度：O(m\sqrt{m})。

代码

直接拿的旧试题的代码，加了亿点常数优化，凑合着看吧。

#include<stdio.h>
#include<vector>
#include<string.h>
#define int long long
using namespace std;
const int maxn=1000005,maxm=1000005,mod=1073741824;
int i,j,k,a,b,c,e,mx,mn,T,cnt,tot,ans;
int lst[maxn],d[maxn],p[maxn],mu[maxn],ok[maxn],ord[maxn],deg[maxn],f[maxn],from[maxm],to[maxm],lcm[maxm],mrk[maxn];
vector<int>v[maxn],w[maxn];
inline int min(int a,int b){
    return a<b? a:b;
}
inline int max(int a,int b){
    return a>b? a:b;
}
inline int gcd(int a,int b){
    return b==0? a:gcd(b,a%b);
}
signed main(){
    p[1]=mu[1]=1;
    for(int i=2;i<maxn;i++){
        if(p[i]==0)
            lst[++cnt]=i,mu[i]=-1;
        for(int j=1;j<=cnt;j++){
            if(i*lst[j]>=maxn)
                break;
            p[i*lst[j]]=1;
            if(i%lst[j]==0){
                mu[i*lst[j]]=0;
                break;
            }
            mu[i*lst[j]]=-mu[i];
        }
    }
    for(int i=1;i<maxn;i++)
        if(mu[i]!=0)
            ok[++tot]=i,ord[i]=tot;
    for(int i=1;i<maxn;i++){
        for(int j=1;i*j<maxn;j++)
            d[i*j]++;
        f[i]=(f[i-1]+d[i])%mod;
    }
    scanf("%lld%lld%lld",&a,&b,&c);
    mn=min(min(a,b),c),mx=max(max(a,b),c);
    e=ans=0;
    for(int i=1;i<=tot;i++){
        if(ok[i]>mx)
            break;
        for(int j=1;j<=tot;j++){
            if(ok[i]*ok[j]>mx)
                break;
            if(mu[ok[i]*ok[j]]==0)
                continue;
                for(int k=j+1;k<=tot;k++){
                if(ok[i]*ok[j]*ok[k]>mx)
                    break;
                if(mu[ok[i]*ok[k]]==0||gcd(ok[j],ok[k])>1)
                    continue;
                from[++e]=ord[ok[i]*ok[j]],to[e]=ord[ok[i]*ok[k]],lcm[e]=ok[i]*ok[j]*ok[k];
                deg[from[e]]++,deg[to[e]]++;
            } 
        }
    }
    for(int i=1;i<=tot;i++){
        if(ok[i]>mn)
            break;
        ans+=mu[ok[i]]*mu[ok[i]]*mu[ok[i]]*f[a/ok[i]]*f[b/ok[i]]*f[c/ok[i]];
    }
    for(int i=1;i<=e;i++){
        v[from[i]].push_back(to[i]),w[from[i]].push_back(lcm[i]);
        v[to[i]].push_back(from[i]),w[to[i]].push_back(lcm[i]);
    }
    for(int i=1;i<=tot;i++){
        if(ok[i]>min(a,b))
            break;
        for(int j=0;j<v[i].size();j++){
            int x=ok[i],y=ok[v[i][j]],z=w[i][j];
            ans=(ans+mu[x]*mu[y]*mu[y]*f[a/z]*f[b/z]*f[c/y]%mod+mod)%mod;
            ans=(ans+mu[x]*mu[x]*mu[y]*f[a/x]*f[b/z]*f[c/z]%mod+mod)%mod;
            ans=(ans+mu[x]*mu[x]*mu[y]*f[a/z]*f[b/x]*f[c/z]%mod+mod)%mod;
        }
    }
    memset(v,0,sizeof(v));
    memset(w,0,sizeof(w));
    for(int i=1;i<=e;i++){
        if(deg[from[i]]>=deg[to[i]])
            v[from[i]].push_back(to[i]),w[from[i]].push_back(lcm[i]);
        else v[to[i]].push_back(from[i]),w[to[i]].push_back(lcm[i]);
    }
    for(int i=1;i<=tot;i++){
        if(ok[i]>mx)
            break;
        for(int j=0;j<v[i].size();j++)
            mrk[v[i][j]]=w[i][j];
        for(int j=0;j<v[i].size();j++){
            int x=v[i][j];
            for(int k=0;k<v[x].size();k++){
                int y=v[x][k],p=mrk[y],q=w[i][j],r=w[x][k];
                if(mrk[y]==0)
                    continue;
                int st1,st2,st3,st4,st5,st6;
                st1=f[a/p]*f[b/q]*f[c/r]%mod;
                st2=f[a/p]*f[b/r]*f[c/q]%mod;
                st3=f[a/q]*f[b/p]*f[c/r]%mod;
                st4=f[a/q]*f[b/r]*f[c/p]%mod;
                st5=f[a/r]*f[b/p]*f[c/q]%mod;
                st6=f[a/r]*f[b/q]*f[c/p]%mod;
                ans=(ans+mu[ok[i]]*mu[ok[x]]*mu[ok[y]]*(st1+st2+st3+st4+st5+st6)%mod+mod)%mod;
            }
        }
        for(int j=0;j<v[i].size();j++)
            mrk[v[i][j]]=0;
    }
    printf("%lld\n",(ans+mod)%mod);
    return 0;
}