【2022省选模拟】分解——Min_25筛

好题!可惜没有原题链接

题目描述

在这里插入图片描述
在这里插入图片描述
注:唯一分解中的质因子按从小到大排序

题解

这题的数据范围看起来就很 Min_25,于是我们试看看这个函数是不是积性的——

完蛋,不是。那做不了了。

其实 Min_25 筛的作用绝不仅仅是筛积性函数的前缀和,而是解决一类和数的质因子分解有关的数论DP问题。这题的函数不是积性函数,但仍然可以用 Min_25 筛来做。设 p α p^{\alpha} pα n n n 的最小质因子及其指数,那么由于 φ \varphi φ σ 1 \sigma_1 σ1 都是积性函数,所以有这样的DP式子:
W ( n ) = φ ( p α ) ⋅ W ( n p α ) + p ⋅ σ 1 ( n p α ) W(n)=\varphi(p^{\alpha})\cdot W(\frac{n}{p^{\alpha}})+p\cdot \sigma_1(\frac{n}{p^{\alpha}}) W(n)=φ(pα)W(pαn)+pσ1(pαn)
我们只需要把 Min_25 筛的第二部分(枚举合数最小质因子并计算)做一些改动,改成维护两个值 ( W ( n ) , σ 1 ( n ) ) (W(n),\sigma_1(n)) (W(n),σ1(n)) 即可。

代码

#include<bits/stdc++.h>//JZM yyds!!
#define ll long long
#define uns unsigned
#define END putchar('\n')
#define fi first
#define se second
#define IF (it->fi)
#define IS (it->se)
#define lowbit(x) ((x)&-(x))
using namespace std;
const int MAXN=1000005;
const ll INF=1e18;
inline ll read(){
	ll x=0;bool f=1;char s=getchar();
	while((s<'0'||s>'9')&&s>0){if(s=='-')f^=1;s=getchar();}
	while(s>='0'&&s<='9')x=(x<<1)+(x<<3)+(s^48),s=getchar();
	return f?x:-x;
}
int ptf[50],lpt;
inline void print(ll x,char c='\n'){
	if(x<0)putchar('-'),x=-x;
	ptf[lpt=1]=x%10;
	while(x>9)x/=10,ptf[++lpt]=x%10;
	while(lpt)putchar(ptf[lpt--]^48);
	if(c>0)putchar(c);
}

const ll MOD=1e9+7,iv2=(MOD+1)>>1;
ll n;
bool nop[MAXN];
int pr[MAXN/9],le,B;
ll pre[MAXN/9],w[MAXN],g[MAXN],g_[MAXN];
int m,id[MAXN],di[MAXN];
inline int ID(ll x){return x<=B?id[x]:di[n/x];}
inline void prepare_min25(){
	B=sqrt(n),nop[0]=nop[1]=1,le=0;
	for(int a=2;a<=B;a++){
		if(!nop[a])pr[++le]=a,pre[le]=(pre[le-1]+a)%MOD;
		for(int i=1,u;i<=le&&(u=pr[i]*a)<=B;i++){
			nop[u]=1;
			if(a%pr[i]==0)break;
		}
	}
	for(ll i=1,j,la;i<=n;i=la+1){
		j=n/i,la=n/j;
		if(j<=B)id[j]=++m;
		else di[la]=++m;
		w[m]=j,g[m]=(j+2)%MOD*((j-1ll)%MOD)%MOD*iv2%MOD;
		g_[m]=j-1;
	}
	for(int i=1;i<=le;i++){
		ll v=1ll*pr[i]*pr[i];
		for(int j=1;j<=m;j++){
			if(w[j]<v)break;
			g[j]=(g[j]+MOD-(g[ID(w[j]/pr[i])]-pre[i-1]+MOD)*pr[i]%MOD)%MOD;
			g_[j]-=g_[ID(w[j]/pr[i])]-i+1;
		}
	}
}
#define pll pair<ll,ll>
inline pll operator+=(pll&a,const pll&b){
	a.fi+=b.fi,a.se+=b.se;
	if(a.fi>=MOD)a.fi-=MOD;
	if(a.se>=MOD)a.se-=MOD;
	return a;
}
inline pll F(ll x,int i){
	if(x<2||i>le+1||x<pr[i])return pll(0,0);
	pll res=pll((g[ID(x)]+MOD-pre[i-1])%MOD,233);
	res.se=(res.fi+g_[ID(x)]-i+1)%MOD;
	for(;i<=le;i++){
		ll e=pr[i],cg=pr[i]+1;
		if(e*pr[i]>x)break;
		for(;e*pr[i]<=x;cg=(cg*pr[i]+1)%MOD,e*=pr[i]){
			pll d=F(x/e,i+1);
			res+=pll(((e-e/pr[i])%MOD*d.fi+pr[i]*d.se)%MOD,cg*d.se%MOD);
			res+=pll(pr[i],(cg*pr[i]+1)%MOD);
		}
	}return res;
}
int main()
{
	freopen("b.in","r",stdin);
	freopen("b.out","w",stdout);
	n=read();
	prepare_min25();
	print((F(n,1).fi+1)%MOD);
	return 0;
}
  • 0
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值