莫比乌斯反演,杜教筛

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\n",ans);
  }

 ------------------------------------------------------------------------------

Codechef DEC16 BOUNCE

求$\sum_{i=1}^n\sum_{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_g\sum_{i=1}^{\lfloor \frac{n}{g} \rfloor}\sum_{j=1}^{\lfloor \frac{n}{g} \rfloor}[gcd(i,j)=g]\cdot i+\sum_g\sum_{i=1}^{\lfloor \frac{n}{g} \rfloor}\sum_{j=1}^{\lfloor \frac{n}{g} \rfloor}[gcd(i,j)=g]\cdot j+\sum_{i=1}^n\sum_{j=1}^n-2$

三项分别考虑,以第一项为例反演。

$\sum_g\sum_d\mu(d) \cdot d\sum_{i=1}^{\lfloor \frac{n}{d*g}\rfloor}\sum_{j=1}^{\lfloor \frac{n}{d*g}\rfloor}i$

上面式子中的前半部分就是$1\ast\mu\cdot n$

证明:$\varphi=\mu \ast n ,1\ast\mu \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}^n\varphi(i)\sum_{j=1}^{\lfloor \frac{n}{i} \rfloor}f(j)=1$

将i=1的项提出可得

$\sum_{i=1}^nf(i)=1-\sum_{i=2}^n\varphi(i)\sum_{j=1}^{\lfloor \frac{n}{i} \rfloor}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\n");
          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\n");else{
        ans=0;
        for (int i=0;i<=10000;i++) tab[i]=-1e9,tabinv[i]=-1e9;
        solve(m);
        printf("%lld\n",ans);    
      }
    }
  } 

 

posted @ 2017-07-29 15:31  z1j1n1  阅读(228)  评论(0编辑  收藏  举报