odeint 的非线性 ODE 解决方案
Non linear ODE solution by odeint
我要解下面的微分方程:
或
没有 F_1 项,代码很简单。但是我没能用包含的 F_1 项解决它,尽管我知道解决方案应该看起来像阻尼谐波振荡。
from scipy.integrate import odeint
import numpy as np
import matplotlib.pyplot as plt
def harmonic_motion(X,t,k,m,tau):
x,v=X
F=-(k/m)*x
F_1=tau/v*(np.gradient(x,t)**2) # This line doesn't work.
dx_dt=v
dv_dt=F-F_1
dX_dt=[dx_dt,dv_dt]
return dX_dt
m=1
k=1
tau=0.1
X0=-3
V0=0
t = np.linspace(0, 15, 250)
sol = odeint(harmonic_motion, [X0,V0], t,args=(k,m,tau))
x=sol[:,0]
v=sol[:,1]
plt.plot(t,x,label='x')
plt.plot(t,v,label='v')
plt.legend()
plt.xlabel('t (s)')
plt.ylabel('x,v (m,m/s)')
plt.show()
带有空气摩擦的 spring 的正确版本是
m*x'' + tau*abs(x')*x' + k*x = 0
这个的第一个订单系统是
def harmonic_motion(X,t,k,m,tau):
x,v = X
return v, -(tau*abs(v)*v + k*x)/m
这个相位图现在应该给出一个很好的原点螺旋。
我要解下面的微分方程:
或
没有 F_1 项,代码很简单。但是我没能用包含的 F_1 项解决它,尽管我知道解决方案应该看起来像阻尼谐波振荡。
from scipy.integrate import odeint
import numpy as np
import matplotlib.pyplot as plt
def harmonic_motion(X,t,k,m,tau):
x,v=X
F=-(k/m)*x
F_1=tau/v*(np.gradient(x,t)**2) # This line doesn't work.
dx_dt=v
dv_dt=F-F_1
dX_dt=[dx_dt,dv_dt]
return dX_dt
m=1
k=1
tau=0.1
X0=-3
V0=0
t = np.linspace(0, 15, 250)
sol = odeint(harmonic_motion, [X0,V0], t,args=(k,m,tau))
x=sol[:,0]
v=sol[:,1]
plt.plot(t,x,label='x')
plt.plot(t,v,label='v')
plt.legend()
plt.xlabel('t (s)')
plt.ylabel('x,v (m,m/s)')
plt.show()
带有空气摩擦的 spring 的正确版本是
m*x'' + tau*abs(x')*x' + k*x = 0
这个的第一个订单系统是
def harmonic_motion(X,t,k,m,tau):
x,v = X
return v, -(tau*abs(v)*v + k*x)/m
这个相位图现在应该给出一个很好的原点螺旋。