[BZOJ 3601]一个人的数论

[BZOJ 3601] 一个人的数论

题意

给定 \(K\) 和 \(n\) 对 \(p_i,r_i\), 令 \(N=\prod p_i^{r_i}\), 求下式对 \(10^9+7\) 取模后的值:
\[ \sum_{k=1,k\perp N}^{N-1} k^K \]
其中 \(p_i,r_i\le 10^9, K\le 100, n\le 1000\).

题解

看见有学弟做了于是也做一做...

推了推发现可能要用Bernoulli数, 他们好神仙啊qaq...

下面是我的画柿子过程, 其中 \(B_i\) 是伯努利数. 感觉这题还是有点意思的.

首先先进行一些套路性的捯饬求和号的过程:
\[ \begin{aligned} \text{Ans}&=\sum_{i=1}^N[i\perp N]i^K \\ &=\sum_{i=1}^N\sum_{d|i,d|N} \mu(d)i^K\\ &=\sum_{d|N}\sum_{k=1}^{\frac N d}\mu(d)(kd)^K \\ &=\sum_{d|N}\mu(d)d^K\sum_{k=1}^{\frac nd}k^K \end{aligned} \]
然后发现最后面是个自然数幂和, 可以尝试用伯努利数怼进去化一化.

PS: 后来问Cage发现其实斯特林数也是能以比较简单的形式代入自然数幂和的, 它的形式是这样的:
\[ n^m=\sum_{k=0}^m \begin{Bmatrix}m \\ k\end{Bmatrix} n^{\underline k} \\ \begin{aligned} \sum_{x=0}^nx^m&=\sum_{x=0}^{n}\sum_{k=0}^m \begin{Bmatrix}m \\ k\end{Bmatrix} x^{\underline k} \\ &=\sum_{k=0}^m\begin{Bmatrix}m \\ k\end{Bmatrix} \sum_{x=0}^n x^{\underline k} \end{aligned} \]
然后内层和式就是混凝土数学里面那个经典的离散微积分(差分/求和)的式子, 也就是:
\[ \Delta(x^\underline n)=nx^\underline{n-1} \\ \sum x^\underline n\delta x=\frac {x^\underline{n+1}}{n+1} \]
所以就有:
\[ \sum_{x=0}^nx^m=\sum_{k=0}^m\begin{Bmatrix}m \\ k\end{Bmatrix} \frac{x^\underline{m+1}}{m+1} \]
然而斯特林数会搞出下降幂来, 在这个题的柿子里并不是很协调的样子, 我们用伯努利数丢进去变换一下和式:
\[ \begin{aligned} \text{Ans}&=\sum_{d|N}\mu(d)d^K \sum_{i=0}^K\frac 1{K+1}B_i{K+1\choose i}\left (\frac N d\right)^{K+1-i}\\ &=\sum_{d|N}\sum_{i=0}^K\mu(d)d^K\frac 1{K+1}B_i{K+1\choose i}\left (\frac N d\right)^{K+1-i}\\ &=\sum_{i=0}^K\sum_{d|N}\mu(d)d^K\frac 1{K+1}B_i{K+1\choose i}\left (\frac N d\right)^{K+1-i}\\ &=\frac 1{K+1}\sum_{i=0}^KB_i{K+1\choose i}\sum_{d|N}\mu(d)d^K\left (\frac N d\right)^{K+1-i}\\ \end{aligned} \]
里层和式里有个 \(d^K\left (\frac N d\right)^{K+1-i}\), 显然它就是 \(N^K\left (\frac N d\right)^{1-i}\). 这样的话 \(N^K\) 就与和式无关了, 丢到外面:
\[ \begin{aligned} \text{Ans}&=\frac 1{K+1}\sum_{i=0}^KB_i{K+1\choose i}\sum_{d|N}\mu(d)N^K\left (\frac N d\right)^{1-i}\\ &=\frac {N^K}{K+1}\sum_{i=0}^KB_i{K+1\choose i}\sum_{d|N}\mu(d)\left (\frac N d\right)^{1-i} \end{aligned} \]
然后根据 \(\mu\) 函数的性质, 只要 \(x\) 中含有完全平方因子则 \(\mu(x)=0\). 于是我们考虑只枚举质因子集合. 然而这次 \(N\) 可能很大没法直接枚举, 尝试DP.

枚举每个质因子 \(p_k\), 如果 \(p_k\not \mid d\) 的话相当于给和式中的每一项都乘上一个 \((p_k^{r_k})^{1-i}\). 否则 \(p_k\mid d\), 相当于和式中的每一项都乘上一个 \(-(p_k^{r_k-1})^{1-i}\) (负号是 \(\mu\) 中来的). 两种方式构造出的值求和就是新的值. 实际上就相当于:
\[ \sum_{d|N}\mu(d)\left (\frac N d\right)^{1-i}=\prod_{k=1}^n \left (\left(p_k^{r_k}\right)^{1-i}-\left (p_k^{r_k-1}\right)^{1-i}\right) \]
于是最后总的答案就是:
\[ \frac {N^K}{K+1}\sum_{i=0}^KB_i{K+1\choose i}\prod_{k=1}^n \left (\left(p_k^{r_k}\right)^{1-i}-\left (p_k^{r_k-1}\right)^{1-i}\right) \]
算上朴素递推Bernoulli数和快速幂, 总时间复杂度 \(O\left((K^2+nK)\log p\right)\). 这题数据数据范围没出满

好像因为只和一个幂次的自然数幂和有关于是直接Lagrange插值插出自然数幂和多项式的系数也是可以的qwq...复杂度一样

然而这题幂次的范围超小甚至可以裸高斯消元解出自然数幂和多项式系数

参考代码

#include <bits/stdc++.h>

const int MAXN=1010;
const int MOD=1e9+7;
const int PHI=MOD-1;

int n;
int K;
int p[MAXN];
int r[MAXN];
int B[MAXN]; // Bernoulli
int inv[MAXN];
int fact[MAXN];

int C(int,int);
int Pow(int,int,int);

int main(){
    scanf("%d%d",&K,&n);
    int N=1;
    for(int i=1;i<=n;i++){
        scanf("%d%d",p+i,r+i);
        N=1ll*N*Pow(p[i],r[i],MOD)%MOD;
    }
    fact[0]=1;
    for(int i=1;i<=K+1;i++)
        fact[i]=1ll*fact[i-1]*i%MOD;
    inv[K+1]=Pow(fact[K+1],MOD-2,MOD);
    for(int i=K+1;i>=1;i--)
        inv[i-1]=1ll*inv[i]*i%MOD;
    B[0]=1;
    for(int i=1;i<=K;i++){
        int sum=0;
        for(int j=0;j<i;j++)
            (sum+=1ll*C(i+1,j)*B[j]%MOD)%=MOD;
        B[i]=1ll*(MOD-sum)*Pow(C(i+1,i),MOD-2,MOD)%MOD;
    }
    int ans=0;
    for(int i=0;i<=K;i++){
        int prod=1;
        for(int j=1;j<=n;j++)
            prod=1ll*prod*(Pow(p[j],1ll*r[j]*(MOD-i)%PHI,MOD)+MOD-Pow(p[j],1ll*(r[j]-1)*(MOD-i)%PHI,MOD))%MOD;
        (ans+=1ll*B[i]*C(K+1,i)%MOD*prod%MOD)%=MOD;
    }
    ans=1ll*ans*Pow(K+1,MOD-2,MOD)%MOD*Pow(N,K,MOD)%MOD;
    printf("%d\n",ans);
    return 0;
}

int C(int n,int m){
    return n<0||m<0||n<m?0:1ll*fact[n]*inv[m]%MOD*inv[n-m]%MOD;
}

inline int Pow(int a,int n,int p){
    int ans=1;
    while(n>0){
        if(n&1)
            ans=1ll*a*ans%p;
        a=1ll*a*a%p;
        n>>=1;
    }
    return ans;
}

转载于:https://www.cnblogs.com/rvalue/p/11032613.html