「CQOI2017」小Q的表格 題解
八萬年沒寫過莫反了……隨便寫點。
先簡單的列出性質：
- \(f(a,b)=f(b,a)\)；
- \(b f(a,a+b) = (a+b)f(a,b)\)。
先考慮這個 \(f\) 兩個恆等式的形式。將 \(f\) 視作一個狀態，那麼有轉移 \((a,a+b) \to (a,b)\)，且 \((a,b) \to (b,a)\)。
不難聯想到輾轉相減求最大公約數的過程，最終 \((a,b)\) 一定可以轉移到 \((\gcd(a,b),\gcd(a,b))\) 上。
那麼，表格上的所有數都可以用對應的一個對角線上的數乘上一個係數表示。修改一個數的話，我們只用管它的對角線元素就行了。
注意到 \(b f(a,a+b) = (a+b)f(a,b)\)
下面的 \(i/j\) 均表示 \(\left\lfloor \dfrac{i}{j} \right\rfloor\)。
然後求和的工作話，我們先記 \(g(i) = f(i,i)\)，直接推式子：
\[\begin{aligned} \sum_{i=1}^k \sum_{j=1}^k f(i,j) &= \sum_{d=1}^n g(d) \sum_{i=1}^n \dfrac{i}{d} \sum_{j=1}^n \dfrac{j}{d} [\gcd(i,j)=d] \\ &= \sum_{d=1}^n g(d) \sum_{i=1}^{n/d} i \sum_{j=1}^{n/d} j[\gcd(i,j)=1] \end{aligned} \]注意到後面的 \(\displaystyle \sum_{j=1}^{\lfloor \frac{n}{d} \rfloor} j[\gcd(i,j)=1]\)
首先這個上界不優美。注意到在原式中，\(i=j \neq 1\) 的時候不會造成貢獻，並且 \(i \neq j\) 的貢獻會算兩次。我們強制 \(i>j\)，然後乘 \(2\) 即可。（下面用 \(\frac{i}{D}\) 的原因是其一定是一個整數）
\[\begin{aligned} 2\sum_{j=1}^i j[\gcd(i,j)=1] &= 2\sum_{j=1}^i j \sum_{D|\gcd(i,j)} \mu(D) \\ &= 2\sum_{D|i} \mu(D) \sum_{j=1}^{\frac{i}{D}} jD \\ &= 2\sum_{D|i} \mu(D) D \sum_{j=1}^{\frac{i}{D}} j \\ &= 2\sum_{D|i} \mu(D) D \dfrac{\left(\frac{i}{D}+1\right)\times \frac{i}{D}}{2} \\ &= \sum_{D|i} \mu(D) D \left(\frac{i}{D}+1\right)\times \frac{i}{D} \\ &= i\sum_{D|i} \mu(D) \left(\frac{i}{D}+1\right) \end{aligned} \]推到這裡基本可以結束了。注意到和式是 \(\mu * (id+1)\)
那麼繼續推之前的式子：
\[\begin{aligned} \sum_{i=1}^k \sum_{j=1}^k f(i,j) &= \sum_{d=1}^n g(d) \sum_{i=1}^n \dfrac{i}{d} \sum_{j=1}^n \dfrac{j}{d} [\gcd(i,j)=d] \\ &= \sum_{d=1}^n g(d) \sum_{i=1}^{n/d} i^2 \varphi(i) \end{aligned} \]後面的和式是個字首和形式。然後預處理 \(i^2 \varphi(i)\) 的字首和非常簡單，外層直接整除分塊做就可以了。
但是還有個修改。這是個方程形式，採用之前的那個關於 \(f(a,b),g(\gcd(a,b))\) 的等式，採用樹狀陣列求出新的 \(g(\gcd(a,b))\)，直接改就行了。求 \(g\) 的區間和也可以用樹狀陣列。
但是有點慢……還湊合，就這樣吧。時間複雜度 \(O(n \sqrt n \log n)\)。不懂這個題的神性在哪裡。
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
char buf[1<<18],*p1=buf,*p2=buf;
#define getchar() (p1==p2 && (p2=(p1=buf)+fread(buf,1,1<<18,stdin),p1==p2)?EOF:*p1++)
LL read()
{
LL x=0;
char c=getchar();
while(c<'0' || c>'9') c=getchar();
while(c>='0' && c<='9') x=(x<<1)+(x<<3)+(c^'0'),c=getchar();
return x;
}
void write(LL x)
{
if(x>9) write(x/10);
putchar(x%10+'0');
}
const LL MOD=1e9+7;
inline LL Add(LL x,LL y){return x+y>=MOD?x+y-MOD:x+y;}
inline LL Sub(LL x,LL y){return x<y?x-y+MOD:x-y;}
inline LL Mul(LL x,LL y){return 1ll*x*y%MOD;}
LL QuickPow(LL x,LL p)
{
LL ans=1,base=x;
while(p)
{
if(p&1) ans=Mul(ans,base);
base=Mul(base,base);
p>>=1;
}
return ans;
}
LL N,m;
inline LL lowbit(LL x){return x&(-x);}
struct BinaryIndexedTree{
LL tr[4000005];
void modify(LL x,LL val){for(LL i=x;i<=N;i+=lowbit(i)) tr[i]=Add(tr[i],val);}
LL query(LL x){LL ans=0;for(LL i=x;i;i^=lowbit(i)) ans=Add(ans,tr[i]);return ans;}
LL query(LL l,LL r){return Sub(query(r),query(l-1));}
}bit;
LL cnt,prime[4000005];
bool vis[4000005];
LL phi[4000005],f[4000005];
void shai(LL up)
{
vis[phi[1]=1]=true;
for(LL i=2;i<=up;++i)
{
if(!vis[i]) phi[i]=i-1,prime[++cnt]=i;
for(LL j=1;j<=cnt && i*prime[j]<=up;++j)
{
vis[i*prime[j]]=true;
if(i%prime[j]==0)
{
phi[i*prime[j]]=Mul(phi[i],prime[j]);
break;
}
phi[i*prime[j]]=Mul(phi[i],prime[j]-1);
}
}
for(LL i=1;i<=up;++i) f[i]=Mul(Mul(i,i),phi[i]),f[i]=Add(f[i],f[i-1]);
}
LL gcd(LL a,LL b){return !b?a:gcd(b,a%b);}
LL g[4000005];
int main(){
shai(4000000);
m=read(),N=read();
for(LL i=1;i<=N;++i) bit.modify(i,g[i]=Mul(i,i));
while(m-->0)
{
LL a=read(),b=read(),x=read(),n=read();
LL d=gcd(a,b);
bit.modify(d,MOD-g[d]);
long long tmp=x;
tmp/=a/d;
tmp/=b/d;
tmp%=MOD;
g[d]=tmp;
bit.modify(d,g[d]);
LL ans=0;
for(LL l=1,r;l<=n;l=r+1) r=n/(n/l),ans=Add(ans,Mul(f[n/l],bit.query(l,r)));
write(ans),puts("");
}
return 0;
}