おはようございます．
今回もシルヴァーマン『はじめての数論』の演習問題を解いていきます．
前回は互除法をコードで書きました．
互除法は $a,b$ の最大公約数 $g$ を求めるアルゴリズムですが，少し工夫すれば
$$ ax + by = g $$
を満たすような $x,y$ まで同時に求めるアルゴリズムに変形することができます．
今回は，そういう問題です．
問題文
(a) から (c) まであるのですが，今回は (a) だけです．
図 6.1 に述べるアルゴリズムは，正の整数 $a$ と $b$ に対する最大公約数 $g$ と同時に方程式 $ax + by = \mathrm{gcd}(a,b)$ の解 $(x,y)$ も計算することを示せ．
図 6.1 というのは，次のようなものです．
解答
わからない
(7) のステップ (2) に戻れという命令は goto文ですが，わかりづらいので while 文で書き直します．
……うん．
書き直したはいいのですけど，これ何なんですかね．
眺めていれば何かひらめくかと思ってしばらく睨みつけていましたが，意味なかったです．なんもおもいつかん．
自分で一から考えてみる
いったん，このアルゴリズムを解釈することは後回しにして，自分でそういうアルゴリズムを考えてみることにします．
さて，$ax + by = g$ なる $x,y$ まで求めるにはどうしたらいいでしょうか．
$v,w$ という変数を $\mathrm{gcd}(a,b) = \mathrm{gcd}(v,w)$ がループごとに常に満たされるような変数とするときに，$v,w$ は $a,b$ の線形結合で表されているはずです．そこで，その係数を保持しておけばいいはずです．
というわけで，次のようなアルゴリズムを考えました．
うっかり$a,b$が正という条件を忘れてますが，まあいいでしょう．
青字はコメントです．invariant と書いてあるのは，ループ開始時に常に成り立つことを意味しています．progress と書いてあるのは，逆にループを繰り返すごとに狭義単調に変化することを指摘しています．
別にコメントは疑似コードに必ずしも書かなくてもよいのですが，（当社比で）格段に見やすくなるので付けておきました．
行列 $A$ を登場させました．これは，もともとの本に載っていたアルゴリズムで変数の数が多すぎてみづらかったので，変数の数を削減しようという意図でこういう書き方をしています．
再び本題のアルゴリズムへ
えー，それでは本題のアルゴリズムに戻ります．
相変わらず意味がわかりませんが，しかし $s = x - qv$ などの変数の設定の意図がわかるようになってきました．
あと，このアルゴリズムの意味がわからない理由がわかってきました．
各変数が何のための変数なのかがわかりにくいんですよ．while ループで保たれる条件がわかれば，アルゴリズムで何をやっているのかがもう少し明確になるでしょう．
さっきの自作のアルゴリズムと比較しつつぐっと睨むと，どうやら変数 $g,w$ は
$$ \mathrm{gcd}(a,b) = \mathrm{gcd}(g,w) $$
を常に（ループ開始時に）満たしていることがわかります．
そこで，さっきの自作アルゴリズムからの類推から $g,w$ は別扱いでよいことがわかります．$t$ も別扱いでよいです．こうした変数は，最大公約数を求める手続きに由来するものだからです．
残りの変数は $x,v,s,q$ です．$q$ は割った時の商なので，まあ良いでしょう．
$s$ ですが，これは一時的な変数ですね．わざわざ定義する必要のない変数なので，カットします．
さっきの自作のアルゴリズムと同様に，変数を行列にまとめて簡潔にしたいですね．そこでいろいろと試行錯誤したところ，
とすれば良い感じになることがわかりました．
最後に，invariant 条件(ループ開始時に常に成りたつ命題) を見つけます．
これは，演繹的に見つけることはたぶん不可能で，ぐっと睨むしかないと思うんですが，さっきの自作アルゴリズムがあるのでそれと比較することにしました．
大きな違いは，$y$ を最後に決めているところですね．
「$g - ax \in b \mathbb{Z}$ であるような $x$ さえ求まっていれば，$b$ で割れば $y$ が求まるでしょう？だから $x$ と $y$ を両方保持する必要はないよね」という思想で作られてるっぽいです．
ということは，invariant 条件も $b$ で割れるという条件である可能性が高いです．
そうやってアタリをつけて試したところ，どうも正解だったようです．
ここまでループ不変な条件がわかっていれば，アルゴリズムの正当性はもう明らかですね．ループごとに $w$ が減少していって，最終的に $g$ が最大公約数に一致します．そのときに $g- ax $ は $b$ の倍数だから，$y$ も自動的に求まりますというだけのことですね．
これで解けました．
感想とまとめ
ユークリッドの互除法を適用した後，後退代入によって $ax + by = g$ なる $x,y $ が求まるというのは有名な話です．
しかし，それをコードで書くというのはやったことがなかったので今回は（珍しく）いい勉強になりました．
いつもこのブログの演習問題記事は自明な問題ばっかりやってたんですけど，今回はそれなりに非自明な問題に取り組めたんじゃないかと思います．
（演習問題は，非自明なやつは解けなさすぎてブログの更新が滞りそうで心配なんですよね．それでつい自明なやつばっかりやってしまうんですよね…．）
あと，今回のテーマである「疑似コードを見て，どんなアルゴリズムか推測する」というのはなかなか愉しい問題だと思いました．
ループ不変な条件とかループの進捗測定器が書いてあれば一発ですが，書かれていなかった場合に推測するのは直観勝負になりそうでおもしろいです．
それこそ Youtuber とかでもできそうな企画だと思いました．
このブログでもまたやってみようかな？
