数值方法导论/线性方程组

目标

• 定义向量和矩阵
• 加、减、乘矩阵
• 求方阵的转置
• 定义矩阵的逆
• 将联立线性方程组设置为矩阵形式

资源

• 书本章节
• Python编程

向量和矩阵

一个矩阵 是一个矩形的数组，例如数字、符号或表达式。矩阵通常用来表示线性变换 和线性方程组。

一个三角矩阵 是一种特殊的方阵。如果 A 中主对角线以下的所有元素都为零，则 A 称为上三角矩阵。

$a_{ij}=0$ 对所有 $i>j$

类似地，如果 A 中主对角线以上的所有元素都为零，则 A 称为下三角矩阵。

$a_{ij}=0$ 对所有 $i<j$

如果主对角线以外的所有元素都为零，则 A 称为对角矩阵。

$a_{ij}=0$ 对所有 $i\neq j$

这个矩阵

${\begin{bmatrix}1&4&2\\0&3&4\\0&0&1\\\end{bmatrix}}$

是上三角，这个矩阵

${\begin{bmatrix}1&0&0\\2&8&0\\4&9&7\\\end{bmatrix}}$

是下三角。

大小为 n 的单位矩阵 是一个 n×n 矩阵，其中主对角线上的所有元素都等于 1，其他所有元素都等于 0。

矩阵乘法

矩阵乘法 接受两个矩阵并生成另一个矩阵。该操作的规则如下所示

${\overset {4\times 2{\text{ matrix}}}{\begin{bmatrix}{\color {Brown}{a_{11}}}&{\color {Brown}{a_{12}}}\\\cdot &\cdot \\{\color {Orange}{a_{31}}}&{\color {Orange}{a_{32}}}\\\cdot &\cdot \\\end{bmatrix}}}{\overset {2\times 3{\text{ matrix}}}{\begin{bmatrix}\cdot &{\color {Plum}{b_{12}}}&{\color {Violet}{b_{13}}}\\\cdot &{\color {Plum}{b_{22}}}&{\color {Violet}{b_{23}}}\\\end{bmatrix}}}={\overset {4\times 3{\text{ matrix}}}{\begin{bmatrix}\cdot &x_{12}&x_{13}\\\cdot &\cdot &\cdot \\\cdot &x_{32}&x_{33}\\\cdot &\cdot &\cdot \\\end{bmatrix}}}$

用圆圈标记的交点处的值为

${\begin{aligned}x_{12}&={\color {Brown}{a_{11}}}{\color {Plum}{b_{12}}}+{\color {Brown}{a_{12}}}{\color {Plum}{b_{22}}}\\x_{13}&={\color {Brown}{a_{11}}}{\color {Violet}{b_{13}}}+{\color {Brown}{a_{12}}}{\color {Violet}{b_{23}}}\\x_{32}&={\color {Orange}{a_{31}}}{\color {Plum}{b_{12}}}+{\color {Orange}{a_{32}}}{\color {Plum}{b_{22}}}\\x_{33}&={\color {Orange}{a_{31}}}{\color {Violet}{b_{13}}}+{\color {Orange}{a_{32}}}{\color {Violet}{b_{23}}}\end{aligned}}$

矩阵的转置

矩阵的转置定义如下：矩阵A的转置是另一个矩阵A^T，可以通过以下任何一种等效操作创建

• 沿其主对角线（从左上角到右下角）反射A以获得A^T
• 将A的行写为A^T的列
• 将A的列写为A^T的行

下图说明了这个概念

矩阵的逆

如果存在一个n×n方阵B，使得

$\mathbf {AB} =\mathbf {BA} =\mathbf {I} _{n}\$

其中I_n表示n×n单位矩阵，使用的乘法是普通的矩阵乘法。B被称为A的逆，记为A⁻¹。

将线性方程组表示为矩阵形式

${\begin{alignedat}{7}2x&&\;+\;&&y&&\;-\;&&z&&\;=\;&&8&\qquad (L_{1})\\-3x&&\;-\;&&y&&\;+\;&&2z&&\;=\;&&-11&\qquad (L_{2})\\-2x&&\;+\;&&y&&\;+\;&&2z&&\;=\;&&-3&\qquad (L_{3})\end{alignedat}}$

我们可以将这个线性方程组表示成矩阵形式，如下所示

${\begin{bmatrix}2x&+&y&-&z\\-3x&-&y&+&2z\\-2x&+&y&+&2z\end{bmatrix}}={\begin{bmatrix}8\\-11\\-3\end{bmatrix}}$

然后我们可以使用矩阵乘法来分离变量。

${\begin{bmatrix}2&1&-1\\-3&-1&2\\-2&1&2\end{bmatrix}}{\begin{bmatrix}x\\y\\z\end{bmatrix}}={\begin{bmatrix}8\\-11\\-3\end{bmatrix}}$

${\begin{bmatrix}25&5&1\\64&8&1\\144&12&1\end{bmatrix}}{\begin{bmatrix}a_{1}\\a_{2}\\a_{3}\end{bmatrix}}={\begin{bmatrix}106.8\\177.2\\279.2\end{bmatrix}}$

线性方程组的一般矩阵形式如下所示

$A_{n\times n}X_{n\times 1}=C_{n\times 1}$

矩阵A是系数矩阵。X是解向量（矩阵），C是右侧向量。如果我们在两边乘以A的逆矩阵，我们可以看到解与A的逆矩阵密切相关。

${\begin{alignedat}{3}A^{-1}AX&&\;=\;&&A^{-1}C\\IX&&\;=\;&&A^{-1}C\\X&&\;=\;&&A^{-1}C\end{alignedat}}$

高斯消元法

来源

• 书本章节

高斯消元法 是一个用来求解线性方程组的算法，类似于求一个可逆方阵的逆矩阵。该算法包含一系列对系数矩阵进行的行变换操作。共有三种基本行变换：

• 交换两行
• 将一行乘以一个非零常数
• 将一行乘以一个常数加到另一行上。

例如，第一个线性方程 (L₁，主元方程) 可以用来消去 $x$ 从接下来的两个方程中：

${\begin{alignedat}{7}2x&&\;+\;&&y&&\;-\;&&z&&\;=\;&&8&\qquad (L_{1})\\-3x&&\;-\;&&y&&\;+\;&&2z&&\;=\;&&-11&\qquad (L_{2})\\-2x&&\;+\;&&y&&\;+\;&&2z&&\;=\;&&-3&\qquad (L_{3})\end{alignedat}}$

然后 L₂ (主元方程) 可以用来消去 $y$ 从 L₃ 中。这个过程被称为**前向消元**。现在 L₃ 中只有一个未知量 $z$，可以用来将 $z$ 代入 L₂ 中求解 $y$。这个过程被称为**后向代入**。

高斯消元法的算法可以这样实现：

''' 
x = gauss_elimin(a, b)
Solves [a][x] = [b] by Gauss elimination.
'''
from numpy import dot, array

def gauss_elimin(a, b):
  (rows, cols) = a.shape
  # elimination phase
  for row in range(0, rows-1): # pivot equation/row
    for i in range(row+1, rows):
      if a[i, row] != 0.0:
        factor = a[i, row]/a[row, row]
        a[i, row+1:rows] = a[i, row+1:rows] - factor*a[row, row+1:rows]
        b[i] = b[i] - factor*b[row]
  # back substitution 
  for k in range(rows-1,-1,-1):
    b[k] = (b[k] - dot(a[k, k+1:rows],b[k+1:rows]))/a[k, k]
  return b

a = array([[3, 2.0], [-6, 6]])
b = array([7, 6])
print gauss_elimin(a, b)

带部分主元选择的高斯消元法

高斯消元法容易受到舍入误差的影响，特别是当方程数量很大时，误差会累积。另一个问题是除零错误。例如，以下矩阵会导致该方法失效，因为第二个主元元素为零，而该方法要求主元元素非零。

$\left[{\begin{array}{ccc|c}1&-1&2&8\\0&0&-1&-11\\0&2&-1&-3\end{array}}\right]$

带部分主元选择的高斯消元法类似于高斯消元法，只是在 $i$ 次前向消元步骤中，找出 $|a_{i,i}|,|a_{i+1,i}|,...,|a_{n,i}|$ 的最大值，并将包含最大值的所在行与主元行进行交换。交换后，之前的系数矩阵如下所示。

$\left[{\begin{array}{ccc|c}1&-1&2&8\\0&2&-1&-3\\0&0&-1&-11\end{array}}\right]$

高斯-赛德尔法

高斯消元法是一种直接（直接）方法，它将原始方程转化为等价的更容易求解的方程。一些方程组没有解，例如，方程数量少于未知量，或者一个方程与另一个方程矛盾。

高斯-赛德尔法 是一种迭代（或间接）方法，它从对解的猜测开始，并反复改进猜测，直到它收敛（收敛条件得到满足）。这种方法的优点是，当系数矩阵很大且稀疏时，它在计算上更有效率，并且可以通过调整误差容限来控制舍入误差。

该算法包含以下步骤：

1. 将每个方程改写为对应未知量的形式，例如，将第一个方程写成 $x_{1}$ 在左侧，将第二个方程写成 $x_{2}$ 在左侧

$x_{i}={\frac {c_{i}-\sum _{j=1,j\neq i}^{n}a_{ij}x_{j}}{a_{ii}}},i=1,2,3,...,n.$

1. 找到 $x_{i}$ 的初始猜测。
2. 使用重写的方程计算新的估计值 - 始终使用最新的估计值。
3. 重复上一步，直到所有 $x_{i}$ 中最大的绝对相对近似误差小于指定的容差。

大多数迭代方法的缺点是，当方程组的系数矩阵不是 对角占优 时，它可能不会收敛。系数矩阵对角占优的方程组总是收敛的。矩阵 A 对角占优，如果

$|a_{ii}|\geq \sum _{j\neq i}|a_{ij}|\quad {\text{for all }}i,$

以及

$|a_{ii}|>\sum _{j\neq i}|a_{ij}|\quad {\text{for at least one }}i.$

一个 使用 numpy 的 Python 实现 只需迭代即可生成解向量。

LU 分解

LU 分解 将系数矩阵 A 分解为下三角矩阵和上三角矩阵的乘积：A = LU。从 A 中导出 L 和 U 的过程称为 LU 分解或 LU 分解，类似于 高斯消元法。

完成 LU 分解后，我们可以使用前向和后向替换来求解由 A 表示的线性方程组。

${\begin{alignedat}{3}AX&&\;=\;&&C\\LUX&&\;=\;&&C\\L^{-1}LUX&&\;=\;&&L^{-1}C\\UX&&\;=\;&&L^{-1}C\\\end{alignedat}}$

令 $L^{-1}C=Z$，这意味着 $LZ=C$ 和 $Z$ 可以通过前向消元法求解。从 $UX=Z$ 中，我们可以使用后向消元法求解 X。

LU 分解的一种方法类似于高斯消元法。对于一个 $3\times 3$

$A={\begin{bmatrix}a_{11}&a_{12}&a_{13}\\a_{21}&a_{22}&a_{23}\\a_{31}&a_{32}&a_{33}\\\end{bmatrix}}={\begin{bmatrix}1&0&0\\l_{21}&1&0\\l_{31}&l_{32}&1\\\end{bmatrix}}{\begin{bmatrix}u_{11}&u_{12}&u_{13}\\0&u_{22}&u_{23}\\0&0&u_{33}\\\end{bmatrix}}.$

$A={\begin{bmatrix}u_{11}&u_{12}&u_{13}\\u_{11}l_{21}&u_{12}l_{21}+u_{22}&u_{13}l_{21}+u_{23}\\u_{11}l_{31}&u_{12}l_{31}+u_{22}l_{32}&u_{13}l_{31}+u_{23}l_{32}+u_{33}\\\end{bmatrix}}$

从最后一个矩阵中，我们可以看到为了消除 $a_{21}$，我们需要将 A 的第一行乘以 $l_{21}$，类似地，为了消除 $a_{31}$，需要将第一行乘以 $l_{31}$。为了计算 L，我们可以简单地记录正向消元步骤中使用的乘数，矩阵 U 与正向消元得到的最终上三角矩阵相同。

Python 示例代码

求矩阵的逆

求矩阵的逆是一个与求解线性方程组相关的过程。如果找到了系数矩阵的逆，那么由系数矩阵表示的方程组也就被解出来了。

${\begin{alignedat}{3}AX&&\;=\;&&C\\A^{-1}AX&&\;=\;&&A^{-1}C\\X&&\;=\;&&A^{-1}C\\\end{alignedat}}$

求矩阵的逆的一种方法是通过解一组线性方程组分别求出逆矩阵的每一列。如果 $AB=I$ 并且

$A={\begin{bmatrix}a_{11}&a_{12}&a_{13}\\a_{21}&a_{22}&a_{23}\\a_{31}&a_{32}&a_{33}\end{bmatrix}}\ B={\begin{bmatrix}b_{11}&b_{12}&b_{13}\\b_{21}&b_{22}&b_{23}\\b_{31}&b_{32}&b_{33}\end{bmatrix}}$

那么

$A{\begin{bmatrix}b_{11}\\b_{21}\\b_{31}\end{bmatrix}}=LU{\begin{bmatrix}b_{11}\\b_{21}\\b_{31}\end{bmatrix}}={\begin{bmatrix}1\\0\\0\end{bmatrix}}$.

我们可以通过使用 LU 分解方法求解方程组来计算 $B$ 的第一列。其余的列可以以类似的方式计算。正如你所看到的，LU 分解方法的优势在于 LU 分解只执行一次，结果被多次使用，这降低了整个解决方案的计算复杂度。