概要
公式解説の、$\Theta(H\log H)$ よりも計算量がよいので、書きます。Project Euler でよく見る手法のひとつですね。正方形では何度も経験があるけれど、長方形では初めて実装しました。
解法
ひとまず真っすぐ横・縦のものは別処理。右下がりのものは右上がりのものと同数ある。などをやって、右上がりのものの計算に帰着します。$(0,0)$ と $(x,y)$ を結ぶ線分の平行移動たちをまとめて数えることで、次のような問題に帰着されます。
$A(H,W) = \{(x,y)\mid 1\leqq x\leqq H, 1\leqq y\leqq W, \gcd(x,y) = 1\}$ とするとき、$\sum_{x,y\in A(H,W)} (H-x)(W-y)$ を求めよ。
展開することで、だいたい次の問題だということが分かります。
$$\sum_{x,y\in A(H,W)} 1,\qquad \sum_{x,y\in A(H,W)} x, \qquad \sum_{x,y\in A(H,W)} y, \qquad \sum_{x,y\in A(H,W)} xy$$を求めよ。
この $4$ つの計算は、だいたい同じように解くことができます。$\sum xy$ のみ説明しましょう。
$H \geqq W$ を仮定しておきます。（$\max(H,W)$ などが無駄に出てきて邪魔なので）
$f(H,W) = \sum_{x,y\in A(H,W)}xy$ とおきます。$A(H,W)$ に関する集計は、$\gcd(x,y) = 1$ の条件が邪魔で行いにくいです。そこで、$B(H,W) = \{(x,y)\mid 1\leqq x\leqq H, 1\leqq y\leqq W\}$ とおき、$F(H,W) = \sum_{x,y\in B(H,W)}xy$ としましょう。
矩形領域 $B(H,W)$ はとても扱いが簡単で、$F(H,W) = (\sum_{x=1}^H x)(\sum_{y=1}^Wy)$ により $F$ の値は $1$ 件あたり $\Theta(1)$ で計算できます。
一方 $(x,y)\in B(H,W)$ のうち、$\gcd(x,y) = d$ を満たすものをまとめて集計すると見なすことで、次が成り立つことが分かります：
$$F(H,W) = \sum_{d=1}^{H}e^2f\biggl(\biggl\lfloor\frac{H}{d}\biggr\rfloor, \biggl\lfloor\frac{W}{d}\biggr\rfloor\biggr)$$
ここから、$f(H,W) = F(H,W) – \sum_{d=2}^{H} d^2f\bigl(\bigl\lfloor\frac{H}{d}\bigr\rfloor, \bigl\lfloor\frac{W}{d}\bigr\rfloor\bigr)$ により、$f$ を小さい側から計算していくことができます。
これだけでは計算量が減らないように思いますが、次の $2$ 点から、このような計算は全件合わせて たとえば $\Theta(H^{3/4})$ で行うことができます：
・ひとつの $f(h,w)$ を求める計算を考える。$h$ 件の $f\bigl(\bigl\lfloor\frac{h}{d}\bigr\rfloor, \bigl\lfloor\frac{w}{d}\bigr\rfloor\bigr)$ の和を求める計算がある。しかしここに現れるペア $\mathrm{pair}(d) = \bigl(\bigl\lfloor\frac{H}{d}\bigr\rfloor, \bigl\lfloor\frac{W}{d}\bigr\rfloor\bigr)$ は、$d\leqq \sqrt{h}$ または、$\bigl\lfloor\frac{H}{d}\bigr\rfloor + \bigl\lfloor \frac{W}{d} \bigr\rfloor \leqq 2\sqrt{h}$ を満たすことなどから、$\Theta(\sqrt{h})$ 件しか存在しない。
・必要なすべての $f(h,w)$ の計算を考える。これも上と同様に、必要な $(h,w)$ の計算が $\Theta(\sqrt{h})$ 件しか存在しないことが分かる。全体の計算量は、$\int_{x=1}^{\sqrt{H}} \sqrt{H/x}dx$ の定数倍により評価できて、$\Theta(H^{3/4})$ となる。
