Skip to content

Instantly share code, notes, and snippets.

@Yunaka12
Created May 12, 2019 04:10
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save Yunaka12/f5ea9a04119b5c77d24bbd9bc1614e95 to your computer and use it in GitHub Desktop.
Save Yunaka12/f5ea9a04119b5c77d24bbd9bc1614e95 to your computer and use it in GitHub Desktop.
2変数関数の極大・極小・鞍点をsympyを利用して求める
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]]))
@variable-funciton
Copy link

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

何が原因かよくわからないのですが,教えていただけませんか?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment