[LOJ#6053]奇怪的函数

最新推荐文章于 2023-11-21 12:28:31 发布

alan_cty

最新推荐文章于 2023-11-21 12:28:31 发布

阅读量635

点赞数

CC 4.0 BY-SA版权

分类专栏： min-25筛数论文章标签： Min-25筛数论 LOJ6053 奇怪的函数

本文链接：https://blog.youkuaiyun.com/alan_cty/article/details/80050163

数论同时被 2 个专栏收录

48 篇文章

订阅专栏

min-25筛

2 篇文章

订阅专栏

本文详细介绍了Min_25筛算法的应用场景及其具体实现过程。针对特定数学函数f(x)，通过高效的筛法计算∑ni=1f(i)并对1e9+7取模，适用于n≤1e10的情况。文章提供了完整的C++代码实现，通过质数筛选及分段处理实现了快速求解。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

Description

给出一个函数f(x)满足：
f(1)=1,
f(p^c)=p xor c,p是质数
f(ab)=f(a)*f(b),(a,b)=1
求 $\sum_{i=1}^{n}f(i)$ 对1e9+7取模
n<=1e10

Solution

min_25筛例题
不会可以看一下神仙zzq的博客
我们只需要知道所有质数的f(x)之和，显然指数都为1，那么直接用质数之和-质数个数即可。
注意2是唯一的偶质数，随手特判一下就好了。

Code

#include <cmath>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define fd(i,a,b) for(int i=a;i>=b;i--)
using namespace std;

typedef long long ll;

const int N=2e5+5,Mo=1e9+7;

int pwr(int x,int y) {
    int z=1;
    for(;y;y>>=1,x=(ll)x*x%Mo)
        if (y&1) z=(ll)z*x%Mo;
    return z;
}

int p[N],g[N],h[N],s[N],id1[N],id2[N],m,inv2,tot,cnt;
bool bz[N];
ll n,w[N];

void get_prime(int n) {
    fo(i,2,n) {
        if (!bz[i]) p[++tot]=i;
        fo(j,1,tot) {
            int k=i*p[j];if (k>n) break;
            bz[k]=1;if (!(i%p[j])) break;
        }
    }
}

int ID(ll x) {return x<=m?id1[x]:id2[n/x];}

int F(ll x,int j) {
    if (x<=1||p[j]>x) return 0;
    int k=ID(x),ans=(g[k]-s[j-1]-(h[k]-(j-1)))%Mo;
    if (j==1) ans+=2;
    fo(i,j,tot) {
        if ((ll)p[i]*p[i]>x) break;
        ll t1=p[i],t2=(ll)p[i]*p[i];
        for(int e=1;t2<=x;t1=t2,t2=t1*p[i],e++)
            (ans+=(ll)F(x/t1,i+1)*(p[i]^e)%Mo+(p[i]^(e+1)))%=Mo;
    }
    return ans;
}

int main() {
    scanf("%lld",&n);
    m=sqrt(n);get_prime(m);
    inv2=pwr(2,Mo-2);
    for(ll l=1,r=0;l<=n;l=r+1) {
        r=n/(n/l);
        w[++cnt]=n/l;
        int res=w[cnt]%Mo;
        g[cnt]=(ll)res*(res+1)%Mo*inv2%Mo;g[cnt]--;
        h[cnt]=(w[cnt]-1)%Mo;
        if (w[cnt]<=m) id1[w[cnt]]=cnt;
        else id2[n/w[cnt]]=cnt;
    }
    fo(i,1,tot) s[i]=(s[i-1]+p[i])%Mo;
    fo(j,1,tot) 
        fo(i,1,cnt) {
            if ((ll)p[j]*p[j]>w[i]) break;
            int k=ID(w[i]/p[j]);
            (g[i]-=(ll)p[j]*(g[k]-s[j-1])%Mo)%=Mo;(g[i]+=Mo)%=Mo;
            (h[i]-=h[k]-(j-1))%=Mo;(h[i]+=Mo)%=Mo;
        }
    printf("%d\n",F(n,1)+1);
    return 0;   
}