今回はSymPyを用いた方程式の解や連立方程式の解を求める方法について紹介する。
1. 方程式の定義方法
SymPyで方程式を定義するにはEq(左辺,右辺)により定義する。例えば
という方程式を定義する場合、左辺にa*x+b,右辺にcを入力する。
Eq(a*x+b,c) >Eq(a*x+b,c)
上は文字式を使用したが、数値を入れてもよい。
Eq(x+1.1,-1) >Eq(x+1.1,-1)
2. 方程式の解法
方程式を解くにはsolveset(方程式,解く変数)を用い、方程式はEq()を用いて定義する。例えば上の方程式を解くには
solveset(Eq(a*x+b,c),x)
>{-(b - c)/a}
またsolveset(方程式=0,x)としてもよい。
solveset(a*x+b-c,x)
>{-(b - c)/a}
またsolve(方程式,解く変数)で解くことができる。solveを用いる場合は方程式は左辺=0となるように定義する。上と同様に解くには
solve (a*x+b-c, x)
>{-(b - c)/a}
公式ドキュメント: https://docs.sympy.org/latest/tutorial/solvers.htmlには
Please note that there is another function called solve which can also be used to solve equations. The syntax is solve(equations, variables) However, it is recommended to use solveset instead.
とあり、solvesetの使用を推奨している。
2次方程式の場合、複数の解が表示される。
solveset(Eq(x**2,1),x) >{-1, 1}
文字式であっても、2次方程式の解の公式が正しく表示される。
solveset(Eq(a*x**2+b*x+c,0),x) >{-b/(2*a) - sqrt(-4*a*c + b**2)/(2*a), -b/(2*a) + sqrt(-4*a*c + b**2)/(2*a)}
解が複素数の場合では
solveset(x**2+1,x) >{-I, I}
解が存在しない方程式では、EmptySet()が返される。例えば
solveset(exp(x),x) >EmptySet()
三角関数を含む方程式で周期的な解が存在する場合は
solveset(sin(x)-1,x) ImageSet(Lambda(_n, 2*_n*pi + pi/2), S.Integers)
pi/2に対して2npiの周期で解があることを示している。
3. 連立方程式の解
連立方程式で線形方程式の場合はlinsolve([式1,式2],[変数1,変数2])などとすることで解くことができる。例えば
をx,yについて解くには
linsolve([a*x+b*y-e,c*x+d*y-f],[x,y])
>{((-b*f + d*e)/(a*d - b*c), (a*f - c*e)/(a*d - b*c))}
一方、連立非線形方程式の場合はnonlinsolve([式1,式2],[変数1,変数2])を用いる。例えば
nonlinsolve([x**2 + x, x - y], [x, y]) >{(-1, -1), (0, 0)}
解が虚数になるようなケースも適切に解が得られる。
nonlinsolve([x**2 + 1, y**2 + 1], [x, y]) >{(-I, -I), (-I, I), (I, -I), (I, I)}
周期関数など含む複雑な非線形連立方程式も解が存在すれば解ける
nonlinsolve([exp(x) - sin(y), 1/y - 3], (x, y)) >{(log(sin(1/3)), 1/3), (ImageSet(Lambda(_n, 2*_n*I*pi + Mod(log(sin(1/3)), 2*I*pi)), S.Integers), 1/3)}
なお非線形方程式に対してlinsolveを用いると解が存在してもEmptyで返してくるので注意が必要だ。
linsolve([x**2 + 1, y**2 + 1], [x, y]) >EmptySet()
4. まとめ
今回はSymPyを用いた方程式の解や連立方程式の解を求める方法について紹介した。
