BZOJ4176
#include <cstdio> #include <map> #define LL long long using namespace std; map <LL,LL> mpa; map <LL,LL> mpb; const LL mo=1e9+7; int b[10000001],phi[10000001],ss[10000001]; int miu[10000001],summiu[10000001],cnt,n; void euler(int n=1e7){ b[1]=1;miu[1]=1;summiu[1]=1; for (LL i=2;i<=n;i++){ if (b[i]==0) {ss[++cnt]=i;miu[i]=-1;} summiu[i]=summiu[i-1]+miu[i]; for (LL j=1;j<=cnt;j++){ if (i*ss[j]>n) break; b[i*ss[j]]=1; if (i%ss[j]==0) {miu[i*ss[j]]=0;break;} else miu[i*ss[j]]=miu[i]*(-1); } } } LL getsummiu(LL num){ if (num<=1e7) return(summiu[num]); if (mpb[num]!=0) return(mpa[num]);//可以不使用map，可能出现的值为n/i的根号n个值 LL res=1;//筛phi时为n*(n+1)/2 for (LL i=2,last;i<=num;i=last+1){ last=min((num/(num/i)),num); res-=(last-i+1)*getsummiu(num/i)%mo; res%=mo; } mpb[num]=1;mpa[num]=res; return(res); } LL f(LL n){ LL ret=0; for (int i=1,last;i<=n;i=last+1){ last=min(n/(n/i),n); ret+=(last-i+1)*(n/i)%mo; ret%=mo; } return(ret*ret%mo); } int main(){ scanf("%d",&n); euler(); LL ans=0; for (int d=1,last;d<=n;d=last+1){ last=min(n/(n/d),n); ans+=(getsummiu(last)-getsummiu(d-1))*f(n/d)%mo; ans%=mo; } ans%=mo;ans+=mo;ans%=mo; printf("%lld ",ans); }
------------------------------------------------------------------------------
Codechef DEC16 BOUNCE
求$sum_{i=1}^nsum_{j=1}^n[j>frac{lob1}{lob2}]cdot[ j<frac{upb1}{upb2}]cdot(frac{lcm(i,j)}{i}+frac{lcm(i,j)}{j}-2)$
先不考虑两条直线带来的贡献，考虑化简
$sum_gsum_{i=1}^{lfloor frac{n}{g} floor}sum_{j=1}^{lfloor frac{n}{g} floor}[gcd(i,j)=g]cdot i+sum_gsum_{i=1}^{lfloor frac{n}{g} floor}sum_{j=1}^{lfloor frac{n}{g} floor}[gcd(i,j)=g]cdot j+sum_{i=1}^nsum_{j=1}^n-2$
三项分别考虑，以第一项为例反演。
$sum_gsum_dmu(d) cdot dsum_{i=1}^{lfloor frac{n}{d*g} floor}sum_{j=1}^{lfloor frac{n}{d*g} floor}i$
上面式子中的前半部分就是$1astmucdot n$
证明:$varphi=mu ast n ,1astmu cdot n ast mu ast n=mu cdot n ast n ast 1 ast mu=varepsilon$
我们称该函数为f，则$varphi ast f=varepsilon$
考虑杜教筛，要求一个函数的前缀和，可以考虑构造另一个函数，使得我们能求出另一个函数的前缀和与两个函数狄利克雷卷积的前缀和。
这里我们找到了$varphi$
$sum_{i=1}^nvarphi(i)sum_{j=1}^{lfloor frac{n}{i} floor}f(j)=1$
将i=1的项提出可得
$sum_{i=1}^nf(i)=1-sum_{i=2}^nvarphi(i)sum_{j=1}^{lfloor frac{n}{i} floor}f(j)$
递归求解即可。
对于给出的限制可以用类欧几里得算法求解。
#include <bits/stdc++.h> #define LL long long using namespace std; const int lim=1.5e7; const LL mo=1e9+7; const LL inv2=(mo+1)/2; const LL inv6=(mo+1)/6; struct data{ LL f,g,h; }; LL upb1,upb2,lob1,lob2,ans; int T,n,a[2000001],bt[lim+1],miu[lim+1],ss[lim+1],cnt; LL tab[200001],sumphi[lim+1],tabinv[200001],suminvphi[lim+1],phi[lim+1],m; LL conv[lim+1]; char st[2000001]; void updupb(LL x,LL y){ if (upb1*y>x*upb2) upb1=x,upb2=y; } void updlob(LL x,LL y){ if (lob1*y<x*lob2) lob1=x,lob2=y; } data likegcd(LL a,LL b,LL c,LL n){ data ret;ret.f=ret.g=ret.h=0; LL N=n%mo,N_N1=N*(N+1)%mo,N_N1_2N1=N%mo*(N+1)%mo*(2*N+1)%mo; if (a==0) return(ret); if (a>=c||b>=c){ LL N_N1_inv2=N_N1*inv2%mo,N_N1_2N1_inv6=N_N1_2N1*inv6%mo; data t=likegcd(a%c,b%c,c,n); ret.f+=t.f+a/c*N_N1_inv2%mo+b/c*(N+1)%mo;ret.f%=mo; ret.g+=t.g+a/c*N_N1_2N1_inv6%mo+(b/c)*N_N1_inv2%mo;ret.g%=mo; ret.h+=(a/c)*(a/c)%mo*N_N1_2N1_inv6%mo+ (b/c)*(b/c)%mo*(N+1)%mo+ (a/c)*(b/c)%mo*N_N1%mo+ t.h+ 2*(a/c)%mo*t.g%mo+ 2*(b/c)%mo*t.f%mo; ret.h%=mo; return(ret); }else{ LL m=(a*n+b)/c,M=m%mo; data t=likegcd(c,c-b-1,a,m-1); ret.f=N*M%mo-t.f;ret.f%=mo; ret.g=(N_N1%mo*M%mo-t.f-t.h)%mo*inv2%mo;ret.g%=mo; ret.h=N*M%mo*(M+1)%mo-2*t.g-2*t.f-ret.f;ret.h%=mo; return(ret); } } LL calc(LL n){ LL N=n%mo; if (lob1!=0){ LL sumi=0,sumj=0; data t=likegcd(upb1,0,upb2,n); sumi+=t.g; sumi-=(upb2+n/upb2*upb2)%mo*((n/upb2)%mo)%mo*inv2%mo; t=likegcd(lob1,0,lob2,n); sumi-=t.g; t=likegcd(lob2,0,lob1,n*lob1/lob2); sumj+=t.g; LL tupb=(n*lob1/lob2)/lob1; sumj-=(lob1+tupb%mo*lob1)%mo*((tupb)%mo)%mo*inv2%mo; t=likegcd(upb2,0,upb1,n*upb1/upb2); sumj-=t.g; LL ran=n*upb1/upb2-n*lob1/lob2,rr=n*upb1/upb2,rl=n*lob1/lob2+1; sumj+=N*((rl+rr)%mo)%mo*(ran%mo)%mo*inv2%mo; return((sumi+sumj)%mo); }else{ LL sumi=0,sumj=0; data t=likegcd(upb1,0,upb2,n); sumi+=t.g; sumi-=(upb2+n/upb2*upb2)%mo*((n/upb2)%mo)%mo*inv2%mo; t=likegcd(upb2,0,upb1,n*upb1/upb2); sumj-=t.g; LL ran=n*upb1/upb2,rr=n*upb1/upb2,rl=1; sumj+=N*((rl+rr)%mo)%mo*(ran%mo)%mo*inv2%mo; return((sumi+sumj)%mo); } } LL calcnod(LL n){ LL N=n%mo; if (lob1!=0){ LL nod=0; data t=likegcd(upb1,0,upb2,n); nod+=t.f;nod%=mo; nod-=n/upb2; t=likegcd(lob1,0,lob2,n); nod-=t.f;nod%=mo; return(nod); }else{ LL nod=0; data t=likegcd(upb1,0,upb2,n); nod+=t.f;nod%=mo; nod-=n/upb2;nod%=mo; return(nod); } } LL getphi(LL n,LL po){ if (po<=lim) return(sumphi[po]); if (tab[n/po]!=-1e9) return(tab[n/po]); LL ret=po%mo*((po+1)%mo)%mo*inv2%mo; for (LL i=2,nxt;i<=po;i=nxt+1){ nxt=po/(po/i); ret-=(nxt-i+1)%mo*getphi(n,po/i)%mo; ret%=mo; } tab[n/po]=ret; return(ret); } LL getinvphi(LL n,LL po){ if (po<=lim) return(suminvphi[po]); if (tabinv[n/po]!=-1e9) return(tabinv[n/po]); LL ret=1; for (LL i=2,nxt;i<=po;i=nxt+1){ nxt=po/(po/i); LL sm=getphi(n,nxt)-getphi(n,i-1); ret-=sm%mo*getinvphi(n,po/i)%mo; ret%=mo; } tabinv[n/po]=ret; return(ret); } void solve(LL n){ for (LL i=1,nxt;i<=n;i=nxt+1){ nxt=n/(n/i); LL sm=getinvphi(n,nxt)-getinvphi(n,i-1); ans+=sm*calc(n/i)%mo; ans%=mo; } LL t=calcnod(n); ans-=2*t%mo; ans%=mo;ans+=mo;ans%=mo; } void euler(){ bt[1]=1;miu[1]=1;phi[1]=1;conv[1]=1; for (int i=1;i<=lim;i++){ if (!bt[i]) {ss[++cnt]=i;miu[i]=-1;phi[i]=i-1;conv[i]=1-i;} for (int j=1;j<=cnt&&ss[j]*i<=lim;j++){ bt[ss[j]*i]=1; if (i%ss[j]==0){ miu[i*ss[j]]=0; phi[i*ss[j]]=phi[i]*ss[j]; conv[i*ss[j]]=conv[i]; break; }else miu[i*ss[j]]=(-1)*miu[i],phi[i*ss[j]]=(ss[j]-1)*phi[i],conv[i*ss[j]]=conv[i]*conv[ss[j]]; } } } void init(){ euler(); for (int i=1;i<=lim;i++) sumphi[i]=sumphi[i-1]+phi[i],sumphi[i]%=mo, suminvphi[i]=suminvphi[i-1]+conv[i],suminvphi[i]%=mo; } int gcd(int x,int y){ if (x%y==0) return(y);else return(gcd(y,x%y)); } int main(){ init(); scanf("%d",&T); while (T--){ scanf("%lld%d",&m,&n); scanf("%s",&st); int lstx=0,lsty=0,flag=1; for (int i=1;i<=n;i++) if (st[i-1]=='R'||st[i-1]=='L'){ a[i]=0; if (st[i-1]=='R'&&lstx) flag=0; if (st[i-1]=='L'&&!lstx) flag=0; lstx^=1; }else{ a[i]=1; if (st[i-1]=='T'&&lsty) flag=0; if (st[i-1]=='D'&&!lsty) flag=0; lsty^=1; } if (!flag){ printf("0 "); continue; } if (a[1]==0) for (int i=1;i<=n;i++) a[i]^=1; int cnt0=0,cnt1=0; upb1=1e9,upb2=1,lob1=0,lob2=1; for (int i=1;i<=n;i++) if (a[i]){ updupb(cnt0+1,cnt1+1); cnt1++; }else{ updlob(cnt0+1,cnt1+1); cnt0++; } int gc=gcd(upb1,upb2);upb1/=gc;upb2/=gc; gc=gcd(lob1,lob2);lob1/=gc;lob2/=gc; if (upb1*lob2<=lob1*upb2) printf("0 ");else{ ans=0; for (int i=0;i<=10000;i++) tab[i]=-1e9,tabinv[i]=-1e9; solve(m); printf("%lld ",ans); } } }