6306. Sequence

22 篇文章 0 订阅
3 篇文章 0 订阅

Description

Input

Output

Sample Input

2 5 4
 

Sample Output

29

Data Constraint

Solution

我们知道f(x)是一个积性函数,因此我们只需要求出f(p^c)的值就可以推出所有的f()。

如何求出所有的f(p^c)?

对于质数 p,令 d 为满足 p x |B 的最大的 x,

那么 f(p^c ) = { p^d · (c − d + 1)^n + ∑(i=0,i<=d-1) p^i · [(c − i + 1)^n − (c − i)^n ] (c > d)

                       ∑(i=0,i<=c) p^i · [(c − i + 1)^n − (c − i)^n ] (c ≤ d)

这相当于是枚举A[]和B的gcd为p^i次方,当i<=d时我们可以A[]的每个数都取p^i,p^i+1...p^c的数,这样每个数的指数都有(c-i+1)种选发,但是要考虑每个数的指数都没有取到i的情况,这是我们要减去p^i+1,p^i+2...p^c的方案,即每个数有(c-i)种取值。

对于i>d时,因为与B求gcd后p的指数依然为d,所以直接贡献为p^d,这是A[]可以随便选d~c的任意一个数。

至于怎么推,我们设s[ i ]表示 i 除去它最小的质因子的数。例如s[100]=100/2/2=25,s[9]=1。

1.那么对于i是质数时

      s[ i ]=1

2.i不是质数且线筛时 i % p[ j ]不为0,即p[ j ]不是i的质因子,那么p[ j ]就是i*p[ j ]的最小质因子,故

     s[ i*p[ j ] ]=i

3.i不是质数且p[ j ]是 i 的最小质因子,那么s的值不变。

    s[ i*p[ j ] ]=s[ i ]

之后我们对于每个已经求出来的f(p^c)。

1.p[ j ]不是i的最小质因子,即i与p[ j ]互质

      f[ i*p[ j ] ]=f[ i ]*f[ p[ j ] ]

2.p[ j ]是i的最小质因子,那么这时i与p[ j ]不互质,这时我们要将所有i的最小质因子除去,再和所有最小质因子的若干次方相乘即可。

  f[ i*p[ j ] ]=f[ s[ i ] ]*f[ i/s[ i ]*p[ j ] ]

 至此,此题得到解决。

Code

#include<cstdio> 
#include<algorithm>
#include<cstring>
#include<cstdlib>
#define I int
#define F(i,a,b) for(I i=a;i<=b;i++)
#define Fd(i,a,b) for(I i=a;i>=b;i--)
#define mem(a,b) memset(a,b,sizeof(a))
#define N 20000010
#define M 998244353
#define ll long long
#define P(x) freopen(x".in","r",stdin);freopen(x".out","w",stdout)
#define rt return
using namespace std;
void rd(ll &x){
	x=0;char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
}
ll n,m,f[N],g[40],B,x,y,d,c,ans=1;
I bz[N],s[N],p[10000000];
ll ksm(ll x,ll k){
	if(k==0) rt 1;
	if(k==1) rt x;
	ll st=ksm(x,k/2);st=(st*st)%M;
	if(k&1) rt (st*x)%M;
	rt st;
}
I main(){
	P("sequence");
	rd(n);rd(m);rd(B);
	n=n%(M-1);
	f[0]=f[1]=1;
	F(i,1,30) g[i]=ksm(i,n);
	F(i,2,m){
		if(!bz[i]){
			p[++p[0]]=i;
			s[i]=1;
			x=B;d=0;
			while(x%i==0){x/=i;d++;}
			x=i;c=0;
			while(x<=m){c++;x=x*i;}
			x=1;
			F(j,1,c){
				x*=i;
				if(j<=d){
					y=1;
					F(k,0,j){
						f[x]=(f[x]+y*((g[j-k+1]-g[j-k]+M)%M)%M)%M;
						y*=i;
					}
				}
				else{
					y=1;
					F(k,0,d-1){
						f[x]=(f[x]+y*((g[j-k+1]-g[j-k]+M)%M)%M)%M;
						y*=i;
					}
					f[x]=(f[x]+(ksm(i,d)*g[j-d+1])%M)%M;
				}
			}
		}
		ans=(ans+f[i])%M;
		F(j,1,p[0]){
			if(i*p[j]>m) break;
			bz[i*p[j]]=1;
			if(i%p[j]==0){
				s[i*p[j]]=s[i];
				f[i*p[j]]=(f[s[i]]*f[i/s[i]*p[j]])%M;
				break;
			}
			s[i*p[j]]=i;
			f[i*p[j]]=(f[i]*f[p[j]])%M;
		}
	}
//	F(i,1,m) ans=(ans+f[i])%M;
	printf("%lld\n",ans);
	rt 0;
}


作者:zsjzliziyang 
QQ:1634151125 
转载及修改请注明 
本文地址:

  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

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

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

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值