CF1097D 题解

· · 题解

题目链接

本题解已全面更新。

首先有一个重要观察：答案关于 n 是积性的。原因很简单，设 a_k(n) 为答案，就有一个显然的递归关系

a_k(n)=\frac{1}{\sigma_0(n)}\sum_{d|n}a_{k-1}(d)

其中 \sigma_0(n) 表示 n 的约数个数。根据积性函数的基本性质，两个积性函数的 Dirichlet 卷积（这里就是 Dirichlet 前缀和）仍然是积性的；两个积性函数对应点值相乘，也是积性的。

所以只要 a_{k-1}(n) 是积性的，a_k(n) 就是积性的。显然 a_0(n)=n 是积性函数，故答案是积性的。

利用这个性质，就只要对 n 做质因数分解，对每个质数幂求出答案然后相乘即可。也就是说，我们只用考虑算 a_k(p^m)，其中 p 为质数。

此时直接用原本的递归关系来计算，复杂度就大大降低了，足以通过本题。

但是这样的 DP 还是太朴素，不符合我鱼鱼之风格，一定要想办法优化一下。
下文中，我们都考虑只对一个固定的 p 来求答案。所以为了方便设 f_{i,j}=a_i(p^j)，在 f 中不写出 p。

在质数幂上，原本的递归关系就变为（有边界值 f_{0,j}=p^j）：

f_{i,j}=\frac{1}{j+1}\sum_{k=0}^jf_{i-1,k}

这个递推式比较特别，第 i 行的信息可以只从第 i-1 行推过来，这是解迭代列可以做的问题。我们设

F_i(x)=\sum_{j=0}^\infty f_{i,j}x^j

递推式就可以写为（对等式两边分别提取 x^j 项系数，对应的就是原式的等式两边）

(xF_i(x))'=\frac{F_{i-1}(x)}{1-x}

这样似乎不太好看，可以设 G_i(x)=xF_i(x)，就改写为

(1-x)xG_i'(x)=G_{i-1}(x)

但这样还是不便于处理，就要考虑基变换换元了。具体来说，我们希望找到合适的 u（也是一个关于 x 的式子），满足 H_i(u)=G_i(x)，同时有

u \frac{\text d H_i(u)}{\text d u}=H_{i-1}(u)

根据链式法则：

\frac{\text dG_i(x)}{\text dx}= \frac{\text dH_i(u)}{\text d u} \frac{\text du}{\text dx}

在等式两边同乘 (1-x)x，左边就是 G_{i-1}(x)，同时也是 H_{i-1}(u)。这样就得到

H_{i-1}(u)=u \frac{\text dH_i(u)}{\text du}=x(1-x)\frac{\text dH_i(u)}{\text d u} \frac{\text du}{\text dx}

化简一下，就得到关于 u 的微分方程：

u=x(1-x) \frac{\text du}{\text dx}

解得 u=x/(1-x)。根据前面的边界条件可知 F_0(x)=(1-px)^{-1}，也就可以解出

H_0(u)=\frac{u}{1-(p-1)u}=\sum_{i\geq 0}(p-1)^iu^{i+1}

关于 uH_i'(u)=H_{i-1}(u) 这个性质，将等式两边同时提取 u^i 系数就会发现，每次迭代就是在 u^i 系数上除以 i 而已。这样就可以轻易地写出 H_k(u)：

H_k(u)=\sum_{i\geq 0}\frac{(p-1)^i}{(i+1)^k}u^{i+1}

接下来的工作就很简单了，再把 H_k(u) 变回 F_k(x) 即可：

\begin{aligned}f_{k,m}&=[x^m]\frac 1xH_k\left( \frac{x}{1-x}\right) \\ &= \sum_{i\geq 0} \frac{(p-1)^i}{(i+1)^k} [x^{m+1}] \left( \frac{x}{1-x}\right)^i \\ &= \sum_{i=0}^m \binom mi \frac{(p-1)^i}{(i+1)^k}\end{aligned}

直接计算即可，(i+1)^{-k} 可以先用线性筛求幂，然后再线性求逆元，单次时间复杂度为 \Theta(m+\pi(m) \log k)，其中 \pi(m) 表示不超过 m 的质数个数。

最后我们来分析一下总时间复杂度，首先有一个分解质因数的复杂度，一般为 \sqrt n 或 \sqrt[4]{n}。
然后是算 f_{k,m} 的时间复杂度，由于 m 之和有一个显然的上界 \log_2 n，且质因数多于一个时，\pi(m) 之和总会减少。也就是说，这部分的时间复杂度不会超过 \mathcal O(\log_2 n +\pi(\log_2 n)\log k)。

这里用的是朴素分解质因数，总时间复杂可以写为

\mathcal O\left(\sqrt n+ \frac{ \log k}{\log \log n}\log n\right) 
#include<cstdio>
#define N 10003
#define ll long long
#define p 1000000007
using namespace std;

inline int power(int a,int t){
    int res = 1;
    while(t){
        if(t&1) res = (ll)res*a%p;
        a = (ll)a*a%p;
        t >>= 1;
    }
    return res;
}

int fac[N],ifac[N];

void init(int n){
    fac[0] = fac[1] = ifac[0] = ifac[1] = 1;
    for(int i=2;i<=n;++i) fac[i] = (ll)fac[i-1]*i%p;
    ifac[n] = power(fac[n],p-2);
    for(int i=n-1;i;--i) ifac[i] = (ll)ifac[i+1]*(i+1)%p;
}

inline int binom(int n,int m){
    if(n<m) return 0;
    return (ll)fac[n]*ifac[m]%p*ifac[n-m]%p;
}

inline int solve(int n,int m,int q){
    int res = 0,tmp = 1;
    for(int i=0;i<=m;++i){
        res = (res+(ll)binom(m,i)*tmp%p*power(i+1,p-n-1))%p;
        tmp = (ll)tmp*(q-1)%p;
    }
    return res;
}

ll n;
int k,ans = 1;

int main(){
    scanf("%lld%d",&n,&k);
    init(233);
    for(ll i=2;i*i<=n;++i){
        if(n%i!=0) continue;
        int cnt = 0;
        while(n%i==0) n /= i,++cnt;
        ans = (ll)ans*solve(k,cnt,i%p)%p;
    }
    if(n>1) ans = (ll)ans*solve(k,1,n%p)%p;
    printf("%d",ans);
    return 0; 
}