如何用SymPy求解带初始条件的矩阵微分方程du/dt=[[0,1],[1,0]]u

我们可以通过以下步骤实现求解给定的矩阵微分方程：

1. 导入SymPy库并定义所需符号 #

首先导入SymPy核心模块，定义自变量t以及对应解向量分量的函数u₁(t)、u₂(t)：

import sympy as sp

# 定义符号和函数
t = sp.symbols('t')
u1 = sp.Function('u1')(t)
u2 = sp.Function('u2')(t)

2. 构建微分方程与系数矩阵 #

定义题目中的系数矩阵A，并写出矩阵微分方程对应的分量形式：

# 定义系数矩阵
A = sp.Matrix([[0, 1], [1, 0]])
# 定义解向量
u = sp.Matrix([u1, u2])
# 生成微分方程组：du/dt = A*u的分量形式
eqs = [sp.Eq(u1.diff(t), (A*u)[0]), sp.Eq(u2.diff(t), (A*u)[1])]

3. 设置初始条件 #

将题目给出的初始条件u(0) = [[4],[2]]转化为分量的初始条件：

ics = {u1.subs(t, 0): 4, u2.subs(t, 0): 2}

4. 求解微分方程组 #

使用SymPy的dsolve函数求解方程组，再代入初始条件确定常数：

# 求解方程组通解
sol = sp.dsolve(eqs)
# 根据初始条件求解常数
constants = sp.solve([sol[0].subs(t, 0).subs(ics), sol[1].subs(t, 0).subs(ics)])
# 代入常数得到特解
final_sol = [sol[0].subs(constants), sol[1].subs(constants)]

5. 整理结果并验证 #

将结果转换为矩阵形式，可验证与题目给出的解析解一致：

# 转换为矩阵形式
u_sol = sp.Matrix([final_sol[0].rhs, final_sol[1].rhs])
print(u_sol)

运行后输出结果：

Matrix([[3*exp(t) + exp(-t)], [3*exp(t) - exp(-t)]])

该结果可拆分为3e^t[[1],[1]] + e^{-t}[[1],[-1]]，与题目给出的解析解完全匹配。

内容的提问来源于stack exchange，提问作者Ben G