LOJ#6053 简单的函数

题面

题目描述

某一天，你发现了一个神奇的函数(f(x))，它满足很多神奇的性质：

1. (f(1)=1)。

2. (f(p^c)=p oplus c)（(p) 为质数，(oplus)表示异或）。

3. (f(ab)=f(a)cdot f(b))（(a) 与 (b) 互质）。

你看到这个函数之后十分高兴，于是就想要求出 (sumlimits_{i=1}^n f(i))。

由于这个数比较大，你只需要输出 (sumlimits_{i=1}^n f(i) mod (10^9+7))。

输入格式

一行一个整数 (n)。

输出格式

一行一个整数 (sumlimits_{i=1}^n f(i) mod 1000000007)。

样例

样例输入 1

6

样例输出 1

16

样例输入 2

233333

样例输出 2

179004642

样例输入3

9876543210

样例输出3

895670833

数据范围与提示

对于(30\%)的数据，(n leq 100)。

对于(60\%)的数据，(n leq 10^6)。

对于(100\%)的数据，(1 leq n leq 10^{10})。

题目分析

要求(f(p^c)=poplus c)，

对除(2)以外的质数(p)而言，等价于求(f(p)=p-1)。

我们可以用该函数来求(g)，

然而，由于该函数并不是一个完全积性函数，

我们需要把它拆成两个部分，分为(f_1(p)=p,f_2(p)=1)，

并预处理出对应的(g(n,j),h(n,j))。

此时，(S(n,j))的初值为

[g(n,|P|)-h(n,|P|)-sum_{i=1}^{j-1}f(p) ]

需要注意的是，若(j=1)，(f(2))会被计算在答案之中，

但是在上面预处理的时候我们算的是(f(2)=1)，需额外加(2)。

代码实现

#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<cstdio>
#include<iomanip>
#include<cstdlib>
#define MAXN 0x7fffffff
typedef long long LL;
const int N=250005,mod=1e9+7;
const int inv2=5e8+4;
using namespace std;
inline LL Getint(){register LL x=0,f=1;register char ch=getchar();while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}while(isdigit(ch)){x=x*10+ch-'0';ch=getchar();}return x*f;}
bool vis[N];
int prime[N],tot;
LL n,sqr,w[N],h[N],g[N],sp[N];
int id1[N],id2[N],m;
void Pre(int n){
	vis[1]=1;
	for(int i=2;i<=n;i++){
		if(!vis[i])prime[++tot]=i,sp[tot]=sp[tot-1]+i;
		for(int j=1;j<=tot&&1ll*i*prime[j]<=n;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j]==0)break;
		}
	}
}
int S(LL x,int y){
	if(x<=1||prime[y]>x)return 0;
	int k=(x<=sqr)?id1[x]:id2[n/x],ret=(g[k]-sp[y-1]-h[k]+y-1)%mod;
	if(y==1)ret+=2;
	for(int i=y;i<=tot&&1ll*prime[i]*prime[i]<=x;i++){
		LL t1=prime[i],t2=1ll*prime[i]*prime[i];
		for(int e=1;t2<=x;++e,t1=t2,t2*=prime[i])
			(ret+=1ll*S(x/t1,i+1)*(prime[i]^e)%mod+((prime[i]^(e+1))%mod))%=mod;
	}
	return ret;
}
int main(){
	n=Getint(),sqr=sqrt(n);
	Pre(sqr);
	for(LL i=1,j;i<=n;i=j+1){
		j=n/(n/i),w[++m]=n/i;
		g[m]=(w[m]%mod)*((w[m]+1)%mod)%mod*inv2%mod-1;
		h[m]=(w[m]-1)%mod;
		if(w[m]<=sqr)id1[w[m]]=m;else id2[j]=m;
	}
	for(int j=1;j<=tot;j++){
		for(int i=1;i<=m&&1ll*prime[j]*prime[j]<=w[i];i++){
			int k=(w[i]/prime[j]<=sqr)?id1[w[i]/prime[j]]:id2[n/(w[i]/prime[j])];
			(g[i]-=prime[j]*(g[k]-sp[j-1]))%=mod;
			(h[i]-=h[k]-j+1)%=mod;
		}
	}
	int ans=S(n,1)+1;
	cout<<(ans+mod)%mod;
	return 0;
}