## bzoj 4815 [Cqoi2017]小Q的表格——反演+分块

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define ll long long
using namespace std;
const int N=4e6+5,M=2005,mod=1e9+7;
int n,g[N],phi[N],f[N],pri[N];bool vis[N];
int base,bh[N],s[M],si[M],fl[N];
void upd(int &x){while(x>=mod)x-=mod;while(x<0)x+=mod;}
int gcd(int a,int b){return b?gcd(b,a%b):a;}
int pw(int x,int k)
{int ret=1;while(k){if(k&1)ret=(ll)ret*x%mod;x=(ll)x*x%mod;k>>=1;}return ret;}
void init()
{
phi[1]=g[1]=1; int cnt=0;
for(int i=2;i<=n;i++)
{
if(!vis[i])pri[++cnt]=i,phi[i]=i-1;
for(int j=1;j<=cnt&&(ll)i*pri[j]<=n;j++)
{
vis[i*pri[j]]=1;
if(i%pri[j]==0){phi[i*pri[j]]=(ll)phi[i]*pri[j]%mod;break;}
else phi[i*pri[j]]=(ll)phi[i]*phi[pri[j]]%mod;
}
g[i]=(g[i-1]+(ll)i*i%mod*phi[i])%mod;//presum!!
}
base=sqrt(n);
for(int i=1;i<=n;i++)f[i]=(ll)i*i%mod,fl[i]=(fl[i-1]+f[i])%mod;
for(int i=1,j=1,k=1;i<=n;i++,k++)
{
bh[i]=j;if(k==base)k=0,j++;
}
/*
for(int i=1,j=1,k=base;i<=n;i++)
{
f[i]=(ll)i*i%mod;bh[i]=j;
si[j]+=f[i]; upd(si[j]); fl[i]=si[j];
if(i==k)s[j]=s[j-1]+si[j],j++,k+=base;
}
*/
}
int calc(int x){int ret=fl[x]+s[bh[x]];upd(ret);return ret;}
int main()
{
int T;scanf("%d%d",&T,&n);init();
int x,y,tn; ll w;
while(T--)
{
scanf("%d%d%lld%d",&x,&y,&w,&tn);//w not %mod!!!
int u=gcd(x,y),d=bh[u];
int tf=w/(x/u)/(y/u)%mod;//not inv

/*//also ok
int chg=tf+calc(u-1)-calc(u);upd(chg);
for(int i=u;bh[i]==bh[u];i++)
fl[i]+=chg,upd(fl[i]);
for(int i=bh[u]+1;i<=bh[n];i++)//n not tn!!!
s[i]+=chg,upd(s[i]);
*/
int pl=tf-f[u];upd(pl);f[u]=tf;
for(int i=u;bh[i]==bh[u];i++)
fl[i]+=pl,upd(fl[i]);
for(int i=bh[u]+1;i<=bh[n];i++)
s[i]+=pl,upd(s[i]);
/*
si[d]=si[d]-f[u]+tf;upd(si[d]);
for(int i=d;i<=bh[n];i++)
s[i]=s[i-1]+si[i],upd(s[i]);
f[u]=tf;

fl[u]=(bh[u-1]==bh[u]?fl[u-1]:0)+f[u];
for(int i=u+1,j=d*base;i<=j;i++)
fl[i]=fl[i-1]+f[i],upd(fl[i]);
*/
int ans=0;
for(int i=1,j;i<=tn;i=j+1)
{
int d=tn/i,sm=0; j=tn/d;
/*
if(bh[j]-bh[i]<=1)
for(int l=i;l<=j;l++)
sm+=f[l],upd(sm);
else
{
sm=s[bh[j]-1]-s[bh[i]-1]+fl[j];
upd(sm);
if(bh[i-1]==bh[i])sm-=fl[i-1],upd(sm);
}
ans=(ans+(ll)sm*g[d])%mod;
*/
ans=(ans+(ll)(calc(j)-calc(i-1))*g[d])%mod;upd(ans);
}
printf("%d\n",ans);
}
return 0;
}

posted on 2018-12-14 09:03  Narh  阅读(193)  评论(0编辑  收藏  举报