Created
May 12, 2019 04:10
-
-
Save Yunaka12/f5ea9a04119b5c77d24bbd9bc1614e95 to your computer and use it in GitHub Desktop.
2変数関数の極大・極小・鞍点をsympyを利用して求める
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
from sympy import Symbol, solve, Derivative, Matrix, simplify | |
x = Symbol('x') | |
y = Symbol('y') | |
f = 2*x**3 + 4*x*y**2 - 10*x*y + y**2 | |
f_x = Derivative(f,x).doit() #xで1階偏微分 | |
f_y = Derivative(f,y).doit() #yで1階偏微分 | |
stationary_points = solve([f_x,f_y],[x,y]) #関数の停留点 | |
stationary_points = list(map(lambda x:simplify(x),stationary_points)) #簡単化 | |
f_xx = Derivative(f, x, 2).doit() #xで2階偏微分 | |
f_yy = Derivative(f, y, 2).doit() #yで2階偏微分 | |
f_xy = Derivative(f_x,y).doit() #xで偏微分してからyで偏微分 | |
f_yx = Derivative(f_y,x).doit() #yで偏微分してからxで偏微分 | |
#ヘッセ行列 | |
hesse = Matrix([ | |
[f_xx,f_xy], | |
[f_yx,f_yy] | |
]) | |
#ヘッセ行列の行列式 | |
det_hessian = hesse.det() | |
#とりあえず出力 | |
print('f_x: {0}'.format(f_x)) | |
print('f_y: {0}'.format(f_y)) | |
print('f_xx: {0}'.format(f_xx)) | |
print('f_yy: {0}'.format(f_yy)) | |
print('f_xy: {0}'.format(f_xy)) | |
print('f_yx: {0}'.format(f_yx)) | |
print('全ての停留点は {0}'.format(stationary_points)) | |
#判定 | |
for i in range(len(stationary_points)): | |
sign_det_hessian = det_hessian.subs([(x,stationary_points[i][0]),(y,stationary_points[i][1])]).evalf() | |
if sign_det_hessian > 0: | |
sign_f_xx = f_xx.subs([(x,stationary_points[i][0]),(y,stationary_points[i][1])]).evalf() | |
if sign_f_xx > 0: | |
print('停留点:{0}はf(x,y)の極小点'.format([stationary_points[i][0],stationary_points[i][1]])) | |
elif sign_f_xx < 0: | |
print('停留点:{0}はf(x,y)の極大点'.format([stationary_points[i][0],stationary_points[i][1]])) | |
else: | |
print('この方法では判定不可能') | |
elif sign_det_hessian < 0: | |
print('停留点:{0}はf(x,y)の鞍点'.format([stationary_points[i][0],stationary_points[i][1]])) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2変数関数 f=x2+x*y+y2として, 上のコードをjupiter notebookで走らせたところ,次のエラーが発生しました:
TypeError Traceback (most recent call last)
in
37 #判定
38 for i in range(len(stationary_points)):
---> 39 sign_det_hessian = det_hessian.subs([(x,stationary_points[i][0]),(y,stationary_points[i][1])]).evalf()
40 if sign_det_hessian > 0:
41 sign_f_xx = f_xx.subs([(x,stationary_points[i][0]),(y,stationary_points[i][1])]).evalf()
TypeError: 'Symbol' object is not subscriptable
何が原因かよくわからないのですが,教えていただけませんか?