python图示李亚普洛夫稳定性

 系统状态轨迹图，图一是不稳定系统，图二是稳定系统。

python代码

from scipy.integrate import odeint
from pylab import *

def deriv(y,t):        # 返回值是y和y的导数组成的数组
    # return array([y[0]+3*y[1], -5*y[0]+2*y[1]])
    # return array([4*y[0] - 2 * y[1], y[0] - 3 * y[1]])
    # return array([y[1], -y[0] - y[1]])
    return array([y[1]+y[0]*(2-y[0]**2-y[1]**2),-y[0]+y[1]*(2-y[0]**2-y[1]**2)])

time = linspace(0.0,50.0,1000)

figure()
for i in range(20):
    theta = i * pi / 10
    yinit = array([3*cos(theta), 3*sin(theta)])
    y = odeint(deriv,yinit,time)
    plot(y[:, 0], y[:, 1], linewidth=0.3)
    # quiver(y[:,0],y[:,1],gradient(y[:,0]),gradient(y[:,1]))

for i in range(20):
    theta = i * pi / 10
    yinit = array([0.03*cos(theta), 0.03*sin(theta)])
    y = odeint(deriv,yinit,time)from scipy.integrate import odeint
from pylab import *
 
def deriv(y,t):        # 返回值是y和y的导数组成的数组
    # return array([y[0]+3*y[1], -5*y[0]+2*y[1]])
    # return array([4*y[0] - 2 * y[1], y[0] - 3 * y[1]])
    # return array([y[1], -y[0] - y[1]])
    return array([y[1]+y[0]*(2-y[0]**2-y[1]**2),-y[0]+y[1]*(2-y[0]**2-y[1]**2)])
 
time = linspace(0.0,50.0,1000)
 
figure()
for i in range(20):
    theta = i * pi / 10
    yinit = array([3*cos(theta), 3*sin(theta)])
    y = odeint(deriv,yinit,time)
    plot(y[:, 0], y[:, 1], linewidth=0.3)
    # quiver(y[:,0],y[:,1],gradient(y[:,0]),gradient(y[:,1]))
 
for i in range(20):
    theta = i * pi / 10
    yinit = array([0.03*cos(theta), 3*sin(theta)])
    y = odeint(deriv,yinit,time)
    plot(y[:, 0], y[:, 1], linewidth=0.3)
    # quiver(y[:,0],y[:,1],gradient(y[:,0]),gradient(y[:,1]))
 
show()
    plot(y[:, 0], y[:, 1], linewidth=0.3)
    # quiver(y[:,0],y[:,1],gradient(y[:,0]),gradient(y[:,1]))

show()

总结

以前做控制理论仿真，首选都是matlab，现在工作用pythono，发现python真是无所不能，科学计算的能力也不遑多让matlab。不过就事论事，matlab轨迹图画的还是比python漂亮一点。