2019 徐州网络赛H min25+阶乘贡献

一开始发现f(xy)=f(x)+f(y),这不是积性函数能求?

思路参考:https://www.cnblogs.com/hua-dong/p/11508790.html

质数p对sum{i=1~n}f(i)的贡献是sum{k=1,p^k<=n}(n/p^k) 

就像阶乘算多少个0,

5倍数多1个0

25的倍数多两个0 

但算25的倍数时5的倍数已经对25的倍数做出过贡献,所以减掉之后25的贡献为1

 

这里f(i)也是一样,p^1的贡献是1 p^2贡献是2 

但算p的倍数时p的倍数已经对p^2的倍数做出过贡献,所以减掉之后p^2的贡献也是1

 

也可以这样想,假设n=16 ,p=2 

2对至少有1个质因子2的数贡献 =n/2

2对至少为2个质因子2的数贡献 =n/4

所以2对仅有有1个质因子2的数 (f(2),f(6),f(10),f(14)) 的贡献 =(n/2-n/4)*1=4

4的贡献 f(4)+f(12)  =(n/4-n/8)*2 

8的贡献 f(8) =(n/8-n/16)*3  

16的贡献f(16)=n/16*4;

4个加起来 发现:ans= n/2+n/4+n/8+n/16 

对于i*f(i)  

2的贡献 f(2)+f(6)+f(10)+f(14) =  1*( (2+4+6+8+...16)-(4+8+12+16) )  

=1*(   2*(1+2+3.....8)  -  4*(1+2+3+4)  ) =1*(2sum(8)-4sum(4)  )

4的贡献  f(4)+f(12)=2*(4sum(4)-8sum(2))

8的贡献           f(8)=3*(8sum(4)-16sum(1))

16的贡献       f(16)=4*(16sum(1))

4个加起来 发现:ans= 2*sum(n/2)+4*sum(n/4)+8*sum(n/8)+16*sum(n/16);

 

所以我们答案= n以内每个质数p的贡献和

 \sum_{i=1}^{n}(n+1-i)f(i)=\sum_{k=1}^{p_k<=n}\sum_{t=1}^{p_k^t<=n}(n+1)\frac{n}{p_k^t}-p_k^tsum(\frac{n}{p_k^t})

 

对于p<=sqrt(n)的质数,我们枚举质数k,质数的倍数t,暴力把p,p^2,p^3的贡献把柿子求出来

对于p>sqrt(n)的质数,n/p<sqrt(n),显然有t=1,那么Sum=\sum_{k=?}^{p_k<=n}(n+1)\frac{n}{p_k}-p_ksum(\frac{n}{p_k})

对于前半部分,我们设g(x)=\sum_{i=1}^{x}[\text{i is prime}]   

众所周知,在min25筛中,对于所有可能的n/i=m的值 g(m)是可以套模板维护出来的(整数分块)

因此 [w[m],w[m+1] ]区间函数和也是很快能维护出来的=g(m)-g(m+1)  

比如n=10  那么g(10) g(5) g(3) g(2) g(1)能求出来  w[1]=10  w[2]=5  w[3]=3....

比如n=10  那么[6,10]=g[1]-g[2]  ,  [4,5]=g(2)-g(3) 的和都能维护出来 

那么我们枚举所有(n/i)=t<sqrt(n)的值 对于每个t,看有多少个n/质数在[t,t+1]内,这区间内每个质数的贡献是t =n/质数

因此ans+=(n+1)*(g[t]-g[t+1])*t

对于后半部分:我们设g(x)=\sum_{i=1}^{x}i[\text{i is prime}]  

同理min25 可求g(m) 

因此ans-=  (g[t]-g[t+1])*sum(t)=(1+t)*t/2 * (g[t]-g[t+1])

 

436ms

//min25筛求f前缀和 
//要求:积性函数+f(p)是多项式(或者是完全积性函数的和 )+ f(p^c)非常好求
//51nod 1239 欧拉函数前缀和 
#include<stdio.h>
#include<math.h>
#include<iostream>
#include<string.h>
using namespace std;
#define ll long long

const ll mod = 998244353;
const ll inv2=(1+mod)/2;
ll qmod(ll a,ll b){ll res=1;while(b){if(b&1)res=res*a%mod;a=a*a%mod;b=b>>1;}return res;}

namespace Min25{

    const int N = 1e6+5;
    ll Sqr,isp[N],pri[N],sump[N],tot=0;//质数上限
    
    void init(int n){//线筛预处理质数 与质数的前缀和(f(p)的前缀和)
        isp[1]=1;
        for (ll i=2;i<=n;++i){
            if (!isp[i]) pri[++tot]=i,sump[tot]=1ll*(sump[tot-1]+i)%mod;
            for (int j=1;i*pri[j]<=n;++j){
                isp[i*pri[j]]=1;
                if (i%pri[j]==0) break;
            }
        }
    }

    ll w[N],id1[N],id2[N],m=0,g[N],h[N]; 

    void getg0(ll n){//对每个n/i 求 g(n/i,0),h(n/i,0);  
        //定义: g(x,p)=sum{i=1~x}[i的最小质因子>pri[p] or i是质数]g0(i) ,要求g0(p)是完全积性函数,不是的话要分开凑。 
        //g(n/i,0)代表 所有数当成质数带入  本来g0(p)=? 仅对质数有效
        //g[1]=g(n/1,0)....g[m]=g(n/m,0) m是整数分块的上限
        memset(h,0,sizeof(h));
        m=0;Sqr=sqrt(n);
        for (ll i=1,j;i<=n;i=j+1){
            j=n/(n/i);w[++m]=n/i;//整数分块 ; 下标映射 ;
            if (w[m]<=Sqr) id1[w[m]]=m;
            else id2[n/w[m]]=m;
 
            //合数也当成质数求和,忽略h0(1) h0(i)=1 g0(i)=i  需要改
            h[m]=(w[m]-1+mod)%mod;
            g[m]=((w[m]+2)%mod)*((w[m]-1)%mod)%mod;
            if (g[m]&1) g[m]+=mod;g[m]/=2;
        }
    }
    void getgp(ll n){//对每个n/i 求g(n/i,|P|) h(n/i,|P|)
 
        // 定义:g(x,p)=sum{i=1~x}[i的最小质因子>p or i是质数]g0(i)  ,其中g0(i)=i; 维护时省略第二维   f(p)=g0(p)-h0(p)凑完全积性 
        // 定义:h(x,p)=sum{i=1~x}[i的最小质因子>p or i是质数]h0(i)  ,其中h0(i)=1; 维护时省略第二维
        // 所以我们这个函数求的是  g(x,|P|)=sum{i=1~x}[i是质数]g0(i)  
        // g[1]=g(n/1,|P|) g[2]=g(n/2,|P|) ...g[m]=g(n/m,|P|) m是整数分块的上限
        for (int j=1;j<=tot;++j)
            for (int i=1;i<=m&&1ll*pri[j]*pri[j]<=w[i];++i){
                int k=(w[i]/pri[j]<=Sqr)?id1[w[i]/pri[j]]:id2[n/(w[i]/pri[j])];//不用改
 
                g[i]=(g[i]-1ll*pri[j]*(g[k]-sump[j-1])%mod)%mod;g[i]=(g[i]+mod)%mod;//需要改 g(n,j)= g(n,j-1) - g0(pri[j])* [ g(k,j-1)-sum_{i=1~j-1}g0(pri[i]) ];
                h[i]=(h[i]-h[k]+j-1);h[i]=(h[i]+mod)%mod;                           //需要改 同上
            }
        //g[i]-g[i+1] 代表[n/(i+1),n/i]区间内的质数p的g0(p)的和
    }

    int fpc(int p,int c){//f(p^c)  需要改
        return (qmod(p,c)-qmod(p,c-1)+mod)%mod;
    }

    int S(ll x,int y,ll n){//定义:S(x,y)=sum{i=1~x}[i的最小质因子>=pri[y] ]*f(i); 求积性函数前缀和ans=S(n,1,n)+f(1) 

        //在此之前,f(n/i,|P|) (n/i)<=m 的各项已经用getgp()维护出来了
        if (x<=1||pri[y]>x) return 0;
        int k=(x<=Sqr)?id1[x]:id2[n/x];
        ll res=(1ll*g[k]-h[k]-sump[y-1]+y-1);(res=res+mod)%mod;//  res= g(n,|P|) - sum_{i=1~j}f(Pi)  需要改

        //if (y==1) res+=2;//特判 f(2)=f(p)+2 其他题可删

        for (int i=y;i<=tot&&1ll*pri[i]*pri[i]<=x;++i){//S的递推式 不改
            ll p1=pri[i],p2=1ll*pri[i]*pri[i];
            for (int e=1;p2<=x;++e,p1=p2,p2*=pri[i])
                (res+=(1ll*S(x/p1,i+1,n)*fpc(pri[i],e)%mod+fpc(pri[i],(e+1))))%=mod;//pri[i]^e 是f(pri[i]^e)的值,需要改
        }
        return res%mod;
    }

}

int main(){
    ll n;
    using namespace Min25;
    scanf("%lld",&n);
    Sqr=sqrt(n);init(Sqr);
    
    getg0(n);getgp(n);
    /*
    for(int i=1;i<=n;i++){
        printf("%d %d %d\n",n/i,i,g[i]);
    }*/
    ll ans=0;
    for(int i=1;i<=tot;i++){//p<sqrt(n)的贡献
        for(ll j=pri[i];j<=n;j*=pri[i]){
            ll t=(n/j)%mod;
            (ans+=t*((n+1)%mod)%mod )%=mod;
            (ans-=(j%mod)*(mod+(1+t)%mod*t%mod*inv2%mod)%mod)%=mod;
            ans=(ans+mod)%mod;
        }
    }
    cout<<ans<<endl;
    for(ll k=1;k<Sqr;k++){//p>sqrt(n) 枚举n/p=k的值
        (ans+=((n+1)%mod)*k%mod*(h[k]-h[k+1]+mod)%mod)%=mod;
        (ans-=mod+(1+k)*k%mod*inv2%mod*(g[k]-g[k+1]+mod)%mod)%=mod;
        ans=(ans+mod)%mod;
    }
    cout<<(ans+mod)%mod<<endl;
    //system("pause");
    return 0;
}

 

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值