【51nod1847】奇怪的数学题

最新推荐文章于 2021-07-24 11:00:15 发布

weixin_30322405

最新推荐文章于 2021-07-24 11:00:15 发布

阅读量218

点赞数

CC 4.0 BY-SA版权

文章标签： ui 数据结构与算法

原文链接：http://www.cnblogs.com/Judge/p/10620755.html

本文介绍了一种求解两个数的次大约数之和的方法，通过杜教筛、第二类斯特林数等数学工具，实现了高效求解。文章详细解释了如何将原问题转化为数论求和，并给出了具体实现代码。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

就当我是 A 了此题吧...

首先预备知识有点多（因为题目要处理的东西都挺毒瘤）：

杜教筛运用（当然你可以用其他筛？）
第二类斯特林数相关定理
下降阶乘幂相关定理
min25 筛运用

~~好了可以关掉本题解了~~

咳咳上面我除了杜教筛都不会，于是熬了好久才抄做掉的

题意

我们首先考虑题目让我们求什么：

\[ANS=\sum_{i=1}^n\sum_{j=1}^n sgcd(i,j)\]

上面的 \(sgcd(i,j)\) 表示 \(gcd(i,j)\) 的次大约数

emmm？有问题么？虽然题目让我们求两个数的次大公约数，但是这个次大公约数必然能被他们的 gcd 整除，那么这个次大公约数也就是 gcd 的次大约数了

noteskey

Part 1

我们考虑用 f(i) 表示一个数的次大约数：

\[ANS=\sum_{i=1}^n\sum_{j=1}^n f(gcd(i,j))\]

然后这里有 gcd 啊，立马开始化式子：

\[ANS=\sum_{d=1}^nf(d)\sum_{i=1}^{n/d}\sum_{j=1}^{n/d} [gcd(i,j)=1]\]

这时候我们就有两个选择了，要么反演变成 μ ，要么转成 φ ，个人感觉可能 φ 方便点（毕竟大家都写的 φ ），但是 μ 也没有差到哪里去

求 φ

\[ANS=\sum_{d=1}^nf(d)\sum_{i=1}^{n/d}\sum_{j=1}^{n/d} \sum_[gcd(i,j)=1]\]

\[ANS=\sum_{d=1}^nf(d)~ (2·\sum_{i=1}^{n/d}φ(i)-1)\]

然后我们发现这里可以数论分块，上杜教筛！

求 μ

\[ANS=\sum_{d=1}^nf(d)\sum_{i=1}^{n/d}\sum_{j=1}^{n/d} \sum_[gcd(i,j)=1]\]

\[ANS=\sum_{d=1}^nf(d) \sum_{k=1}^{n/d}μ(k) \sum_{i=1}^{n/(kd)}\sum_{j=1}^{n/(kd)} \]

\[ANS=\sum_{d=1}^nf(d) \sum_{k=1}^{n/d}μ(k) {n \over kd}^{~2} \]

然后我们发现还是可以数论分块，上杜教筛！

~~然后这篇题解到此结束，完结撒花~~

有什么不对的样子...emmm 哦！前面的 f 还没说怎么求啊！

Part 2

我们重新考虑 f 是个啥：

如果说 x 是个质数，那么 f(x) 显然是 1 ，然后我们考虑一个合数的 f 值其实就是其本身除去最小值因子

总的来说，一个数 x 的 f 值就是 x 除去其最小质因子（当然 1 是有特例，算 0 的）

那么按照套路我们肯定是要求 f 的前缀和的（不然怎么数论分块啊！）

所以我们肯定也是要求出类似 \(\sum_{i=1}^n i^k\) 这样的东西（我们考虑除去一个最小值因子之后的贡献时肯定要用到的嘛）

然后我们看看这玩意儿怎么求先

k 次幂和的快速计算

可能看到这玩意儿你会想到拉格朗日插值...

但是很遗憾这道题的模数很不友好（就是那个自然溢出的玩意儿），我们只能另寻他路

所以首先我们要用的就是第二类斯特林数的性质了：

\[n^k = \sum_{i=1}^k \left\{ ~ ^k_i \right\}i! (~ ^n_i)\]

那么前缀和就是：

\[\sum_{i=1}^n i^k = \sum _{j=1}^n \sum_{i=1}^k \left\{ ~ ^k_i \right\}i! (~ ^j_i )\]

\[\sum_{i=1}^n i^k =\sum_{i=1}^k \left\{ ~ ^k_i \right\} \sum _{j=1}^n i! ( ~^j_i)\]

然后看到后面的东西其实是一个下降阶乘幂

\[\sum_{i=1}^n i^k =\sum_{i=1}^k \left\{ ~ ^k_i \right\} \sum _{j=1}^n j^\underline{i}\]

分析这个下降阶乘幂的前缀和：

\[\sum _{j=1}^n j^\underline{i}= {(n+1)^\underline{i+1}\over i+1}\]

至于证明...要用积分（我们看到这玩意儿挺像多项式求积分的形式的，就是感觉底数 n+1 有点不对劲，不过毕竟这不是简单的多项式而是下降阶乘幂，所以就默认它存在吧）

那么原来的式子带进去：

\[\sum_{i=1}^n i^k =\sum_{i=1}^k \left\{ ~ ^k_i \right\} {(n+1)^\underline{i+1}\over i+1}\]

这样的话我们就可以预处理斯特林数然后 \(O(k)\) 跑出 \(i^k\) 的前缀和了（注意这里的 k 上限 50 ，所以不会出事）

那么我们考虑质数合数贡献分开来算，就是，就是说我们只要让合数的 f 值前缀和加上质数个数（每个 f 贡献为 1 ，相当于求质数个数）就是要求的 f 前缀和了

质数个数可能还比较好做，但是还要考虑合数的次大约数前缀和怎么求：

我们考虑之前说的一个数除去它的最小质因子，这有点像欧拉筛哈，但是欧拉筛是线性的，对于这样的问题无能为力

于是我们考虑 min25 筛中其实也有每个数被最小质因子筛去的这一性质，于是我们欢乐地套上 min25 筛就过了此题

code

这里用的是第一种化简（也就是求 φ 的），但其实这题关键在于求出 f 的前缀和所以求 φ 还是 μ 无所谓了啦

//by Judge
#include<bits/stdc++.h>
#define fp(i,a,b) for(register int i=(a),I=(b)+1;i<I;++i)
#define ui unsigned int
using namespace std;
const int M=1e6+3;
typedef ui arr[M<<1];
ui n,m,mm,K,cnt,cntval,sq; ui ans,S[60][60];
arr v,p,loc1,loc2,g,rec,phi,pcnt,pw,G; map<ui,ui> mp;
inline ui qpow(ui x,int p){ ui s=1;
    for(;p;p>>=1,x=x*x)
        if(p&1) s=s*x; return s;
}
inline ui Pos(ui x){ //判断返回两种 id 中哪种的 val 
    return x<=sq?loc1[x]:loc2[n/x];
}
inline ui f(ui x){  // 求 f 前缀和 
    return G[Pos(x)]+pcnt[Pos(x)];
}
inline void prep(int n){ v[1]=phi[1]=1; //预处理各种变量 
    for(int i=2;i<=n;++i){ if(!v[i]) p[++cnt]=i,phi[i]=i-1;
        for(int j=1;j<=cnt&&i*p[j]<=n;++j){ v[i*p[j]]=1;
            if(!(i%p[j])){phi[i*p[j]]=phi[i]*p[j];break;}
            phi[i*p[j]]=phi[i]*(p[j]-1);
        } phi[i]+=phi[i-1]; //线性筛 phi  
    } S[0][0]=1;
    fp(i,1,K){ S[i][1]=1;
        fp(j,2,i) S[i][j]=S[i-1][j]*(ui)j+S[i-1][j-1]; //计算第二类斯特林数 
    }
    fp(i,1,cnt) pw[i]=qpow(p[i],K); //计算 p[i] 的 K 次幂 
}
inline void init_loc(ui n){ cntval=0; //预处理出要处理的 n/i 的 val 
    for(ui l=1,r;l<=n;l=r+1)
        r=n/(n/l),rec[++cntval]=n/l;
    reverse(rec+1,rec+1+cntval);
    fp(i,1,cntval)
        if(rec[i]<=sq) loc1[rec[i]]=i;
        else loc2[n/rec[i]]=i;
}
inline ui sum(ui l,ui r){ //计算 l~ r 的和 
    ui tmp1=l+r,tmp2=r-l+1;
    if(tmp1&1) tmp2>>=1;
    else tmp1>>=1;
    return tmp1*tmp2;
}
inline ui get_phi(ui n){ if(n<=mm) return phi[n]; //杜教筛求 phi 前缀和 
    if(mp[n]) return mp[n]; ui ans=sum(1,n);
    for(ui l=2,r;l<=n;l=r+1) r=n/(n/l),
        ans-=get_phi(n/l)*(r-l+1); return mp[n]=ans;
}
inline ui calc_sumk(long long n){  //计算 n^K 的前缀和 
    ui ans=0,tmpval,tmp;
    fp(i,1,K){ // K 次幂 ，循环 K 次 
        tmpval=1,tmp=(n-i+1)%(i+1);
        for(long long j=n-i+1;j<=n+1;++j,++tmp){
            if(tmp>=i+1) tmp-=i+1;
            if(!tmp) tmpval*=(ui)j/(i+1);
            else tmpval*=(ui)j; //计算下降幂 
        }
        ans+=S[K][i]*tmpval;
    } return ans;
}
inline void get_g(){ // g 存 rec[i] 范围内 质数 ^ k 次幂 的前缀和 ， pcnt 存 rec[i] 内质数个数（也就是 所有质数 f 的前缀和） 
    fp(i,2,cntval) g[i]=calc_sumk(rec[i])-1,pcnt[i]=rec[i]-1;
    fp(j,1,m) for(int i=cntval;i&&1ll*p[j]*p[j]<=rec[i];--i)
        g[i]-=pw[j]*(g[Pos(rec[i]/p[j])]-g[Pos(p[j]-1)]),  // g 当前为最小质因子大于 p[j] 或者为质数 的数的 k 次幂前缀和 
        pcnt[i]-=pcnt[Pos(rec[i]/p[j])]-pcnt[Pos(p[j]-1)],
        G[i]+=g[Pos(rec[i]/p[j])]-g[Pos(p[j]-1)]; // G 中记录合数的 f 前缀和，因为是从大到小筛所以这里的 g 表示的是 g[rec[i]/p[j]][j-1] 
}
int main(){ cin>>n>>K;
    prep(mm=pow(n,2/3.0)),sq=sqrt(n);
    m=upper_bound(p+1,p+1+cnt,sq)-1-p;
    init_loc(n),get_g(),ans=0; ui pre=0,now;
    for(ui l=2,r;l<=n;l=r+1,pre=now) r=n/(n/l),
        now=f(l),ans+=(get_phi(n/l)*2-1)*(now-pre);
    return !printf("%u\n",ans);
}

转载于:https://www.cnblogs.com/Judge/p/10620755.html