当我使用scipy求解常微分方程时,获得的结果与sympy的结果不同。具体来说,它们在小数点后各不相同,并且没有虚部。
由于我要求解的常微分方程是sympy获得的第一个符号解,因此我首先使用sympy.lambdify将其转换,然后将其引入scipy.odeint方法。
这是代码的一部分
t = sy.symbols('t')
y_0_t_func = sy.Function('y_0')
tspan = np.linspace(0,10,10)
def y_0_Q_x(t):
return sy.nsimplify(
(lambda_1 * pow(A_0,2) * F_0 * P_theta *
(k_0 + beta * (function_B_0(t) / function_B_1(t)) - F_0 *P_theta)
- (lambda_2 / 2) * pow(function_B_0(t) / function_B_1(t),2) +
y_1_t.rhs * (alpha - F_0 * P_theta + x_1_t.rhs)))
f_y_0 = sy.lambdify((y_0_t_func(t),t),r * y_0_t_func(t) - y_0_Q_x(t))
y_0_numerical = np.flipud(
integrate.odeint(f_y_0,g_2 * S_Ba,np.flipud(tspan)))
y_0_t = sy.dsolve(
sy.nsimplify((sy.diff(y_0_t_func(t),t,1) - r * y_0_t_func(t)
+y_0_Q_x(t))),y_0_t_func(t),ics={y_0_t_func(10): g_2 * S_Ba})
sy.pprint(y_0_t.rhs.evalf(subs={t: 0}))
print(y_0_numerical[0])
程序的结果是
-951.483500299349 - 2.98155597433514e-19⋅ⅈ
-951.48350659