【SymPy】（七）方程求解（八）矩阵

【SymPy】（一）SymPy简介

【SymPy】（二）使用SymPy需要避开的坑

【SymPy】（三）基本操作（四）打印

【SymPy】（五）简化

【SymPy】（六）微积分

（七）方程求解

from sympy import *
x, y, z = symbols('x y z')

1 等式

从前面的文章我们知道，SymPy中的符号方程不是用=或==表示的，而是用Eq表示的。

然而，还有一个更简单的方法：我们知道a=b相当于a−b=0。在SymPy中，任何不在Eq中的表达式都被SynPy的求解函数（如solveset()）假定等于0。这意味着我们可以不使用x==y，而使用x-y（被求解函数认为是x-y==0）。

如果要求解的方程已经等于0，则这一点特别有用。我们在编程时不用键入Eq：solveset(Eq(expr,0),x)，只需使用solveset(expr,x)即可求解方程。
例如：

2 求解代数方程

求解代数方程的主要函数是solveset。其语法为：solveset(equation, variable=None, domain=S.Complexes)

式中，等式可以是Eq实例或表达式的形式（假定为零）。

还有一个函数叫做solve，也可以用来解方程。它的语法为solve(equations, variables)，后面会介绍它的用途。

当求解一个等式时，solveset的输出有：有限集，域，映射集

如果无解，则返回一个空集，如果无法找到解决方案，则返回一个条件集合。

在solveset模块中，使用linsolve求解线性方程组。下面是linsolve的语法示例：

• （1） 方程式形式:

• 增广矩阵形式：

• $Ax=b$形式
  

注：解的顺序对应于给定符号的顺序。

在solveset模块中，使用nonlinsolve求解非线性方程组。

下面是nonlinsolve的示例。

• 当只有实数解存在时：

• 当只有复数解存在时：

• 实数解和复数解同时存在时

• 无穷多个解时：

---

注：
（1）解的顺序对应于给定符号的顺序。

（2）nonlinsolve不会返回LambertW形式的解（如果存在LambertW形式解的话）。
solve则可用于这种情况：

（3）nonlinsolve不能很好地求解具有三角函数的方程组。
solve也可用于此类情况（但不能给出所有解决方案）

solveset只给出每个解。要得到包含多重性的多项式的解，使用roots。

roots的输出{0: 1, 3: 2}意思是：0是重数为1的根，3是重数为2的根。

solveset 无法求解由LambertW（超越方程求解器）求解的方程。
而solve可以：

---

3 求解微分方程

要求解微分方程，可使用dsolve。首先，通过将cls=function传递给symbols函数来创建未定义的函数。

f, g = symbols('f g', cls=Function)

f和g现在是未定义的函数。我们可以调用f(x)，它表示一个未知函数。

f(x)的未估值导数：

要求解常微分方程，请将其和要求解的函数传递给dsolve。

如求解微分方程$f^{\prime \prime }(x)-2f^{\prime }(x)+f(x)=\sin (x)$：

dsolve返回Eq的一个实例。这是因为一般情况下，微分方程的解不能为函数显式求解。

dsolve解中的任意常数是C1、C2、C3等形式的符号。

（八）矩阵

from sympy import *

要在SymPy中生成矩阵，使用Matrix对象。通过提供构成矩阵的行向量的列表来构造矩阵。

例如，构造矩阵：
$$\left[\begin{array}{cc}1& -1\\ 3& 4\\ 0& 2\end{array}\right]$$

为了便于生成列向量，单独一个元素列表被认为是列向量

矩阵的操作与SymPy或Python中的任何其他对象一样。

关于SymPy矩阵需要注意的一点是，与SymPy中的其他对象不同，它们是可变的。这意味着可以对它们进行适当的修改。如果需要Matrix的不可变，请使用ImmutableMatrix。

1 基本操作

1.1 获取矩阵形状

要得到矩阵的形状，请使用shape：

1.2 访问行和列

若要获取矩阵的单个行或列，请使用row或col。例如，M.row(0)将获取第一行。M.col(-1)将获取最后一列。

1.3 删除和插入行和列

要删除行或列，请使用row_del或r col_del。这些操作将修改矩阵。

要插入行或列，请使用 row_insert或col_insert。

2 基本方法

简单的加法和乘法操作仅通过使用+、*和**来完成。要求矩阵的逆，则计算-1次方。

求Matrix的转置，用T：

3 矩阵构造函数

要创建单位矩阵，请使用eye。eye(n)将创建一个n×n单位矩阵。

zeros(n, m)创造一个 n × m 的零矩阵

相似地，ones创建全1矩阵

要创建对角矩阵，请使用diag。diag的参数可以是数字或矩阵。其中数字被解释为1×1矩阵。输入的矩阵按对角线排列。其余元素用0填充。

4 先进的方法

4.1 行列式

要计算矩阵的行列式，请使用det。

4.2 简化的阶梯形式

若要将矩阵转换为简化的阶梯形式，请使用rref

rref返回两个元素的元组。第一个元素是简化的行梯队形式，第二个元素是主元的列索引。

注：rref返回的元组的第一个元素是Matrix类型。第二种是tuple类型。

4.3 零空间（Nullspace）

要找到矩阵的零空间，请使用nullspace，得到输出张成矩阵零空间的列向量。

使用ASCII Pretty Printer格式输出：

4.4 列空间（Columnspace）

要找到矩阵的零空间，请使用columnspace，得到输出张成矩阵列空间的列向量。

4.5 特征值、特征向量和对角化（Eigenvalues, Eigenvectors, and Diagonalization）

要求矩阵的特征值，请使用eigenvals。eigenvals返回一个字典：{特征值：代数重数}

表示M特征值有-2、3和5，特征值-2和3的代数重数为1，特征值5的代数重数2。

要找到矩阵的特征向量，请使用eigenvects。返回(特征值：代数重数，[特征向量])。

从结果我们也可以知道特征值5的几何重数为2，因为它有两个特征向量。因为所有特征值的代数重数和几何重数是相同的，所以M是可对角化的。

使用diagonalize实现对角化，返回一个元组$(P,D)$，D为对角矩阵且$M=PD{P}^{-1}$

Tip：要创建一个名为λ的符号，可以对SymPy符号和Python变量使用相同的名称——lamda（不带b）。它可以被打印为λ。

如果只需要特征多项式，请使用charpoly。这比eigenvals更有效，因为有时符号根的计算成本很高。

未完待续：

【SymPy】（九）高级表达式操作