洛咕 P3700 [CQOI2017]小Q的表格

神仙題orz

首先推一下給的兩個式子中的第二個

\(b\cdot F(a,a+b)=(a+b)\cdot F(a,b)\)

先簡單的想，\(F(a,a+b)\)和\(F(a,b)\)會相互影響

可以換一種角度想，\(F(a,b-a)\)和\(F(a,b)\)會相互影響\((b>a)\)

那麼可以從\(F(x,y)\)一路推下去

\(F(x,y)=F(x,y-x)=F(x,y-2x)=\cdots=F(x,y\mod x)\)

（注意這裡的\(\text{mod}\)結果是0的話就沒有辦法再減了）

這時橫座標比縱座標大了，利用題目給的式子1，swap橫縱座標

\(F(x,y)=F(x,y\mod x)=F(y\mod x,x)=F(y\mod x,x\mod(y\mod x))=\cdots\)

總結一下，如果繼續這樣推下去，當橫、縱座標相等時會就不能再減了

剛才是怎麼推的呢，就是當x,y不等的時候每次都把x對y取模然後交換x,y

是不是很熟悉，就是求gcd的過程

那麼可以推出來，\(\gcd(x,y)\)相等的會相互影響

具體影響多少是顯然的，\(F(x,y)=F(\gcd(x,y),\gcd(x,y))\times \frac{xy}{\gcd(x,y)^2}\)

所以只要維護\(F(i,i)\)的值就行了，下面設\(F(i,i)=F(i)\)

現在要算\(\sum_{i=1}^n\sum_{j=1}^nF(i,j)\)

\(ans=\sum_{i=1}^n\sum_{j=1}^nF(\gcd(i,j))\)

考慮列舉\(\gcd(i,j)\)，

\(ans=\sum_{k=1}^nF(k)\sum_{i=1}^n\sum_{j=1}^n\frac{ij}{k^2}[\gcd(i,j)=k]\)

\(ans=\sum_{k=1}^nF(k)\sum_{i=1}^{n/k}\sum_{j=1}^{n/k}ij[\gcd(i,j)=1]\)

設\(g(n)=\sum_{i=1}^{n}\sum_{j=1}^{n}ij[\gcd(i,j)=1]\)

可以只求一半，\(g(n)=2\sum_{i=1}^{n}\sum_{j=1}^{i}ij[\gcd(i,j)=1]-\sum_{i=1}^{n}i^2[\gcd(i,i)=1]\)

右邊那一塊是顯然的

\(g(n)=2\sum_{i=1}^{n}\sum_{j=1}^{i}ij[\gcd(i,j)=1]-1\)

設\(h(n)=\sum_{i=1}^{n}in[gcd(i,n)=1]\)，\(g(n)=2\sum_{i=1}^{n}h(i)-1\)

怎麼求\(h\)呢，顯然可行的\(i\)有\(\varphi(n)\)種，如果\(\gcd(i,n)=1\)，那麼\(\gcd(n-i,n)=1\)

所以\(i\)和\(n-i\)可以配對，加起來是\(n\)，一共\(\varphi(n)/2\)對，\(h(n)=\frac{\varphi(n)\times n^2}{2}\)，注意特判\(h(n)=1\)，因為不能配對

然後求出了\(h\)就可以求出\(g\)，\(F\)每次只會修改一個。要求的\(ans=\sum_{k=1}^nF(k)g(n/k)\)顯然\(n/k\)只有根號種取值，數論分塊即可，\(F\)樹狀陣列字首和

#include<bits/stdc++.h>
#define il inline
#define vd void
#define int ll
#define mod 1000000007
typedef long long ll;
il int gi(){
    int x=0,f=1;
    char ch=getchar();
    while(!isdigit(ch)){
        if(ch=='-')f=-1;
        ch=getchar();
    }
    while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
    return x*f;
}
int g[4000010],h[4000010];
int pri[4000010],pr,phi[4000010];
bool yes[4000010];
int n,m,f[4000010];
il vd update(int x,int p){while(x<=n)f[x]=(f[x]+p)%mod,x+=x&-x;}
il int query(int x){int ret=0;while(x)ret=(ret+f[x])%mod,x-=x&-x;return ret;}
signed main(){
    m=gi(),n=gi();
    int x,y,k,t,G;
    for(int i=1;i<=n;++i)update(i,1ll*i*i%mod);
    phi[1]=1;
    for(int i=2;i<=n;++i){
        if(!yes[i])phi[i]=i-1,pri[++pr]=i;
        for(int j=1;j<=pr&&i*pri[j]<=n;++j){
            yes[i*pri[j]]=1;
            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)h[i]=1ll*phi[i]*i%mod*i%mod;
    g[1]=1;
    for(int i=2;i<=n;++i)g[i]=(g[i-1]+h[i])%mod;
    while(m--){
        x=gi(),y=gi();G=std::__gcd(x,y);
        t=(gi()/(x/G)/(y/G))%mod;
        update(G,(0ll+t-query(G)+query(G-1)+mod)%mod);
        k=gi();
        int ans=0;
        for(int i=1;i<=k;++i){
            int j=k/(k/i);
            ans=(ans+1ll*(query(j)-query(i-1)+mod)*g[k/i])%mod;
            i=j;
        }
        printf("%lld\n",ans);
    }
    return 0;
}