题解 函数

传送门

除了 \(1e13\) 的数据范围以外都很 min_25 筛
然而我 min_25 筛的板子忘了

那么是 powerful number 的板子题
\(f(x)\) 的定义题面已经给出
我们构造一个 \(g(x)=x^k\)，尝试找到一个 \(g*h=f\)
发现 \(f(1)=g(1)h(1)=1\)，那么 \(h(1)=1\)
发现 \(f(p)=g(1)h(p)+g(p)h(1)=1=g(p)h(1)\)
那么 \(h(p)=0\)
然后可以得到 \(h(p^i)=p^k-p^{2k}\)，证明暂时不会

然后推一下式子：

\[\begin{aligned} \sum\limits_{i=1}^nf(i)&=\sum\limits_{i=1}^n\sum\limits_{d\mid i}g(i)h(\frac{i}{d})\\ &=\sum\limits_{i=1}^nh(i)\sum\limits_{j=1}^{\lfloor \frac{n}{j}\rfloor}g(j) \end{aligned}\]

然后前面那个东西只在根号个位置有值

• 关于自然数幂和：
  见这里的总结
  一个他没有提到的做法是利用第二类斯特林数
  发现 \(n^k=\sum\limits_{i=0}^k\begin{Bmatrix}k\\i\end{Bmatrix}i!\binom{n}{i}\)
  然后外面套个 \(\sum\) 可以拆开变成
  \[\sum\limits_{i=0}^k\begin{Bmatrix}k\\i\end{Bmatrix}i!\binom{n+1}{j+1} \]

于是就可以做了

• 关于爆搜 powerful number：
  如果这样写：

  void dfs(int u, ll now, ll val) {
  	if (u>pcnt) return ;
  	if (now*pri[u]*pri[u]>n) return ;
  	dfs(u+1, now, val);
  	now*=pri[u]; val=val*h[u]%mod;
  	while (1) {
  		if ((now*=pri[u])>n) break;
  		ans=(ans+val%mod*G(n/now))%mod;
  		dfs(u+1, now, val);
  	}
  }
  

  貌似是 \(O(n)\) 的，原因未详细分析
  需要这样写：

  ll dfs(ll n, int pos, ll val) {
  	ll ans=0;
  	for (int i=pos; i<=pcnt; ++i) {
  		ll t=n/pri[i]/pri[i];
  		if (!t) break;
  		for (; t; t/=pri[i]) ans=(ans+dfs(t, i+1, val*h[i]%mod))%mod;
  	}
  	return (ans+val*G(n))%mod;
  }
  

整体复杂度 \(O(n\sqrt n)\)

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define INF 0x3f3f3f3f
#define N 100000010
#define ll long long
//#define int long long

char buf[1<<21], *p1=buf, *p2=buf;
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf, 1, 1<<21, stdin)), p1==p2?EOF:*p1++)
inline ll read() {
	ll ans=0, f=1; char c=getchar();
	while (!isdigit(c)) {if (c=='-') f=-f; c=getchar();}
	while (isdigit(c)) {ans=(ans<<3)+(ans<<1)+(c^48); c=getchar();}
	return ans*f;
}

ll n, k;
const ll mod=1e9+7;
inline void md(ll& a, ll b) {a+=b; a=a>=mod?a-mod:a;}
inline ll qpow(ll a, ll b) {ll ans=1; for (; b; a=a*a%mod,b>>=1) if (b&1) ans=ans*a%mod; return ans;}

namespace force{
	ll ans;
	// bool npri[N];
	bitset<N> npri;
	int pri[5761460], f[N], pcnt;
	void solve() {
		f[1]=1;
		// cout<<double(sizeof(pri)+sizeof(f)+sizeof(npri))/1000/1000<<endl;
		for (int i=2; i<=n; ++i) {
			if (!npri[i]) pri[++pcnt]=i, f[i]=qpow(i, k);
			for (int j=1,x; j<=pcnt&&1ll*i*pri[j]<=n; ++j) {
				npri[x=i*pri[j]]=1;
				if (!(i%pri[j])) {f[x]=f[i]; break;}
				else f[x]=1ll*f[i]*f[pri[j]]%mod;
			}
		}
		for (int i=1; i<=n; ++i) md(ans, f[i]);
		printf("%lld\n", ans);
	}
}

namespace task{
	bitset<N> npri;
	int pri[5761460], pcnt;
	ll st[21][21], fac[23], inv[23], h[N], sqr, ans;
	inline ll C(int n, int k) {return fac[n]*inv[k]%mod*inv[n-k]%mod;}
	inline ll qpow(ll a, ll b) {ll ans=1; for (; b; a=a*a%mod,b>>=1) if (b&1) ans=ans*a%mod; return ans;}
	inline ll G(ll n) {n%=mod; ll ans=0, t=n+1; for (int i=1; i<=k; ++i) t=t*(n-i+1)%mod, ans=(ans+st[k][i]*fac[i]%mod*inv[i+1]%mod*t)%mod; return ans;}
	// void dfs(int u, ll now, ll val) {
	// 	if (u>pcnt) return ;
	// 	if (now*pri[u]*pri[u]>n) return ;
	// 	dfs(u+1, now, val);
	// 	now*=pri[u]; val=val*h[u]%mod;
	// 	while (1) {
	// 		if ((now*=pri[u])>n) break;
	// 		// cout<<"now: "<<now<<endl;
	// 		ans=(ans+val%mod*G(n/now))%mod;
	// 		// cout<<"val: "<<(val*G(n/now)%mod+mod)%mod<<endl;
	// 		dfs(u+1, now, val);
	// 	}
	// }
	ll dfs(ll n, int pos, ll val) {
		ll ans=0;
		for (int i=pos; i<=pcnt; ++i) {
			ll t=n/pri[i]/pri[i];
			if (!t) break;
			for (; t; t/=pri[i]) ans=(ans+dfs(t, i+1, val*h[i]%mod))%mod;
		}
		return (ans+val*G(n))%mod;
	}
	void solve() {
		sqr=sqrt(n);
		st[0][0]=1; fac[0]=fac[1]=1; inv[0]=inv[1]=1;
		for (int i=2; i<=k+1; ++i) fac[i]=fac[i-1]*i%mod;
		for (int i=2; i<=k+1; ++i) inv[i]=(mod-mod/i)*inv[mod%i]%mod;
		for (int i=2; i<=k+1; ++i) inv[i]=inv[i-1]*inv[i]%mod;
		for (int i=1; i<=k; ++i) for (int j=1; j<=i; ++j) st[i][j]=(st[i-1][j-1]+j*st[i-1][j])%mod;
		for (int i=2; i<=sqr; ++i) {
			if (!npri[i]) pri[++pcnt]=i, h[pcnt]=(qpow(i, k)-qpow(1ll*i*i%mod, k))%mod; //, cout<<(h[i]+mod)%mod<<endl;
			for (int j=1,x; j<=pcnt&&1ll*i*pri[j]<=sqr; ++j) {
				npri[x=i*pri[j]]=1;
				if (!(i%pri[j])) break;
			}
		}
		// for (int i=1; i<=n; ++i) {
		// 	ll tem=0;
		// 	for (int j=1; j<=i; ++j) tem=(tem+qpow(j, k))%mod;
		// 	cout<<tem<<' ';
		// } cout<<endl;
		// for (int i=1; i<=n; ++i) cout<<G(i)<<' '; cout<<endl;
		printf("%lld\n", (dfs(n, 1, 1)%mod+mod)%mod);
	}
}

signed main()
{
	// freopen("function.in", "r", stdin);
	// freopen("function.out", "w", stdout);
	
	n=read(); k=read();
	// force::solve();
	task::solve();

	return 0;
}