Pre knowledge: Special topic of Mobius inversion (basic part)
Du Jiao sieve
Let \ (f(n) \) be an integral function
Find \ (S(n)=\sum_{i=1}^nf(i)\),\(n\leq 10^{10}\)
We consider prefixing and adding another integral function \ (f(n) \) on the \ (f(n) \) volume
This leads to a sub problem
Next, let's focus on the term \ (d=1 \), that is \ (g(1)\times S(n) \)
In this formula, \ (g(1) \) is easy to find, and the prefix sum of the convolution of \ (g \) and \ (f \) can also be easy to find by constructing an appropriate \ (g \), and the latter can be solved by integer division, blocking and recursion
If we preprocess the division of \ (S \) in advance, it can be proved that the complexity of calculating \ (S(n) \) is about \ (O(n^{\frac{2}{3}) \)
[BZOJ3944]--Sum
Most common
First look at \ (mu(n) \), our task is to construct an appropriate \ (g \)
We know that \ (\ mu \times I = e \) (\ (\ times \) represents Dirichlet convolution)
So we construct \ (g(n)=I(n)=1 \)
Then the prefix sum of \ (g \) and \ (\ mu \) is still \ (e \), so the sum of this block is 1
in other words
It can be solved by Du jiaosieve
Look at \ (\ varphi(n) \), we know \ (\ varphi(n)\times I=id \)
Then we also consider \ (g(n)=I(n)=1 \) on the volume
Then the value of the formula is
On the left is the sequence of equal differences, which is easy to sum
#include<bits/stdc++.h> using namespace std; typedef long long LL; const int N = 4641600; LL s[N]; LL mu[N]; int v[N],prime[N],tot=0; LL phi[N]; LL p[N]; void init(int n) { phi[1]=1; mu[1]=1; for(int i=2;i<=n;i++) { if(!v[i]) { prime[++tot]=i; v[i]=i; mu[i]=-1; phi[i]=i-1; } for(int j=1;j<=tot;j++) { if(prime[j]>v[i]||i*prime[j]>n) break; v[i*prime[j]]=prime[j]; if(i%prime[j]==0) { mu[i*prime[j]]=0; phi[i*prime[j]]=phi[i]*prime[j]; break; } else { mu[i*prime[j]]=mu[i]*mu[prime[j]]; phi[i*prime[j]]=1ll*phi[i]*(prime[j]-1); } } } for(int i=1;i<=n;i++) mu[i]=mu[i-1]+mu[i]; for(int i=1;i<=n;i++) phi[i]=phi[i-1]+phi[i]; } bool vis[N]; LL S(LL n,int k) { if(n<=N-10) return mu[n]; if(vis[k]) return s[k]; LL res=1; LL l=2,r; for(;l<=n;l=r+1) { r=n/(n/l); res=res-1ll*(r-l+1)*S(n/l,1ll*k*l); } vis[k]=1; s[k]=res; return res; } LL P(LL n,int k) { if(n<=N-10) return phi[n]; if(vis[k]) return p[k]; LL res=1ll*(1ll+n)*n/2; LL l=2,r; for(;l<=n;l=r+1) { r=n/(n/l); res=res-1ll*(r-l+1)*P(n/l,k*l); } vis[k]=1; p[k]=res; return res; } int main() { freopen("sum.in","r",stdin); freopen("sum.out","w",stdout); init(N-10); int T; cin>>T; while(T--) { LL n; cin>>n; cout<<P(n,1)<<' '; memset(vis,0,sizeof(vis)); cout<<S(n,1)<<endl; memset(vis,0,sizeof(vis)); } return 0; }
[HDOJ5608]function
The title gives a function \ (F(n)=n^2-3n+2 \)
And he guarantees that \ (F(n)=\sum_{d|n}f(d)\)
Find \ (\ sum {I = 1} ^ NF (I) \)
We first find the recurrence formula of \ (f(n) \)
Mobius can be obtained by inversion
\(f(n)=\sum_{d|n}F(d)\mu(\frac{n}{d})\)
This is obviously an integral function
However, the scope of \ (n \) is very large. Consider Du Jiao sieve
It is found that the \ (f (n) \) given by the topic can be regarded as the convolution of \ (f (n) \) and \ (I(n) \)
So the answer is
The three items in the left half can be calculated separately, of which
\(1^2+2^2+3^2+......n^2=n(n+1)(2n+1)/6\)
#include<bits/stdc++.h> using namespace std; typedef long long LL; const int N = 1e6+7; const LL mod = 1e9+7; LL Pow(LL a,LL b) { LL res=1; while(b) { if(b&1) res=1ll*res*a%mod; a=1ll*a*a%mod; b>>=1; } return res; } LL inv(LL n) { return Pow(n,mod-2); } LL inv6=Pow(3,mod-2),inv2=Pow(2,mod-2); LL F(LL n){return (1ll*n*(n+1)/2%mod*(2ll*n+1)%mod*inv6%mod-3*(1+n)*n/2%mod+2*n+mod)%mod;} int mu[N],prime[N],tot=0; int v[N]; LL f[N]; void init(int n) { mu[1]=1; f[1]=0; for(int i=2;i<=n;i++) { if(!v[i]) { v[i]=i; mu[i]=-1; prime[++tot]=i; } for(int j=1;j<=tot;j++) { if(prime[j]>v[i]||i*prime[j]>n) break; v[i*prime[j]]=prime[j]; if(i%prime[j]==0) { mu[i*prime[j]]=0; break; } else mu[i*prime[j]]=mu[i]*mu[prime[j]]; } } for(LL d=1;d<=n;d++) for(LL T=d;T<=n;T+=d) f[T]=(f[T]+(d*d-3*d+2)%mod*mu[T/d]+mod)%mod; for(int i=1;i<=n;i++) f[i]=(f[i-1]+f[i])%mod; } LL bin[N],top=0; int tag[N]; int cnt=0; unordered_map<LL,LL>s,vis; LL S(LL n) { if(n<=1e6) return f[n]; if(vis[n]) return s[n]; LL res=F(n); LL l=2,r; for(;l<=n;l=r+1) { r=n/(n/l); res=res-(1ll*(r-l+1)*S(n/l)%mod); res=(res%mod+mod)%mod; } vis[n]=1; s[n]=res; return res; } inline int read() { int X=0; bool flag=1; char ch=getchar(); while(ch<'0'||ch>'9') {if(ch=='-') flag=0; ch=getchar();} while(ch>='0'&&ch<='9') {X=(X<<1)+(X<<3)+ch-'0'; ch=getchar();} if(flag) return X; return ~(X-1); } void write(LL x) { if(x < 0) { putchar('-'); x = -x; } if(x > 9) write(x/10); putchar(x % 10 + '0'); return; } int main() { freopen("a.in","r",stdin); freopen("a.out","w",stdout); init(1e6); LL T; cin>>T; LL n; while(T--) { n=read(); if(n<=1e6) { write(f[n]); printf("\n"); continue; } write(S(n)); printf("\n"); } return 0; }
[51nod 2026]Gcd and Lcm
Title Requirements
Where \ (f(n)=\sum_{d|n}\mu(d)d\)
At first, I wanted to reverse this formula, but I found that I couldn't do it anymore
Consider analyzing the properties of a seat \ (f(n) \)
This function is obviously an integral function, so we consider each prime factor separately
\(f(p^k)=1-p \), that is to say, for each prime factor, its index has nothing to do with it, so the \ (f \) of a number is only related to the number of types of other prime factors
Then we analyze \ (gcd(i,j) \) and \ (lcm(i,j) \)
1: For the common factors of \ (i,j \), they are multiplied twice, let them be A
2: For the respective factors of \ (i,j \), they are multiplied only once, assuming that they are B and C respectively
That is, A(BCA)
We can get (AB)(AC) by exchanging the order
That is \ (f(gcd(i,j))\times f(lcm(i,j))=f(i)\times f(j) \)
The formula becomes
The problem becomes to find the prefix and sum of \ (f \)
that is
The second half can be divided into blocks
The problem is further transformed into finding the prefix and sum of \ (\ mu(d)d \)
Consider Du Jiao sieve
We consider constructing a \ (g(n)=Id(n)=n \)
That is, $$g\times (\mu(d)d)$$
So the prefix and are
Just take it in
#include<bits/stdc++.h> using namespace std; typedef long long LL; const int mod = 1e9+7; const int N =1e6+7; int mu[N]; LL F[N]; int v[N],prime[N],tot=0; void init(int n) { mu[1]=1; for(int i=2;i<=n;i++) { if(!v[i]) { v[i]=i; prime[++tot]=i; mu[i]=-1; } for(int j=1;j<=tot;j++) { if(v[i]<prime[j]||i*prime[j]>n) break; v[i*prime[j]]=prime[j]; if(i%prime[j]==0) { mu[i*prime[j]]=0; break; } else mu[i*prime[j]]=mu[i]*mu[prime[j]]; } } for(LL i=1;i<=n;i++) F[i]=(F[i-1]+mu[i]*i%mod)%mod; } unordered_map<LL,LL> s,vis; LL S(LL n) { if(n<=1e6) return F[n]; if(vis[n]) return s[n]; LL res=1; LL l=2,r; for(;l<=n;l=r+1) { r=n/(n/l); res=(res-(l+r)*(r-l+1)/2%mod*S(n/l)%mod+mod)%mod; } vis[n]=1; s[n]=res; return res; } LL f(LL n) { LL l=1,r; LL ans=0; for(;l<=n;l=r+1) { r=(n/(n/l)); ans=(ans+(S(r)-S(l-1)+mod)%mod*(n/l)%mod)%mod; } return ans; } int main() { freopen("c.in","r",stdin); freopen("c.out","w",stdout); init(1e6); LL n; cin>>n; cout<<f(n)*f(n)%mod; return 0; }
[51Nod 1238 sum of least common multiples V3
cancer
seek
For more details, see Special topic of Mobius inversion (basic part) [BZOJ2693] jzptab in
Let's use the last formula directly here
Where \ (S(n)=\sum_{i=1}^ni\)
So the key problem is to screen out \ (f (n) = n \ sum_ Prefix and of {d|n} \ mu (d) d \)
Obviously, in the front part, we can roll up a \ (g(n)=n^2 \) to teach the sieve. The whole formula can be divided into blocks
However, don't forget that there is an integer division block outside, so it becomes \ (O(n) \)
So we can only consider other methods
By observing and typing tables, we find that \ (f\times (n^2) = (n) \)
As for the specific proof, I'm sorry I don't understand, so we put a \ (g(n)=n^2 \) on the volume
Prefix and change to
It can be solved
#include<bits/stdc++.h> using namespace std; const int N = 4641600; typedef unsigned long long LL; int mu[N],prime[N],tot=0; int v[N]; LL F[N]; const LL mod = 1e9+7; void init(LL n) { F[1]=1; for(int i=2;i<=n;i++) { if(!v[i]) { v[i]=i; F[i]=1LL*i*((1-i+mod)%mod)%mod%mod; prime[++tot]=i; } for(int j=1;j<=tot;j++) { if(prime[j]>v[i]||i*prime[j]>n) break; v[i*prime[j]]=prime[j]; if(i%prime[j]==0) { F[i*prime[j]]=1LL*F[i]*prime[j]%mod; break; } else F[i*prime[j]]=1LL*F[i]*F[prime[j]]%mod; } } for(int i=1;i<=n;i++) F[i]=(F[i-1]+F[i])%mod; } unordered_map<LL,LL>s,vis; LL Pow(LL a,LL b) { LL res=1; while(b) { if(b&1) res=1LL*res*a%mod; a=1LL*a*a%mod; b>>=1LL; } return res; } LL inv6=Pow(6,mod-2),inv2=Pow(2,mod-2); LL sqr(LL n){n%=mod;return 1LL*n*(n+1)%mod*(2LL*n+1)%mod*inv6%mod;} LL sum(LL n){n%=mod;return 1LL*(1+n)%mod*n%mod*inv2%mod;} bool flag=1; LL Sum(LL n) { if(n<=N-10) return F[n]; if(vis[n]) return s[n]; LL res=sum(n); LL l=2,r; for(;l<=n;l=r+1) { r=n/(n/l); res=(res-1LL*((sqr(r)-(sqr(l-1)))+mod)%mod*Sum(n/l)%mod+mod)%mod; } vis[n]=1; s[n]=res; return res; } int main() { freopen("d.in","r",stdin); freopen("d.out","w",stdout); init(N-10); LL n; cin>>n; LL ans=0; LL l=1,r; for(;l<=n;l=r+1) { r=(n/(n/l)); ans=(ans+1LL*sum(n/l)*sum(n/l)%mod*((Sum(r)-Sum(l-1)+mod)%mod)%mod)%mod; } cout<<ans; return 0; }
51 [Nod 1237] sum of maximum common divisor V3
See the above for details Special topic of Mobius inversion (basic part) [NOI2010day1]T1 energy acquisition in
The final formula is
Just process the prefix and of \ (\ varphi(n) \) directly
#include<bits/stdc++.h> using namespace std; typedef long long LL; const LL mod = 1e9+7; const LL N = 4641600; LL phi[N]; LL v[N],prime[N],tot=0; void init(LL n) { phi[1]=1; for(LL i=2;i<=n;i++) { if(!v[i]) { v[i]=i; phi[i]=i-1; prime[++tot]=i; } for(LL j=1;j<=tot;j++) { if(prime[j]>v[i]||i*prime[j]>n) break; v[i*prime[j]]=prime[j]; if(i%prime[j]==0) { phi[i*prime[j]]=1LL*phi[i]*prime[j]%mod; break; } else phi[i*prime[j]]=1LL*phi[i]*(prime[j]-1)%mod; } } for(LL i=1;i<=n;i++) phi[i]=(phi[i-1]+phi[i])%mod; } LL Pow(LL a,LL b) { LL res=1; while(b) { if(b&1) res=1LL*res*a%mod; a=1LL*a*a%mod; b>>=1; } return res; } LL inv2=Pow(2ll,mod-2); unordered_map<LL,LL>vis,s; LL mul(LL a,LL b) { LL res=0; while(b) { if(b&1) res=(res+a)%mod; a=(a+a)%mod; b>>=1; } return res; } LL Sum(LL n) { if(n<=N-10) return phi[n]; if(vis[n]) return s[n]; LL res=mul(mul(n,n+1),inv2); LL l=2,r; for(;l<=n;l=r+1) { r=(n/(n/l)); res=(res-1LL*(r-l+1)%mod*Sum(n/l)%mod)%mod; } res=(res+mod)%mod; vis[n]=1; s[n]=res; return res; } int main() { freopen("d.in","r",stdin); freopen("d.out","w",stdout); init(N-10); LL n; cin>>n; LL res=0; LL l=1,r; for(;l<=n;l=r+1) { r=n/(n/l); res=(res+mul(Sum(r)-Sum(l-1),mul(n/l,n/l))+mod)%mod; } cout<<res; return 0; }
[bzoj 4176] Lucas's number theory
See the above for details Special topic of Mobius inversion (basic part) The sum of approximate numbers of [SDOI2015] in
\(\ mu \) can be processed by Du Jiao sieve, and the rear part can be divided into blocks
In other words, didn't it become \ (O(n) \)
#include<bits/stdc++.h> using namespace std; typedef long long LL; const int N = 2500005; const int mod = 1e9+7; int mu[N],prime[N],tot=0; LL v[N]; void init(int n) { mu[1]=1; for(int i=2;i<=n;i++) { if(!v[i]) { v[i]=i; mu[i]=-1; prime[++tot]=i; } for(int j=1;j<=tot;j++) { if(prime[j]>v[i]||i*prime[j]>n) break; v[i*prime[j]]=prime[j]; if(i%prime[j]==0) { mu[i*prime[j]]=0; break; } else mu[i*prime[j]]=mu[i]*mu[prime[j]]; } } for(int i=1;i<=n;i++) mu[i]=(mu[i]+mu[i-1]); } unordered_map<LL,LL> s; LL Sum(LL n) { if(n<=N-10) return mu[n]; if(s[n]) return s[n]; LL res=1; LL l=2,r; for(;l<=n;l=r+1) { r=n/(n/l); res=(res-1ll*(r-l+1)*Sum(n/l)%mod)%mod; } res=(res+mod)%mod; s[n]=res; return res; } LL S(LL n) { int l=1,r; LL res=0; for(;l<=n;l=r+1) { r=n/(n/l); res=(res+(LL)(n/l)*(r-l+1)%mod)%mod; } return res; } int main() { freopen("math.in","r",stdin); freopen("math.out","w",stdout); init(N-10); LL n; cin>>n; LL res=0; int l=1,r; for(;l<=n;l=r+1) { r=n/(n/l); LL val=S(n/l); res=(res+1ll*val*val%mod*(Sum(r)-Sum(l-1)+mod)%mod)%mod; } cout<<res; return 0; }