[分块莫比乌斯反演] BZOJ 4815 [Cqoi2017]小Q的表格

最新推荐文章于 2022-01-28 15:38:33 发布

原创最新推荐文章于 2022-01-28 15:38:33 发布 · 616 阅读

1 ·

CC 4.0 BY-SA版权

莫比乌斯反演&杜教筛同时被 2 个专栏收录

38 篇文章

订阅专栏

分块

24 篇文章

订阅专栏

本文介绍了一种解决特定一维数论问题的方法，利用辗转相减法和欧几里得算法，通过预处理得到φ函数的前缀和，并采用分块技巧优化查询和更新操作的时间复杂度。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

那个神奇的关系式其实是辗转相减的形式
稍微发现下就能知道这其实是个一维的东西
$f_{a,b}={ab\over gcd^2(a,b)}*f_{gcd(a,b),gcd(a,b)}$
然后推一推就知道
$ans=\sum_{i=1}^n f_{i,i}*g(\lfloor {n\over i}\rfloor)$ 其中 $g_{n}=\sum_{i=1}^n i*i*\phi(i)$

对于 $f$ 的询问和修改是不均衡的修改是 $O(m)$ 次询问是 $O(m\sqrt n)$ 次
那么分快维护 $O(1)$ 查询 $O(\sqrt n)$ 修改
复杂度 $O(m\sqrt n)$

#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cmath>
using namespace std;
typedef long long ll;

const int P=1e9+7;

inline ll Pow(ll a,int b){
  ll ret=1;
  for (;b;b>>=1,a=a*a%P)
    if (b&1)
      ret=ret*a%P;
  return ret;
}

const int N=4000005;
const int _B=2005;

int prime[N],num;
int vst[N],phi[N];
inline void Pre(int n){
  phi[1]=1;
  for (int i=2;i<=n;i++){
    if (!vst[i]) prime[++num]=i,phi[i]=i-1;
    for (int j=1;j<=num && (ll)prime[j]*i<=n;j++){
      vst[prime[j]*i]=1;
      if (i%prime[j]==0){
    phi[i*prime[j]]=phi[i]*prime[j];
    break;
      }
      phi[i*prime[j]]=phi[i]*(prime[j]-1);
    }
  }
  for (int i=1;i<=n;i++) phi[i]=(phi[i-1]+(ll)i*i%P*phi[i]%P)%P;
}

int f[N];
int lp[_B],rp[_B],B,pos[N],tot;
ll pre[N],suf[N],sum[_B];
int a[N];

inline void Modify(int x,int y){
  int b=(x-1)/B+1;
  a[x]=y;
  pre[lp[b]]=a[lp[b]];
  for (int i=lp[b]+1;i<=rp[b];i++)
    pre[i]=(pre[i-1]+a[i])%P;
  suf[rp[b]]=a[rp[b]];
  for (int i=rp[b]-1;i>=lp[b];i--)
    suf[i]=(suf[i+1]+a[i])%P;
  for (int i=1;i<=tot;i++)
    sum[i]=(sum[i-1]+pre[rp[i]])%P;
}

inline ll Query(int l,int r){
  int lb=(l-1)/B+1,rb=(r-1)/B+1;
  if (lb==rb) return (pre[r]+suf[l]+P-pre[rp[lb]])%P;
  return (sum[rb-1]+P-sum[lb]+pre[r]+suf[l])%P;
}

inline ll Solve(int n){
  ll ret=0;
  for (int i=1,j;i<=n;i=j+1){
    j=n/(n/i);
    ret+=Query(i,j)*phi[n/i]%P;
  }
  return ret%P;
}

int n,m;

int main(){
  freopen("t.in","r",stdin);
  freopen("t.out","w",stdout);
  scanf("%d%d",&m,&n);
  B=sqrt(n); Pre(n);
  for (int i=1;i<=n;i++) pos[i]=(i-1)/B+1;
  tot=pos[n];
  for (int i=1;i<=tot;i++) lp[i]=(i-1)*B+1,rp[i]=i*B; rp[tot]=n;
  for (int i=1;i<=n;i++) a[i]=(ll)i*i%P;
  for (int b=1;b<=tot;b++){
    pre[lp[b]]=a[lp[b]];
    for (int i=lp[b]+1;i<=rp[b];i++)
      pre[i]=(pre[i-1]+a[i])%P;
    suf[rp[b]]=a[rp[b]];
    for (int i=rp[b]-1;i>=lp[b];i--)
      suf[i]=(suf[i+1]+a[i])%P;
    sum[b]=(sum[b-1]+pre[rp[b]])%P;
  }
  int x,y,d,K; ll z;
  while (m--){
    scanf("%d%d%lld%d",&x,&y,&z,&K);
    z%=P;
    d=__gcd(x,y);
    Modify(d,z*d%P*d%P*Pow((ll)x*y%P,P-2)%P);
    int ret=Solve(K);
    printf("%d\n",ret);
  }
  return 0;
}