- 前回のあらすじ
- 最初の不等式が意味すること（part 1)
- 二つ目の不等式が意味すること
- デカルト座標か極座標か？（最初の不等式が意味すること：part 2)
- 動径$r$の範囲の目処
- pythonのコード
- 付録（ちょっとした改良）
前回のあらすじ
2024年の東京大学の入試問題（数学）を解き始めた。小手調べに問題１からとりかかっている。下の２つの不等式で表される条件を満たす(x,y)を求める問題なので、シミュレーションでまずは答えの概略を得ようという方針である。
\begin{equation} -1 \le -\frac{y}{\sqrt{2(x^2 + y^2)}} \le -\frac{1}{2} \end{equation} および \begin{equation} \frac{\sqrt{3}}{2}\le\frac{y+2}{\sqrt{2(x^2+y^2+2y+2)}}\le 1 \end{equation}
最初の不等式が意味すること（part 1)
シミュレーションの前に、まずは最初の不等式から分析を始める。まず、$y<0$の領域が「はじかれる」ことがすぐにわかる。 納得いかない人は、例えば$y=-1$などの具体的な値を代入してみるといいだろう。
実は最初の不等式はもうちょっと有用な情報をもっている。これについては後でもう一度考える。
二つ目の不等式が意味すること
次に二つの不等式は$x\leftrightarrow -x$の変換（パリティ反転のようなもの）に対して対称性をもつ点である。これは点Aがyz平面にあるので、幾何学的にみてもyz平面を対称面として角度$\theta,\phi$に対称性があるからである。つまり$x>0$の領域も$x<0$の領域も、x軸の矢印の向きが変わるだけで幾何学的な関係は不変だということである（言葉で書くと仰々しいが、絵を描いてみれば一発でわかるだろう）。ということで、答えもy軸に対して左右対称の図形になるはずである。
ということで、シミュレーションでは$x>0, y>0$、つまり第一象限だけ計算すればよく、最後にy軸に対してパタンと対称的に折り返せば答えが得られる。
また条件(i)で「点Pは原点Oとは異なる」とあるので$x=y=0$は領域から除くことなる点も覚えておこう。実際、最初の不等式に$x=y=0$を代入すると面倒なことが起こりそうな感じがある（原点は境界を表す二つの直線の交点に対応し、境界線上にある）。一方、二つ目の不等式は満たされてしまう（楕円の内部という条件になるので問題なく満たされてしまう）。
デカルト座標か極座標か？（最初の不等式が意味すること：part 2)
シミュレーションでは、$x>0, y>0$を満たす全ての点$(x,y)$について具体的な値を使って２つの不等式を計算し、条件が満たされる場合は⚪︎、満たされない場合は×でラベルすることになる。そして、最後に⚪︎と判定された座標だけを描画するが完了である。
この「絨毯爆撃」は稠密に行いところであるが、数値計算では離散点を選んで「サンプリング」しないといけない。ある意味、この「サンプリング」はスクラッチカードをコインでこするような作業であるから、このタイプのシミュレーションを「スクラッチカードシミュレーション」と今後呼ぶことにしよう。
さて、計算のサンプリング点は、デカルト座標に沿って「格子状」に行うのが基本であるが、時には極座標に沿って回転するようなサンプリングをした場合がわかりやすい場合もある。果たして、この問題ではどちらがいいだろうか？
答えは「どちらでもいい」であるが、最初の不等式は極座標を使うと簡単になるので、今回は極座標を選んでみることにする。 $x=r\cos\alpha, y=r\sin\alpha$を代入すると、最初の不等式は次のように書き換えられる。 \begin{equation} \frac{1}{\sqrt{2}} \le \sin\alpha \le \sqrt{2} \end{equation} 不等式の右側の「端」は明らかに満たされる(正弦関数は１より小さいから）。したがって、この不等式は \begin{equation} \frac{\pi}{4}\le \alpha \le \frac{\pi}{2} \end{equation} と解釈できる。もちろん、$x\rightarrow -x$の対称性があるので、本当のところは \begin{equation} \frac{\pi}{4} \le \alpha \le \frac{3\pi}{4} \end{equation} という意味である。これはすなわち、$y>x$(ただし$y>0$)とその（y軸に関する）折り返し領域が適用範囲です、ということを言っているのである。
2つ目の不等式に移ろう。極座標で条件を書き直すと、 \begin{equation} \frac{\sqrt{3}}{2}\le \frac{r\sin\alpha+2}{\sqrt{2(r^2+2r\sin\alpha+2)}} \le 1 \end{equation} となる。
動径$r$の範囲の目処
シミュレーションでは、$r$の大きさを次第に大きくしながら、与えた$r$に対して$\alpha$を回していくという感じなる。つまり、一番上のループは$r$に関してであり、二つ目のループが$\alpha$に関してとなる。$r$の大きさをどの程度にするかは、試行錯誤でもよいが、一応の目処はつけておいた方がいいだろう。
まずは$r\rightarrow\infty$の場合を見てみたい。果たして第１象限を目一杯考えなければならないのであろうか？ということである。 \begin{equation} \frac{r\sin\alpha+2}{\sqrt{2(r^2+2r\sin\alpha+2)}} \rightarrow \frac{\sin\alpha}{\sqrt{2}} < \frac{1}{\sqrt{2}}=\frac{\sqrt{2}}{2} \end{equation} なので、これは条件の範囲外である。つまり、あまり遠くにある点(x,y)は考えなくていい、ということだ。
次は$\alpha=\pi/2$の方向について調べてみよう。代入すると \begin{equation} \frac{\sqrt{3}}{2}\le \sqrt{\frac{r+2}{2(r^2+2r+2)}} < 1 \end{equation} となる。このままではよくわからないので、二乗して2次式の形に書き直す。不等式は二つに分割されるが、整理すると \begin{equation} r^2>0 \end{equation} \begin{equation} r^2-2r-2 \le 0 \end{equation} となるため、最初の条件が自明なものとなって不要である。したがって、二つ目の不等式だけが重要となり、$r>0$に注意すると \begin{equation} 0<r\le 1+\sqrt{3} < 3 \end{equation} を得る。
次は$\alpha=\frac{\pi}{4}$の場合について計算してみる。不等式を整理すると \begin{equation} 2r^2-\sqrt{2}r-2\le 0 \rightarrow (r-\sqrt{2})(r+\frac{1}{\sqrt{2}})\le 0 \end{equation} を得るので、 \begin{equation} 0<r\le \sqrt{2} \end{equation} となる。
ということで、だいたい$r=3$程度、まあ$r=5$まで計算しておけば、最初のシミュレーションとしては十分そうである。
pythonのコード
以前、$0<\log_x y<1$を調べた時に使ったコードを再利用しよう。今回は、関数を定義して、不等式の真ん中の複雑な表式の計算に利用する。
import math def f1(r,q): return (r*math.sin(q) + 2.0)/math.sqrt(2*(r**2+2*r*math.sin(q)+2.0)) Nmax = 100 Rmin = 0.0 Rmax = 5.0 Amin = math.pi/4 Amax = math.pi/2 dr = (Rmax - Rmin) / Nmax da = (Amax - Amin) / Nmax for j in range(1, Nmax+1): r = Rmin + j*dr for i in range(1, Nmax+1): alph = Amin + i*da result = f1(r,alph) if (result > math.sqrt(3)/2.) and (result < 1.): print(i,j,r*math.cos(alph),r*math.sin(alph),result)
面倒臭いので、角度と動径に関するメッシュ数は共通(Nmax=100)にしてある。もちろん、問題によっては離散化の仕方は異なっていい。 今回はとにかく計算するだけ。答えを外部ファイルに出力すると、こんな感じになる。
1 1 0.035076571288542795 0.035631925946260266 0.9995501892587716 2 1 0.03479563982961572 0.035906314888159443 0.999552675108376 3 1 0.03451256201172186 0.03617848895942246 0.9995551575191258 4 1 0.03422735529643444 0.03644843137107058 0.9995576359191038 5 1 0.03394003727664709 0.03671612547178428 0.9995601097375942 6 1 0.033650625675488675 0.03698155474893048 0.9995625784052117 7 1 0.03335913834522999 0.0372447028295811 0.9995650413540283 ...
この数値データをgnuplotを使って表示すると次のような図となった。
$y=x$と$y=\sqrt{3}+1$の直線を加えてある。
なんだか指紋のような感じの図形である。答えを知っていれば、まあこれで良さそうなのはわかるのだが、このシミュレーションの結果から楕円の方程式にたどり着くことはできるのであろうか？
付録（ちょっとした改良）
計算点をすべて描画し、適合した点と適合してない点で色分けして描画してみる。条件文の結果を次のように分ける。つまり適合したときは最後の数字を１に、適合しないときは最後の数字を０にして、外部ファイルに書き出す。
if (result > math.sqrt(3)/2.) and (result < 1.): print(i,j,r*math.cos(alph),r*math.sin(alph),result,1) else: print(i,j,r*math.cos(alph),r*math.sin(alph),result,0)
外部ファイルをawkで処理して、適合したもの(tmp1.dat)と適合しないもの(tmp0.dat)にファイルを分けて、gnuplotで描画する。
awk '{if($6 == 0) print $3,$4}' tmp.dat> tmp0.dat
awk '{if($6 == 1) print $3,$4}' tmp.dat> tmp1.dat
gnuplot> plot [0:5][0:5] "tmp0.dat" w p pt 6 ps .5, "tmp1.dat" w p pt 6 ps 0.5,x,1+sqrt(3)
結果は次の通り。
