python求解常微分方程

### Python求解常微分方程 #### 使用`SciPy`库中的`solve_ivp` 对于初值问题，可以利用`SciPy`库里的`solve_ivp`函数来解决。此函数提供了多种可选的积分器，并支持事件检测功能。 ```python from scipy.integrate import solve_ivp import matplotlib.pyplot as plt import numpy as np def model(t, y): return -2 * t * y # dy/dt=-2ty 的例子 sol = solve_ivp(model, [0, 4], [1], dense_output=True) t = np.linspace(0, 4, 100) plt.plot(t, sol.sol(t).T) plt.xlabel('Time') plt.ylabel('Solution y') plt.title('ODE Solution using solve_ivp') plt.show() ``` 这段代码展示了如何设置模型函数并调用`solve_ivp`来进行求解[^1]。 #### 利用`SymPy`解析求解 当面对能够找到闭合形式解的情况时，可以选择采用`SymPy`进行符号运算得到精确的结果。 ```python from sympy import * init_printing(use_unicode=True) f = Function('f')(x) # 定义未知函数f(x) eq = Eq(f.diff(x), f) # 构建一阶线性齐次微分方程df/dx=f solution = dsolve(eq, f) pprint(solution) # 打印漂亮格式化的解答C₁⋅ℯˣ ``` 这里通过定义符号表达式的办法构建了简单的指数增长型的一阶线性齐次微分方程，并得到了它的通解[^3]。 #### 数值方法实现——四阶龙格-库塔法 除了依赖现成工具外，还可以手动编写程序模拟一些经典的数值算法过程，比如下面展示的是四阶Runge-Kutta方法的一个简单版本： ```python def runge_kutta_4th_order(func, init_cond, time_span, step_size): """Four-order Runge Kutta method.""" ts = np.arange(*time_span, step_size) ys = np.zeros_like(ts, dtype=float) ys[0] = init_cond for i in range(len(ts)-1): k1 = func(ts[i], ys[i]) k2 = func(ts[i]+step_size/2., ys[i]+k1*step_size/2.) k3 = func(ts[i]+step_size/2., ys[i]+k2*step_size/2.) k4 = func(ts[i]+step_size, ys[i]+k3*step_size) ys[i+1] = ys[i] + (k1 + 2*k2 + 2*k3 + k4)*step_size / 6. return ts, ys # 测试案例 ts, ys = runge_kutta_4th_order(lambda t, y: -2*t*y, 1., (0., 4.), .1) plt.figure(figsize=(8,5)) plt.plot(ts, ys, label='Approximation by RK4') plt.legend(); plt.grid(True); plt.show(); ``` 该片段实现了四阶龙格-库塔公式以近似特定类型的连续时间动态系统的响应曲线[^4]。