1. 程式人生 > >BZOJ4916 神犇和蒟蒻(尤拉函式+杜教篩)

BZOJ4916 神犇和蒟蒻(尤拉函式+杜教篩)

  第一問是來搞笑的。由尤拉函式的計算公式容易發現φ(i2)=iφ(i)。那麼可以發現φ(n2)*id(n)=Σd*φ(d)*(n/d)=nΣφ(d)=n2 。這樣就有了杜教篩所要求的容易算字首和的兩個函式。一通套路即可。

#include<iostream> 
#include<cstdio>
#include<cmath>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<map>
using namespace std;
#define ll long long #define P 1000000007 #define N 1000010 char getc(){char c=getchar();while ((c<'A'||c>'Z')&&(c<'a'||c>'z')&&(c<'0'||c>'9')) c=getchar();return c;} int gcd(int n,int m){return m==0?n:gcd(m,n%m);} int read() { int x=0,f=1;char c=getchar(); while (c<'0'||c>'
9') {if (c=='-') f=-1;c=getchar();} while (c>='0'&&c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar(); return x*f; } int n,phi[N],iphi[N],prime[N],cnt,inv6=166666668; map<int,int> f; bool flag[N]; inline void inc(int &x,int y){x+=y;if (x>=P) x-=P;} int sumone(int
x){return (1ll*x*(x+1)>>1)%P;} int sumtwo(int x){return 1ll*x*(x+1)%P*(x<<1|1)%P*inv6%P;} int work(int x) { if (x<=min(n,N-10)) return iphi[x]; if (f.find(x)!=f.end()) return f[x]; int s=sumtwo(x); for (int i=2;i<=x;i++) { int t=x/(x/i); inc(s,P-1ll*(sumone(t)-sumone(i-1)+P)*work(x/i)%P); i=t; } f[x]=s;return s; } int main() { #ifndef ONLINE_JUDGE freopen("bzoj4916.in","r",stdin); freopen("bzoj4916.out","w",stdout); const char LL[]="%I64d\n"; #else const char LL[]="%lld\n"; #endif n=read();cout<<1<<endl; flag[1]=1;phi[1]=1; for (int i=2;i<=min(n,N-10);i++) { if (!flag[i]) prime[++cnt]=i,phi[i]=i-1; for (int j=1;j<=cnt&&prime[j]*i<=min(n,N-10);j++) { flag[prime[j]*i]=1; if (i%prime[j]==0) {phi[prime[j]*i]=phi[i]*prime[j];break;} phi[prime[j]*i]=phi[i]*(prime[j]-1); } } for (int i=1;i<=min(n,N-10);i++) iphi[i]=1ll*i*phi[i]%P; for (int i=1;i<=min(n,N-10);i++) inc(iphi[i],iphi[i-1]); cout<<work(n); return 0; }