题解 P3700 【[CQOI2017]小Q的表格】
zhoutb2333
2018-04-12 09:49:12
首先题目里说$f(a,b)=f(b,a)$，说明我们只用维护半个棋盘的信息就够了。
然后$b *f(a,a+b)=(a+b)*f(a,b)$这个式子我们化成$\frac{f(a,b)}{a*b}=\frac{f(a,a+b)}{a*(a+b)}$
那么就有$\frac{f(a,b)}{a*b}=\frac{f(a,a \% b)}{a*(a \% b)}$
发现这个形式很像辗转相除法求$\gcd$的过程，所以我们如果一直这样递归下去的话可以得到$\frac{f(a,b)}{a*b}=\frac{f(\gcd(a,b),\gcd(a,b))}{\gcd^2(a,b)}$。
即$f(a,b)=\frac{a*b}{d^2}*f(d,d)$，令$d=\gcd(a,b)$
发现如果修改值的话，$\gcd$不同的位置的$f$值相互独立。也就是说修改一个位置$(a,b)$的值后有且仅有$\gcd(x,y)=d$的$(x,y)$位置的$f$值会跟着变化。（这句话对求答案没什么用，是为方便理解写的）
然后回去看$ans$，
$$\sum \limits ^{k}_{i=1}\sum \limits ^{k}_{j=1} f(i,j)$$
$$= \sum \limits ^{k}_{d=1} f(d,d) \sum \limits _{d|i} \sum \limits _{d|j} \frac{i*j}{d^2} \varepsilon(\gcd(i,j)==d)$$
$$= \sum \limits ^{k}_{d=1} f(d,d) \sum \limits^{\lfloor \frac{k}{d} \rfloor} _{i=1} \sum \limits^{\lfloor \frac{k}{d} \rfloor} _{j=1} i*j* \varepsilon(\gcd(i,j))$$
所以说我们可以根据$k/d$分块，这样需要动态维护$f(i,i)$的前缀和。观察到$n$为
$4*10^6$而$m$为$10^4$，所以我们大概需要一个数据结构，支持单点修改$O(\sqrt{n})$，查前缀和$O(1)$。分块正好可以做到这个。具体就是块内每个位置维护本块内的前缀和，然后每个块末尾维护一下总的前缀和，这样查询的时候等于上一块的总答案加上这个位置的块答案。修改暴力重构，只会影响到$\sqrt{n}$左右个值。
设$G(x)=\sum \limits^{x} _{i=1} \sum \limits^{x} _{j=1} i*j* \varepsilon(\gcd(i,j))$。那么根据$\sum \limits ^{x}_{i=1}i*\varepsilon(\gcd(i,x))=\frac{x*(\varphi(x)+\varepsilon(x))}{2}$，有$G(x)=\sum \limits ^{x}_{i=1} i^2 *\varphi(i)$（考虑$G(x)$比$G(x-1)$多了什么），这个东西可以线性筛预处理。（然而蒟蒻这一步只会$n \sqrt{n}$的算法导致做不出来了..QAQ）
注意不要爆$int$
代码略丑，维护块的方式比较清奇
``` cpp
%:pragma GCC optimize(2)
#include<bits/stdc++.h>
#define maxn 4000010
#define ll long long
using namespace std;
int pri[maxn/8],phi[maxn],pricnt=0;
int bel[maxn],bpos[2010],epos[2010],size,tn,n,q;
//bel为每个位置属于哪个块，bpos,epos为每个块的初始位置和末尾位置
ll val[maxn],ans[maxn],ANS[2010],f[maxn];
//ans为块内前缀和，ANS为总前缀和
const ll p=1e9+7;
bool notpri[maxn];
inline int gcd(int x,int y){
if(!x||!y)
return x|y;
return x>y?gcd(y,x%y):gcd(x,y%x);
}
inline ll sum(int x){
if(!x)
return 0;
return (ANS[bel[x]-1]+ans[x])%p;
}
inline ll pw(ll a,ll b){
ll ret=1;
while(b){
if(b&1)
ret=ret*a%p;
a=a*a%p;
b>>=1;
}
return ret;
}
void init(){
phi[1]=1;
for(int i=2;i<=n;i++){
if(!notpri[i]){
pri[++pricnt]=i;
phi[i]=i-1;
}
for(int j=1;j<=pricnt;j++){
if(i*pri[j]>n)
break;
notpri[i*pri[j]]=true;
if(i%pri[j]==0){
phi[i*pri[j]]=phi[i]*pri[j];
break;
}
phi[i*pri[j]]=phi[i]*(pri[j]-1);
}
}
for(int i=1;i<=n;i++)
f[i]=(f[i-1]+1LL*phi[i]*i%p*i%p)%p;
}
int main(){
scanf("%d%d",&q,&n),size=sqrt(n),init();
for(int i=1;i<=n;i++){
bel[i]=i/size+1;
if(bel[i]>bel[i-1])
bpos[bel[i]]=i,epos[bel[i]-1]=i-1;
}
tn=bel[n],epos[tn]=n;
for(int i=1;i<=n;i++)
val[i]=1LL*i*i%p;
for(int i=1;i<=tn;i++){
for(int j=bpos[i];j<=epos[i];j++)
ans[j]=(ans[j-1]*(j>bpos[i])+val[j])%p;
ANS[i]=(ANS[i-1]+ans[epos[i]])%p;
}
for(int i=1;i<=q;i++){
int a,b,k,d,pos;
ll x,Ans=0;
scanf("%d%d%lld%d",&a,&b,&x,&k),x%=p,d=gcd(a,b);
val[d]=x*d%p*d%p*pw(1LL*a*b%p,p-2)%p;
for(int j=d;j<=epos[bel[d]];j++)
ans[j]=(ans[j-1]*(j>bpos[bel[d]])+val[j])%p;
for(int j=bel[d];j<=tn;j++)
ANS[j]=(ANS[j-1]+ans[epos[j]])%p;
for(int j=1;j<=k;j=pos+1){
pos=k/(k/j);
Ans=(Ans+f[k/j]*(sum(pos)-sum(j-1)+p)%p)%p;
}
printf("%lld\n",Ans);
}
return 0;
}
```