实践：使用NumPy求解方程组

线性方程组可以表示为 $Ax = b$ 的形式。恒等矩阵和逆矩阵 ($A^{-1}$) 是解决这类方程组的重要概念。此处，将使用NumPy来求得解向量 (vector) $x$。

我们将通过NumPy的两种主要方式来实现这一点：

1. 直接计算逆矩阵 ($A^{-1}$)，然后计算 $x = A^{-1}b$。此方法紧密遵循数学定义，但通常效率较低，并且对于计算机而言，数值稳定性可能较差，特别是对于大型系统。理解此方法很有帮助。
2. 使用NumPy的专用求解函数 (np.linalg.solve)。这是实践中推荐的方法。它通常更快且数值更可靠，因为它使用精密算法（通常基于矩阵分解，如LU分解），这些算法避免显式地形成逆矩阵。

让我们通过一个例子进行学习。考虑以下线性方程组：

$$2x_1 + x_2 - x_3 = 8 \\ -3x_1 - x_2 + 2x_3 = -11 \\ -2x_1 + x_2 + 2x_3 = -3$$

我们的目标是找到同时满足所有这三个方程的 $x_1$、$x_2$ 和 $x_3$ 的值。

在NumPy中设置系统

首先，我们需要将这个系统表示为矩阵形式 $Ax = b$。

矩阵 $A$ 包含变量的系数：

$$A = \begin{bmatrix} 2 & 1 & -1 \\ -3 & -1 & 2 \\ -2 & 1 & 2 \end{bmatrix}$$

向量 (vector) $b$ 包含右侧的常数项：

$$b = \begin{bmatrix} 8 \\ -11 \\ -3 \end{bmatrix}$$

向量 $x$ 包含我们想要求解的未知变量：

$$x = \begin{bmatrix} x_1 \\ x_2 \\ x_3 \end{bmatrix}$$

让我们将 $A$ 和 $b$ 创建为NumPy数组：

import numpy as np

# 定义系数矩阵 A
A = np.array([
    [2, 1, -1],
    [-3, -1, 2],
    [-2, 1, 2]
])

# 定义常数向量 b
# 我们使用 reshape(-1, 1) 使其成为一个显式的列向量（3行，1列）
# 尽管NumPy通常会处理广播，但显式定义是良好的做法。
b = np.array([8, -11, -3]).reshape(-1, 1)

print("矩阵 A:\n", A)
print("\n向量 b:\n", b)

方法1：使用矩阵逆求解

从数学角度看，如果 $A$ 可逆，我们可以通过计算 $x = A^{-1}b$ 来找到 $x$。让我们尝试使用NumPy的 np.linalg.inv() 函数来求 $A$ 的逆，然后使用 @ 运算符执行矩阵乘法。

# 计算 A 的逆矩阵
try:
    A_inv = np.linalg.inv(A)
    print("A 的逆矩阵 (A_inv):\n", A_inv)

    # 计算 x = A_inv * b
    # 使用 @ 运算符进行矩阵乘法（点积）
    x_inverse_method = A_inv @ b

    print("\n解向量 x (使用逆矩阵方法):\n", x_inverse_method)

except np.linalg.LinAlgError:
    print("\n矩阵 A 是奇异的（不可逆）。无法使用逆矩阵方法求解。")

执行此代码应会给出逆矩阵 $A^{-1}$ 和得到的解向量 (vector) $x$。这种方法直接反映了数学公式 $x = A^{-1}b$。然而，计算逆矩阵可能计算开销大，且容易出现浮点不准确。

方法2：使用 np.linalg.solve 求解（推荐）

NumPy 提供了一个更直接、高效且数值更稳定的函数：np.linalg.solve(A, b)。此函数用于求解系统 $Ax = b$ 中的 $x$，而无需显式计算 $A$ 的逆矩阵。

# 使用 np.linalg.solve 直接求解系统 Ax = b
try:
    x_solve_method = np.linalg.solve(A, b)
    print("\n解向量 x (使用 np.linalg.solve):\n", x_solve_method)

except np.linalg.LinAlgError:
    print("\n矩阵 A 是奇异的（不可逆）。np.linalg.solve 无法找到唯一解。")

您会发现 x_solve_method 给出的结果与 x_inverse_method 相同（可能存在微小的浮点差异）。对于大多数实际应用，np.linalg.solve 是更优选的函数。

验证解

检查我们的解是否正确总是一个好主意。我们可以通过将计算得到的 $x$ 代回原始方程 $Ax = b$，然后查看结果是否等于 $b$ 来实现这一点。

让我们使用从 np.linalg.solve 得到的解：

# 使用从 np.linalg.solve 获得的解
x_solution = x_solve_method # Or use x_inverse_method

# 计算 A * x
b_calculated = A @ x_solution

print("\n原始向量 b:\n", b)
print("\n计算得到的 A @ x:\n", b_calculated)

# 检查计算得到的 b 是否与原始 b 接近
# 使用 np.allclose 比较浮点数
# 它检查两个数组在给定容差范围内是否按元素相等。
are_close = np.allclose(b_calculated, b)

print(f"\n计算得到的 Ax 是否接近原始 b？ {are_close}")

如果 are_close 为 True，则确认我们的解向量 (vector) $x$ 在浮点运算的小容差范围内正确满足方程组 $Ax=b$。

处理不可逆（奇异）矩阵

如果矩阵 $A$ 是奇异的（即，它没有逆矩阵）会怎样？这对应于一个方程组，该方程组要么没有解，要么有无限多个解。在这种情况下，np.linalg.inv() 和 np.linalg.solve() 都会引发 LinAlgError。

让我们尝试使用一个奇异矩阵：

# 奇异矩阵示例（第3列是第1列 + 第2列）
A_singular = np.array([
    [1, 2, 3],
    [4, 5, 9],
    [7, 8, 15]
])

# 一个虚拟的 b 向量
b_dummy = np.array([1, 1, 1]).reshape(-1, 1)

print("\n使用奇异矩阵 A_singular 进行测试:\n", A_singular)

try:
    x_singular = np.linalg.solve(A_singular, b_dummy)
    print("\n解（奇异情况）：\n", x_singular) # 这行代码不会被执行
except np.linalg.LinAlgError as e:
    print(f"\n使用奇异矩阵求解时出错: {e}")

try:
    A_singular_inv = np.linalg.inv(A_singular)
    print("\n逆矩阵（奇异情况）：\n", A_singular_inv) # 这行代码不会被执行
except np.linalg.LinAlgError as e:
    print(f"\n求奇异矩阵逆矩阵时出错: {e}")

如预期，NumPy 正确识别出矩阵是奇异的并引发错误，表明使用这些标准方法无法找到唯一解。

在本实践部分，您学习了如何处理线性方程组，使用NumPy数组 $A$ 和 $b$ 表示它，并求解未知向量 (vector) $x$。您练习了两种方法：使用显式逆矩阵（np.linalg.inv）和使用直接求解器（np.linalg.solve），并理解后者通常因效率和数值稳定性而被优选。您还学习了如何验证您的解，并看到了NumPy如何处理因奇异矩阵而导致唯一解不存在的情况。