题目大意
维护一个 \(n\) 阶矩阵 \(A\)，其中最开始 \(A_{i,j}=\gcd(i,j)\)，共有 \(q\) 次操作，每次操作将矩阵某一行或某一列同乘一个数 \(y\),求每一次操作后矩阵的所有元素之和。对 \(10^9+7\) 取模。
\(n,q\leq 10^5\)，且保证数据随机生成。
思路
根据欧拉函数的性质，有
\[\sum_{c|i,c|j}\varphi(c)=\gcd(i,j) \]则考虑维护 \(n\) 个矩阵 \(M\)，\(M_i\) 的大小为 \(\lfloor\frac{n}{i}\rfloor\)，表示的是在其中一个因数为 \(i\) 时另一个要满足的因数的系数，即在上述查询中已经满足 \(c|i\)，下面要查 \(c|j\)，则一开始每个矩阵的初始元素为 \(\varphi(i)\)。不难发现所有 \(M_i\) 的所有元素的和即为答案。
之后因为行与列思路相等，不妨先考虑行。则假设修改第 \(x\) 行。则查询 \(x\) 的每一个因数 \(c_i\)。将 \(M_{c_i}\) 的第 \(\lfloor\frac{x}{c_i}\rfloor\) 行同乘 \(y\)，表示在这个因数下因为要满足另一行的所有的 \(gcd\) 均要乘 \(y\)。之后发现要么全部乘一行要么乘一列，则设 \(r_{i,j},c_{i,j}\) 表示 \(M_i\) 的第 \(j\) 行/列的系数，开始为 \(1\),则 \(M_i\) 的所有元素和为
\[\varphi(i)\times\sum r_{i,j} \times \sum c_{i,j} \]然后维护即可。
发现由于数据随机生成，则在 \([1,n]\) 中随机选取一个数的期望因子数位 \(\ln n^{(1)}\)，则时间复杂度为 \(O(n\ln n)\)，空间复杂度为 \(r,c\) 的空间为 \(O(n\ln n)\)。
(1)在 \([1,n]\) 中随机选取一个数的期望因子数位 \(\ln n\) 的证明：\(\sum_{i=1}^n d(i)=\sum_{i=1}^n\sum_{c|i}1=\sum_{c=1}^n\lfloor\frac{n}{i}\rfloor\approx n\ln n\)
则单一的期望值为 \(\ln n\)
代码
#include <bits/stdc++.h>
#define endl "\n"
using namespace std;
typedef long long ll;
const ll MAXN=1e5+5;
const ll MOD=1e9+7;
bool np[MAXN];
ll pm[MAXN],phi[MAXN];
vector<ll>r[MAXN],c[MAXN];
void pshuffle(){
phi[1]=1;
np[1]=true;
for(int i=2;i<MAXN;++i){
if(!np[i]){
pm[++pm[0]]=i;
phi[i]=i-1;
}
for(int j=1;j<=pm[0]&&i*pm[j]<MAXN;++j){
ll p=pm[j];
np[i*p]=true;
if(i%p==0){
phi[i*p]=phi[i]*p;
break;
}
phi[i*p]=phi[i]*(p-1);
}
}
}
vector<ll>V[MAXN];
ll rs[MAXN],cs[MAXN];
signed main(){
ios::sync_with_stdio(false);
cin.tie(0);cout.tie(0);
pshuffle();
ll n,q;
cin>>n>>q;
ll ans=0;
for(int i=1;i<=n;++i){
r[i].push_back(0);
c[i].push_back(0);
for(int j=1;j<=n/i;++j){
r[i].push_back(1);
c[i].push_back(1);
V[i*j].push_back(i);
}
ans+=phi[i]*(n/i)*(n/i);
rs[i]=cs[i]=n/i;
ans%=MOD;
}
while(q--){
char op;
ll x,y;
cin>>op>>x>>y;
if(op=='R'){
for(auto v:V[x]){
ans=((ans-phi[v]*rs[v]%MOD*cs[v]%MOD)%MOD+MOD)%MOD;
rs[v]+=(r[v][x/v]*(y-1)%MOD);
rs[v]%=MOD;
r[v][x/v]*=y;
r[v][x/v]%=MOD;
ans+=phi[v]*rs[v]%MOD*cs[v]%MOD;
ans%=MOD;
}
}else{
for(auto v:V[x]){
ans=((ans-phi[v]*rs[v]%MOD*cs[v]%MOD)%MOD+MOD)%MOD;
cs[v]+=(c[v][x/v]*(y-1)%MOD);
cs[v]%=MOD;
c[v][x/v]*=y;
c[v][x/v]%=MOD;
ans+=phi[v]*rs[v]%MOD*cs[v]%MOD;
ans%=MOD;
}
}
cout<<ans<<endl;
}
return 0;
}