如果这篇博客帮助到你，可以请我喝一杯咖啡~
CC BY-NC-SA 4.0 （除特别声明或转载文章外）
第一次见这种题，写篇题解。
我们先考虑消元后的结果，$A_{i,j}=ij\gcd(i,j)$，设消元后的矩阵为 $A’$，设 $B_{i,j}=\gcd(i,j)$，则 $A’{i,j}=ijB’{i,j}$，这部分的证明就是对第 $i$ 行做除以 $i$ 的行变换，对第 $j$ 列做除以 $j$ 的列变换，进行消元后变换回去。
那么现在的问题就是求解 $B’$，这里放张 $6\times 6$ 的 $B$ 和 $B’$，供大家理解。
\[B=\begin{bmatrix}1&1&1&1&1&1\\1&2&1&2&1&2\\1&1&3&1&1&3\\1&2&1&4&1&2\\1&1&1&1&5&1\\1&2&3&2&1&6\end{bmatrix}\\ B'=\begin{bmatrix}1&1&1&1&1&1\\0&1&0&1&0&1\\0&0&2&0&0&2\\0&0&0&2&0&0\\0&0&0&0&4&0\\0&0&0&0&0&2\end{bmatrix}\\\]我们注意到一点，在消元时，第 $i$ 行只会让 $i$ 的倍数行减去本行。因为其他行的第 $i$ 列的元素已经为 $0$ 了。
根据打表，我们推测 $B’_{i,j}=\begin{cases} 0 & i\nmid j\\ \varphi(i) & i\mid j \end{cases}$。
设矩阵行数为 $n$，当 $n=1$ 时猜测显然成立，我们假设矩阵行数为 $n-1$ 时猜测成立。
若 $n\nmid k$，显然有 $\gcd(n,k)\mid n$，在消元到 $\gcd(n,k)$ 行时，前面的所有行对 $B_{n,k}$ 和 $B_{\gcd(n,k),k}$ 的影响相同，又因为 $\gcd(n,k)=\gcd(\gcd(n,k),k)$。所以此时 $B_{n,k}$ 和 $B_{\gcd(n,k),k}$ 相等，此行消元完成后 $B’_{n,k}=0$。
若 $n\mid k$：
\[\begin{aligned} B'_{n,k}=&\gcd(n,k)-\sum\limits_{d\mid n,d\ne n}\varphi(d)\\ =&\gcd(n,k)+\varphi(n)-\sum\limits_{d\mid n}\varphi(d)\\ =&n+\varphi(n)-n\\ =&\varphi(n) \end{aligned}\]所以若行数为 $n-1$ 行时成立，则行数为 $n$ 时也成立，又因为行数为 $1$ 时成立，所以对于 $\ge 1$ 的所有行数都成立。
所以 $A’_{i,j}=\begin{cases} 0 & i\nmid j \\ ij\varphi(i) & i\mid j \end{cases}$。
所以操作 $4$ 的答案为 $\prod\limits_{i=1}^xi^2\varphi(i)$。
消元前：
操作 $2$：$xy\gcd(x,y)$。
操作 $3$：$\sum\limits_{i=1}^x\sum\limits_{j=1}^xij\gcd(i,j)$。
\[\begin{aligned} &\sum\limits_{i=1}^x\sum\limits_{j=1}^xij\gcd(i,j)\\ =&\sum\limits_{i=1}^x\sum\limits_{j=1}^xij\sum\limits_{d\mid i,d\mid j}\varphi(d)\\ =&\sum\limits_{d=1}^x\varphi(d)d^2\sum\limits_{i=1}^{\left\lfloor\frac{x}{d}\right\rfloor}i\sum\limits_{j=1}^{\left\lfloor\frac{x}{d}\right\rfloor}j\\ =&\sum\limits_{d=1}^x\varphi(d)d^2\left(\frac{\left\lfloor\frac{x}{d}\right\rfloor(\left\lfloor\frac{x}{d}\right\rfloor+1)}{2}\right)^2 \end{aligned}\]消元后：
操作 $2$：$\begin{cases} 0 & x\nmid y \\ xy\varphi(x) & x\mid y \end{cases}$。
操作 $3$：$\sum\limits_{j=1}^x\sum\limits_{i\mid j}ij\varphi(i)$。
\(\begin{aligned} &\sum\limits_{j=1}^x\sum\limits_{i\mid j}ij\varphi(i)\\ =&\sum\limits_{j=1}^xj\sum\limits_{i\mid j}i\varphi(i)\\ \end{aligned}\) $f(i)=i\varphi(i)$ 的 DGF 为 $\frac{\zeta(x-2)}{\zeta(x-1)}$，$g(j)=\sum\limits_{i\mid j}f(i)$ 的 DGF 为 $\frac{\zeta(x)\zeta(x-2)}{\zeta(x-1)}$，$h(j)=jg(j)$ 的 DGF 为 $\frac{\zeta(x-1)\zeta(x-3)}{\zeta(x-2)}$，线性筛下这个函数就好了，这个函数在 $p^k(p\in \mathrm{Prime})$ 处的取值为 $\frac{p^{3k+1}+p^k}{p+1}$。