自动搬运

查看原文

来自洛谷，原作者为

	云浅知处 喵

搬运于2025-08-24 22:55:29，当前版本为作者最后更新于2024-02-19 20:21:45，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

---

对于 $n=\prod p_i^{c_i}$，有：

$$\sigma_0^t(n^s)=\prod(1+s\times c_i)^{t}$$

然后设：

$$f(n)=\sum_{d\mid n} \mu\left(\dfrac{n}{d}\right)\sigma_0^t(n^s) $$

所以有：

$$\sigma_0^t(n^s)=\sum_{d\mid n}f(d)$$

由于 $\sigma_0^t$ 是积性函数，不妨先假设 $n=p^c$。则：

$$\sigma_0^t(p^{cs})=\sum_{i=0}^{c}f(p^i)$$

也就是说对于 $n=p^c$，$f(p^c)$ 是 $\sigma_0^t$ 关于 $c$ 的高维差分。

类似的，设：

$$f(n)=\sum_{d\mid n} g(d)$$

那个 $g$ 是 $f$ 关于 $c$ 的差分，即 $g$ 是 $\sigma_0^t$ 的二阶高维差分。具体地就有：

$$g(p^c)=(1+c_is)^t-[c_i\ge 1]2(1+(c_i-1)s)^t+[c_i\ge 2] (1+(c_i-2)s)^t $$

那么 $g$ 是积性函数，就有：

$$g(n)=\prod((1+c_is)^t-[c_i\ge 1]2(1+(c_i-1)s)^t+[c_i\ge 2] (1+(c_i-2)s)^t) $$

矩阵可以拆成：

$$\begin{aligned} a_{i,j}&=\sum_{d\mid \gcd(i,j)}g(d)=\sum_{d=1}^{n}[d\mid i][d\mid j]g(d)\\\\ b_{i,j}&=[j\mid i]g(j)\\\\ c_{i,j}&=[i\mid j]\\\\ A&=B*C \end{aligned} $$

$*$ 是矩阵乘法。$B,C$ 都是上、下三角矩阵，于是

• $\text{det}(B)=\prod_{i=1}^n b_{i,i}=\prod_{i=1}^ng(i)$
• $\text{det}(C)=\prod_{i=1}^nc_{i,i}=1$

故答案为：

$$\det(A)=\text{det} (B)\times \det(C)=\prod_{m=1}^{n}\prod((1+c_is)^t-[c_i\ge 1]2(1+(c_i-1)s)^t+[c_i\ge 2] (1+(c_i-2)s)^t) $$

依次考虑每个质数的贡献。

对于 $p\le \sqrt{n}$ 的部分暴力枚举 $p$ 和 $c$，暴力求快速幂，复杂度 $\mathcal{O}(\sqrt{n}\log t)$。

对于 $p>\sqrt{n}$，可知 $c\le 1$，当 $c=0$ 时贡献为 $1$ 不用考虑，当 $c=1$ 时 $g(p)=(1+s)^t-2$，故只需要筛素数个数即可，在本题范围内，使用 min25 筛即可。

#include<bits/stdc++.h>

#define ll long long
#define mk make_pair
#define fi first
#define se second

using namespace std;

const int mod=1e9+7;
int ksm(int x,ll y,int p=mod){
	int ans=1;y%=(p-1);
	for(int i=y;i;i>>=1,x=1ll*x*x%p)if(i&1)ans=1ll*ans*x%p;
	return ans%p;
}
int inv(int x,int p=mod){return ksm(x,p-2,p)%p;}
mt19937 rnd(time(0));
int randint(int l,int r){return rnd()%(r-l+1)+l;}
void add(int &x,int v){x+=v;if(x>=mod)x-=mod;}
void Mod(int &x){if(x>=mod)x-=mod;}
int cmod(int x){if(x>=mod)x-=mod;return x;}

const int N=1e6+5;
const int B=320000;
ll n;
int cp,pri[N],s,t;
bool vis[N];
void Euler(int V){
	for(int i=2;i<=V;i++){
		if(!vis[i])pri[++cp]=i;
		for(int j=1;j<=cp&&pri[j]*i<=V;j++){
			vis[i*pri[j]]=1;
			if(i%pri[j]==0)break;
		}
	}
}
ll f[N<<1],num[N<<1];
int id1[N],id2[N];
int getid(ll x){
	if(x>B)return id2[n/x];
	return id1[x];
}
void calcpri(){
	Euler(B);
	int ncnt=0;
	for(ll l=1,r=0;l<=n;l=r+1){
		ll val=(n/l);
		if(val<=B)num[id1[val]=++ncnt]=val;
		else num[id2[n/val]=++ncnt]=val;
		r=n/(n/l);r=min(r,n);
	}
	for(int i=1;i<=ncnt;i++)f[i]=num[i]-1;
	for(int i=1;i<=cp;i++){
		for(int j=1;j<=ncnt&&num[j]>=1ll*pri[i]*pri[i];j++){
			int k=getid(num[j]/pri[i]);
			f[j]-=f[k]-(i-1);
		}
	}
}

int g[100];
void solve(){
	cin>>n>>s>>t;
	if(s==0)return puts(n==1?"1":"0"),void();

	int ans=1;calcpri();
	g[1]=cmod(ksm(s+1,t)+mod-2);
	for(int i=2;i<=40;i++){
		add(g[i],ksm(1ll*i*s%mod+1,t));
		add(g[i],mod-2ll*ksm(1ll*(i-1)*s%mod+1,t)%mod);
		add(g[i],ksm(1ll*(i-2)*s%mod+1,t));
	}
	int lim=n/B;
	for(int i=1;i<=cp&&n/pri[i]>lim;i++){
		ll val=pri[i];int k=1;
		while(val<=n){
			ll cc=n/val-n/(val*pri[i]);
			ans=1ll*ans*ksm(g[k],cc)%mod,val*=pri[i],k++;
		}
	}
	ll cc=0;
	// floor(n/p) >= i <==> i*p <= n <==> p <= floor(n/i)
	// floor(n/p) == i <==> floor(n/(i+1)) < p <= floor(n/i)
	for(int i=1;i<=lim;i++){ // floor(n/p) = i
		ll cnt=f[getid(n/i)]-f[getid(n/(i+1))];
		cc+=1ll*i*cnt;
	}
	ans=1ll*ans*ksm(g[1],cc)%mod;
	cout<<ans<<endl;
}

signed main(void){

	solve();

	return 0;
}