Maxima で綴る数学の旅

紙と鉛筆の代わりに、数式処理システムMaxima / Macsyma を使って、数学を楽しみましょう

-数学- ハートのえくぼ、代数曲線の孤立点

パピヨンと秋

 

以前掲載した記事:

に、コメントを頂きました(コメントは出来れば普通に書いて頂いて良いのですが、、、)。コメントの内容を要約すると、(%o1)の方程式で表される代数曲線はハートの形を描くが、実は見えない特異点(孤立点)が含まれる。それを、この代数曲線上の点(x,y)のxの最大値、最小値、yの最大値、最小値を求めることで見つけよ。

(%i1) H:x^6 + 3*x^4*y^2 + 6*x^4*y - 2*x^4 + 3*x^2*y^4 - 2*x^2*y^3 - 6*x^2*y^2 - 6*x^2*y + 3*x^2 + y^6 - 3*y^4 + 3*y^2 - 1 = 0;
$$ \tag{%o1} y^6+3\,x^2\,y^4-3\,y^4-2\,x^2\,y^3+3\,x^4\,y^2-6\,x^2\,y^2+3\,y^2+6\,x^4\,y-6\,x^2\,y+x^6-2\,x^4+3\,x^2-1=0 $$

f:id:jurupapa:20161207232942p:plain 

問題が、多変数多項式の方程式を制約とする最大最小問題に定式化されると、当然使う道具はQuantifier Elimination (限定子除去)がぴったりです。

まずxの最大、最小を求めてみます。グラフを見ると大体-1<=x<=1にみえます。実際はどうなのでしょうか。

(%i2) load(qepmax)$

QEを使う場合、yが(実数に)存在するためのxの条件を求めれば良いのです。
(%i3) HX:qe([[E,y]],H);
$$ \tag{%o3} \left(64\,x^6-219\,x^4+192\,x^2-64\leq 0\right) \wedge \left(\left(64\,x^6-219\,x^4+192\,x^2-64=0\right) \lor \left(\left(x-1\leq 0\right) \wedge \left(x+1\geq 0\right)\right)\right) $$

計算時間は2秒くらいでした。予想した答え-1<=x<=1だけではなく、xの6次式が付いています。この論理式はこのままでは分かりにくいので、Aをこのxの6次式と置いて、論理式を簡略化してみます。残念ながらここは手計算で行いました。

(%i4) A=part(HX,1,1);
$$ \tag{%o4} A=64\,x^6-219\,x^4+192\,x^2-64 $$
(%i5) A<=0 %and (A=0 %or (-1<=x %and x<=1));
$$ \tag{%o5} \left(-1\leq x\right) \wedge \left(x\leq 1\right) \lor \left(A=0\right) \wedge \left(A\leq 0\right) $$

丁寧に論理式を変形していくと、結局、(%o6)のようにA=0あるいは-1<=x<=1、となります。
(%i6) A=0 %or (-1<=x %and x<=1);
$$ \tag{%o6} \left(-1\leq x\right) \wedge \left(x\leq 1\right) \lor \left(A=0\right) $$

このA=0を満たすxがどんな値なのか、具体的に求めてみましょう。
(%i7) A:part(HX,1,1);
$$ \tag{%o7} 64\,x^6-219\,x^4+192\,x^2-64 $$
(%i8) HXS:solve(A=0,x);
$$ \tag{%o8} \left[ x=-\frac{\sqrt{\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+1233}}{8\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{6}}} , x=\frac{\sqrt{\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+1233}}{8\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{6}}} , x=-\frac{\sqrt{\sqrt{3}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\,i-137\,3^{\frac{5}{2}}\,i-\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}+146\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}-1233}}{2^{\frac{7}{2}}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{6}}} , x=\frac{\sqrt{\sqrt{3}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\,i-137\,3^{\frac{5}{2}}\,i-\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}+146\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}-1233}}{2^{\frac{7}{2}}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{6}}} , x=-\frac{\sqrt{-\sqrt{3}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\,i+137\,3^{\frac{5}{2}}\,i-\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}+146\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}-1233}}{2^{\frac{7}{2}}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{6}}} , x=\frac{\sqrt{-\sqrt{3}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\,i+137\,3^{\frac{5}{2}}\,i-\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}+146\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}-1233}}{2^{\frac{7}{2}}\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{6}}} \right] $$

非常に複雑形でxが求まりました。これを見易くするために、浮動小数点形式に変換して見ます。

(%i9) %,numer,rectform;
$$ \tag{%o9} \left[ x=-1.520184372286849 , x=1.520184372286849 , x=-0.226227439061609\,i-0.778868468556658 , x=0.226227439061609\,i+0.778868468556658 , x=0.226227439061609\,i-0.778868468556658 , x=0.778868468556658-0.226227439061609\,i \right] $$

最初の2つだけが実数であとは複素数です。

最初のxに対応するyを求めて見ましょう。元の式のxに6次方程式の解の一つを入れて、yのみの方程式を作ります。

(%i10) H,HXS[1];
$$ \tag{%o10} y^6+\frac{3\,\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)\,y^4}{64\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}}-3\,y^4-\frac{\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)\,y^3}{32\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}}+\frac{3\,\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)^2\,y^2}{4096\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}}-\frac{3\,\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)\,y^2}{32\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}}+3\,y^2+\frac{3\,\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)^2\,y}{2048\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}}-\frac{3\,\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)\,y}{32\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}}+\frac{\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)^3}{262144\,\left(13824\,\sqrt{17}+71577\right)}-\frac{\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)^2}{2048\,\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}}+\frac{3\,\left(1233+73\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}+\left(71577+13824\,\sqrt{17}\right)^{\frac{2}{3}}\right)}{64\,\left(71577+13824\,\sqrt{17}\right)^{\frac{1}{3}}}-1=0 $$
(%i11) %,numer;

複雑すぎてよくわからないので、係数を浮動小数点形式に変換して見ます。
$$ \tag{%o11} y^6+3.93288157723548\,y^4-4.62192105149032\,y^3+5.155852500186077\,y^2+18.17746815484311\,y+7.593578252988485=0 $$

この方程式の近似解を求めてみます。
(%i12) allroots(%);
$$ \tag{%o12} \left[ y=-0.6610497889701161 , y=-0.661049816831282 , y=2.289943689914807\,i-0.6610498029006993 , y=-2.289943689914807\,i-0.6610498029006993 , y=1.144971844957403\,i+1.322099605801399 , y=1.322099605801399-1.144971844957403\,i \right] $$

この最初の2つが実数解で、ほとんど同じ値です。グラフを書くとわかるのですが、おそらくこれらは同じ値で重解ではないか、と思います。とすると、これで (x,y)=(-1.520184372286849,-0.6610497889701161)及び(1.520184372286849,-0.6610497889701161)という2つの孤立点が見つかったことになります。

同じようにしてxが存在することを条件としてQEすることで、yの範囲が求まります。あとは同じような計算をすることであと2つほど孤立点を見つけることができます。

このやり方(x, yの最大、最小を求める方法)では全ての孤立点を求められる保証はありません。

少し文献をネットで調べてみたところ、数式処理システムRisa/Asirのグラフ描画ルーチンでは全ての孤立点をきちんと描画することが出来るようです。

学位論文 竹島卓 数式処理システムRisa/Asirの開発と応用 神戸大学、2005年

の第5章が参考になります。