東大数学2024問題１(part 3): 境界を抽出する - 複天一流：どんな手を使ってでも問題解決を図るブログ

前回のあらすじ
境界を抜き出す
データの読み込み
動径と角度でグループ化されたデータ構造を利用して境界を抽出する
境界を抽出する

前回のあらすじ

前回の記事では、２つの不等式を満たす$(x,y)$をしミューレションによって数値的に見つけ出し、その結果をgnuplotで描画した。大体の感じはわかったが、その境界の方程式を算出するところには至らなかった。

今回は不等式を満たす領域の境界を表す方程式を計算して求めるところまでやりたい。

境界を抜き出す

ある領域の境界を抜き出すという作業は、デジタル画像データの処理でよく利用される。画像の輪郭を抜き出して、個人認証に利用したり、ロボットやAIが作業する物体の識別に利用したりするからだ。今回のプログラムはそこまで大げさなものではないが、基本的な考え方は同じであると思われるので、いづれ役に立つ日もくるのではないだろうか？（補償なしの、安っぽい希望で申し訳ない）

前回書いたプログラムをscratch-card.pyと名づけることにすると、このプログラムの出力で6つ目の数字（整数）が1となっている座標（データの３つ目と４つ目の浮動小数の数字）が不等式の条件を満たす点である。この領域の「境界」だけを抜き出すためには、極座標で計算したループの構造を利用することになる。

ループでは動径$r$を与えてから、角度$\alpha$を動かした。したがって、与えられた$r$のうち、$\alpha$の値が最も小さい点に対応するデータが、境界点を与えることになる。この観点からpythonでプログラムを書いてみた。

i=1
while True:
    try:
        data_seq = input().split()
        if int(data_seq[1]) == i and int(data_seq[5]) == 1:
            i += 1
            print(data_seq[0],data_seq[1],data_seq[2], data_seq[3],data_seq[5])

    except EOFError:
        break

このプログラムはread-data.pyと呼ぶことにする。以下では、このプログラムを部分分解して、ポイントを確認してみる。

データの読み込み

読み込むデータはscratch-card.pyで作ったデータである。したがって、次の二つの方法で読み込むことができる。

[1] 外部データに一度scratch-card.pyの結果を保存するやり方
% python3 scratch-card.py > tmp.dat
% python3 read-data.py < tmp.dat > result.dat

[2]パイプを使ってデータをプログラム間でやりとりするやり方
% scratch-card.py | read-data.py > result.dat

実際の計算では、１つ目の方法を採用することが多いだろう。２つ目の方法では、プログラムを実行形式にした場合を想定している。

今回のデータ読み取り方法の利点は、データの型を気にしなくていいという点である。それはPythonがもっているinput()関数のおかげである。この関数はデータを「文字列」と解釈して取り込み、「文字列」として出力する。

以前、赤玉計算の問題を解いたときはC言語を使ったのでscanf()を利用したが、今回は同じような処理をpythonのinput()を使って行っている。すでに上述したとおり、scanfの場合には読み込むデータの型を指定する必要があったのに対し、inputの場合にはどんなデータでも「文字列」と決めつけて読み込むので、変数型について気にする必要がない。その分だけプログラムは簡略化されるが、その後の計算では好みのデータ型へ変換する点について気をつけねばならない。

またinputで読み込んだデータは「区切り」というものがない。scratch-card.pyでは「空白」を用いて複数の数値データを区切っていたが、inputで読み込んでしまうと空白も文字列の一部と見做される。

したがって、たとえば、

@@XOOO,XOOX@

というデータがあるとする（赤玉の時に使った記号である）。2つの記号シーケンスがカンマ(,)で区切られているというデータ構造である。これはいわゆるcsvファイルである（comma separated values)。

このデータを次のようにしてpythonで読み取ったとする。

dseq = input()

するとdseq[0]=@ dseq[1]=@ dseq[2]=X dseq[3]=O .... dseq[6]=, dseq[7]=X ... dseq[11]=@となる。またlen(dseq)は12という値を返す。データの区切りのつもりでつかったカンマ(,)もデータの一部として取り込まれてしまう。

データが2つのシーケンスからなることをプログラムにわからせるために利用するのがsplit()関数である。今、csv形式でデータを区切っているので次のようにプログラムしてデータを２つに分割する。

dseq = input().split(',')

すると, len(dseq)は2となり、dseq[0]=@@XOOO dseq[1]=XOOX@となって、人間が意図した構造にデータが分割された。

splitのデフォルトの区切り文字は「空白」なので、scratch-card.pyのデータを読み込むときはsplit()をそのまま使用することになる。

最後に、scratch-card.pyが出力するデータ行は100である。予め行数がわかっているので、その分だけ読み込むようにプログラムすることは可能だが、C言語のようにEOFを利用し、任意の長さのデータを読み込めるようにしたいときもあるだろう。今回のread-data.pyでは後者を採用している。そこで利用するのがtry... except EOFError...の構文である。コードを見れば使い方はすぐにわかるだろう。

以上で、データの読み込み部分についての説明はおしまいである。

動径と角度でグループ化されたデータ構造を利用して境界を抽出する

さて、プログラムのメインである境界抽出のアルゴリズム（のようなもの）についてみていこう。「アルゴリズム」と呼ぶのがちょっとおこがましい感じもするくらい、今回のプログラムでは簡単に処理できる。

まず、動径$r$と角度$\alpha$によってデータはまとめられており、まずは$r$を固定した上で、$\alpha$を動かしている。つまり、与えられた$r$に対して最も小さな$\alpha$の値をとる点が「境界点」となる。ある$r$で境界点を見つけたら、次の$r$に移り、そこでまた境界点を見つける、といった作業を（$r$の最大値になるまで）繰り返せばよい。

条件を満たす点かどうかは6番目の整数値data_seq[5]で判定できる(1なら条件を満たす点、0なら条件を満たさない点）。

整数の変数iが、着目している動径$r$の値に対応するラベルである。最初に合致したら、すぐに次の値i+1になってしまうため、領域内部の座標は表示されずスルーしていくだけである。データが進んでdata_seq[1]の値が変化したタイミングこそが境界点であり、そのときは座標が表示される....これをデータの終わりまで繰り返す、という感じである。

東大の問題を自力で解けばわかるように、動径が小さいところ$(r\le \sqrt{2})$では境界は$y=x$で表される。この部分はデータを「見る」だけで人間でも判定できる。問題は、動径が大きくなって境界の方程式が楕円のそれに切り替わる領域である。答えは \begin{equation} x^2 + \frac{(y-1)^2}{3} = 1 \end{equation} つまり、 \begin{equation} y=\sqrt{3(1-x^2)}+1 \end{equation} であるのだが、これをシミュレーションで計算するにはどうすればよいかである。