Mathematicaでグレブナー基底を使って連立方程式を解いてみる
Mathematicaでグレブナー基底を使って連立方程式を解いてみる
◼
執筆者:技術部・降旗
◼
掲載日:
◼
使用したバージョン:Mathematica14.0.0
◼
使用した主な関数:GroebnerBasis, Solve
はじめに
はじめに
普段Mathematicaで連立方程式を解きたい場合にはSolveやNSolveなどの関数を使いますが、先日グレブナー基底を使うことでも連立方程式を解くことができることを知りました。そこで今回はグレブナー基底を使って実際に連立方程式を解いてみたいと思います。
Mathematicaでグレブナー基底を求める
Mathematicaでグレブナー基底を求める
グレブナー基底がどのようなものか、また、グレブナー基底を使って連立方程式を解く方法の詳細については、知らない人向けに分かりやすく説明しているブログなどをウェブ上で簡単に見つけることができるためここでは省略しますが、グレブナー基底の発想の原点は「多項式同士の除算計算」の一意性にあり、連立多項式を解く際に解が存在するとしたらグレブナー基底は解に至る手順に一意性を保証するそうです。
Mathematicaでグレブナー基底を求めるにはGroebnerBasis関数を使用します。
Mathematicaでグレブナー基底を求めるにはGroebnerBasis関数を使用します。
2元1次連立方程式
2元1次連立方程式
まずは下記の簡単な2元1次連立方程式を解いてみたいと思います。
x+3y=7 |
2x+y=4 |
この方程式の解は x=1, y=2 です。
In[]:=
ContourPlot[{x+3y==7,2x+y==4},{x,-5,5},{y,-5,5}]
Out[]=
In[]:=
Solve[{x+3y==7,2x+y==4},{x,y}]
Out[]=
{{x1,y2}}
グレブナー基底を求めるために、まずは右辺の値を左辺に移行することにより連立方程式から多項式を求めます。
x+3y-7 |
2x+y-4 |
そして、この多項式を GroebnerBasis関数に入力して実行します。 GroebnerBasis関数の1つ目の引数は多項式のリスト、2つ目の引数は変数のリストです。
In[]:=
GroebnerBasis[{x+3y-7,2x+y-4},{x,y}]
Out[]=
{-2+y,-1+x}
なお、多項式の形にすることは必須ではないようで、連立方程式の形のまま入力しても実行できました。
In[]:=
GroebnerBasis[{x+3y==7,2x+y==4},{x,y}]
Out[]=
{-2+y,-1+x}
今回は1次の方程式のため、得られた結果からすぐに y=2, x=1 ということが分かりますが、求まったグレブナー基底に==0を付けて、これらの解が求められることを確かめてみます。
In[]:=
Solve[-2+y==0]
Out[]=
{{y2}}
In[]:=
Solve[-1+x==0]
Out[]=
{{x1}}
2元2次連立方程式
2元2次連立方程式
続いて2元2次連立方程式についても同様にグレブナー基底を使って解いてみたいと思います。
この連立方程式の解は4つあります。
この連立方程式の解は4つあります。
2 x 2 y |
2 x |
In[]:=
ContourPlot[{x^2-4y^2==4,x^2-2x-2y==16},{x,-10,10},{y,-10,10}]
Out[]=
In[]:=
sol=Solve[{x^2-4y^2==4,x^2-2x-2y==16},{x,y}]
Out[]=
x,y-16-2+,x,y-16-2+,x,y-16-2+,x,y-16-2+
1
2
2
1
2
2
1
2
2
1
2
2
結果はRootの形で返されるため少し分かりづらいです。
おおよその値が知りたい場合にはN関数を使います。
おおよその値が知りたい場合にはN関数を使います。
In[]:=
N[sol]
Out[]=
{{x-3.45113,y1.40626},{x-2.86614,y-1.02649},{x4.58765,y-2.06437},{x5.72961,y2.68461}}
それでは、先ほどと同様に多項式の形にしてからグレブナー基底を求めます。
In[]:=
gb=GroebnerBasis[{x^2-4y^2-4,x^2-2x-2y-16},{x,y}]
Out[]=
{32+12y-27-4+4,6+x+y-2}
2
y
3
y
4
y
2
y
この連立方程式では2つの基底が求まります。
少し話が脱線しますが、求まったグレブナー基底を==0と置いてプロットしてみると、グラフの形はかなり違いますが、解の位置は変わっていないように見えることが分かります。
少し話が脱線しますが、求まったグレブナー基底を==0と置いてプロットしてみると、グラフの形はかなり違いますが、解の位置は変わっていないように見えることが分かります。
In[]:=
ContourPlot[{32+12y-27-4+4==0,6+x+y-2==0},{x,-10,10},{y,-10,10},Epilog->{AbsolutePointSize[5],Point[{x,y}]/.sol},PlotLegends->"Expressions"]
2
y
3
y
4
y
2
y
Out[]=
さらに元の方程式と重ねて表示させてみます。
下の例では元の2つの方程式を赤とマゼンタ、求まった2つの基底を青とシアンで描画します。
下の例では元の2つの方程式を赤とマゼンタ、求まった2つの基底を青とシアンで描画します。
In[]:=
ContourPlot[{x^2-4y^2==4,x^2-2x-2y==16,32+12y-27-4+4==0,6+x+y-2==0},{x,-10,10},{y,-10,10},Epilog->{AbsolutePointSize[5],Point[{x,y}]/.sol},ContourStyle->{Red,Magenta,Blue,Cyan},PlotLegends->"Expressions"]
2
y
3
y
4
y
2
y
Out[]=
本題に戻り、先ほどと同様に求まったグレブナー基底から方程式の解を求めてみます。
求まった1つ目の基底はyのみに関する4次の多項式のため、この式を==0と置いた方程式を解くことでyを求めることができます。
求まった1つ目の基底はyのみに関する4次の多項式のため、この式を==0と置いた方程式を解くことでyを求めることができます。
4つの解が求まりましたがSolveの結果と一致する結果が無いように見えます。
しかし、FullSimplifyを使ってSolveの結果を簡約化すると、順序は異なりますが同じ結果が得られていることを確認できます。
しかし、FullSimplifyを使ってSolveの結果を簡約化すると、順序は異なりますが同じ結果が得られていることを確認できます。
求まった2つ目の基底はxとyに関する式のため、この式を==0と置いた方程式に求まったyを代入してから方程式を解くことで、対応するxの値が求められます。試しに得られた1つ目のyを代入してみます。
こちらも最初のSolveで得られた解と一致するものが無いように見えますが、FullSimplifyを使って簡約化すると同じ解が得られます。
同様に他のyに対応するxの値を求めてみます。
得られた結果をまとめます。
おわりに
おわりに
今回はMathematicaでグレブナー基底を使って連立方程式の解を求めてみました。もちろん最初から Solve を使って直接連立方程式の解を求める方が簡単なため、今後の課題としてはグレブナー基底を用いた場合にメリットがあるか(グレブナー基底を使って解くことで早く解が求まる問題があるのか)などを調べてみたいと思います。