使用 Python 的 SymPy 进行符号计算
Python 的 SymPy 库为用户提供了强大的符号计算功能，使得数学表达式的处理变得简单而高效。本文将详细介绍如何使用 SymPy 进行符号计算，包括基本操作、微积分、方程求解和更多高级功能。
一、 什么是 SymPy？
SymPy 是一个用于符号数学的 Python 库，允许用户进行符号计算、代数运算、微积分、方程求解等。与数值计算不同，符号计算允许用户以符号的形式处理数学表达式，从而获得更精确的结果。
可以通过 pip 安装 SymPy：
pip install sympy
二、 基本操作
import sympy as sp# 定义符号变量
x, y = sp.symbols('x y')
# 创建表达式
expr1 = x ** 2 + 3 * x + 5
expr2 = x ** 2 - 4 * x + 4# 打印表达式
print(f"Expression 1: {expr1}") # 输出: Expression 1: x**2 + 3*x + 5
print(f"Expression 2: {expr2}") # 输出: Expression 2: x**2 - 4*x + 4# 表达式的加法
sum_expr = expr1 + expr2
print(f"Sum: {sum_expr}") # 输出: Sum: 2*x**2 - x + 9# 表达式的减法
diff_expr = expr1 - expr2
print(f"Difference: {diff_expr}") # 输出: Difference: 3*x + 1# 表达式的乘法
prod_expr = expr1 * expr2
print(f"Product: {prod_expr}") # 输出: Product: x**4 - x**3 + 7*x**2 - 12*x + 20# 表达式的除法
quot_expr = expr1 / (x + 1)
print(f"Quotient: {quot_expr}") # 输出: Quotient: (x**2 + 3*x + 5)/(x + 1)# 表达式的简化
simplified_expr = sp.simplify(prod_expr)
print(f"Simplified Product: {simplified_expr}") # 输出: Simplified Product: x**4 - x**3 + 7*x**2 - 12*x + 20# 替换变量
replaced_expr = expr1.subs(x, 2)
print(f"Replaced Expression (x=2): {replaced_expr}") # 输出: Replaced Expression (x=2): 15# 计算表达式的数值
numeric_value = expr1.evalf(subs={x: 2})
print(f"Numeric Value (x=2): {numeric_value}") # 输出: Numeric Value (x=2): 15.0000000000000# 获取 x 的系数
coefficient = expr1.coeff(x, 1) # 获取 x 的一次项系数
print(f"Coefficient of x in Expression 1: {coefficient}") # 输出: Coefficient of x in Expression 1: 3Expression 1: x**2 + 3*x + 5
Expression 2: x**2 - 4*x + 4
Sum: 2*x**2 - x + 9
Difference: 7*x + 1
Product: (x**2 - 4*x + 4)*(x**2 + 3*x + 5)
Quotient: (x**2 + 3*x + 5)/(x + 1)
Simplified Product: (x**2 - 4*x + 4)*(x**2 + 3*x + 5)
Replaced Expression (x=2): 15
Numeric Value (x=2): 15.0000000000000
Coefficient of x in Expression 1: 3
三、微积分计算
import sympy as sp # 导入 SymPy 库，并将其命名为 sp，以便后续使用# 定义符号变量
x = sp.symbols('x') # 创建一个符号变量 x，用于后续的数学表达式# 创建一个表达式
expr = x ** 2 + 3 * x + 5 # 定义一个多项式表达式：x^2 + 3x + 5# 3.1 求导
# 求导
derivative = sp.diff(expr, x) # 使用 diff 函数对 expr 进行求导，关于变量 x
print(f"Derivative: {derivative}") # 输出导数结果，应该是 2*x + 3# 3.2 积分
# 不定积分
indefinite_integral = sp.integrate(expr, x) # 使用 integrate 函数计算 expr 的不定积分，关于变量 x
print(f"Indefinite Integral: {indefinite_integral}") # 输出不定积分结果，应该是 x**3/3 + 3*x**2/2 + 5*x# 定积分
definite_integral = sp.integrate(expr, (x, 0, 1)) # 计算 expr 在区间 [0, 1] 上的定积分
print(f"Definite Integral: {definite_integral}") # 输出定积分结果，应该是 8/3# 使用 doit 计算不定积分
integral_expr = sp.Integral(expr, x) # 创建一个积分对象，表示对 expr 的不定积分
result_doit = integral_expr.doit() # 使用 doit 方法计算该积分
print(f"Indefinite Integral using doit: {result_doit}") # 输出不定积分结果，应该是 x**3/3 + 3*x**2/2 + 5*x# 3.3 极限
# 计算极限
limit_expr = sp.limit(sp.sin(x) / x, x, 0) # 计算 sin(x)/x 当 x 趋近于 0 时的极限
print(f"Limit: {limit_expr}") # 输出极限结果，应该是 1# 使用 evalf 计算数值
numeric_value = definite_integral.evalf() # 计算定积分的数值近似值
print(f"Numeric Value of Definite Integral: {numeric_value}") # 输出定积分的数值结果，应该是 2.66666666666667# 使用 evalf 计算导数在某点的数值
numeric_derivative = derivative.evalf(subs={x: 1}) # 计算导数在 x=1 时的数值
print(f"Numeric Value of Derivative at x=1: {numeric_derivative}") # 输出导数在 x=1 时的数值结果，应该是 5.00000000000000Derivative: 2*x + 3
Indefinite Integral: x**3/3 + 3*x**2/2 + 5*x
Definite Integral: 41/6
Indefinite Integral using doit: x**3/3 + 3*x**2/2 + 5*x
Limit: 1
Numeric Value of Definite Integral: 6.83333333333333
Numeric Value of Derivative at x=1: 5.00000000000000
四、矩阵计算
import sympy as sp# 1. 创建矩阵
A = sp.Matrix([[1, 2], [3, 4]]) # 创建一个 2x2 矩阵
B = sp.Matrix([[5, 6], [7, 8]]) # 创建另一个 2x2 矩阵print(f"Matrix A:\n{A}") # 输出矩阵 A
print(f"Matrix B:\n{B}") # 输出矩阵 B# 2. 矩阵运算
# 2.1 矩阵加法
C = A + B
print(f"Matrix Addition (A + B):\n{C}") # 输出矩阵 A 和 B 的和# 2.2 矩阵减法
D = A - B
print(f"Matrix Subtraction (A - B):\n{D}") # 输出矩阵 A 和 B 的差# 2.3 矩阵乘法
E = A * B
print(f"Matrix Product (A * B):\n{E}") # 输出矩阵 A 和 B 的乘积# 2.4 矩阵转置
F = A.T
print(f"Transpose of Matrix A:\n{F}") # 输出矩阵 A 的转置# 3. 矩阵的行列式和逆
# 3.1 计算行列式
det_A = A.det()
print(f"Determinant of Matrix A: {det_A}") # 输出矩阵 A 的行列式# 3.2 计算逆矩阵
if A.det() != 0: # 确保矩阵 A 是可逆的inv_A = A.inv()print(f"Inverse of Matrix A:\n{inv_A}") # 输出矩阵 A 的逆
else:print("Matrix A is not invertible.")# 4. 计算特征值和特征向量
# 4.1 计算特征值
eigenvals = A.eigenvals() # 特征值
print(f"Eigenvalues of Matrix A: {eigenvals}") # 输出特征值# 4.2 计算特征向量
eigenvects = A.eigenvects() # 特征向量
print(f"Eigenvectors of Matrix A: {eigenvects}") # 输出特征向量Matrix A:
Matrix([[1, 2], [3, 4]])
Matrix B:
Matrix([[5, 6], [7, 8]])
Matrix Addition (A + B):
Matrix([[6, 8], [10, 12]])
Matrix Subtraction (A - B):
Matrix([[-4, -4], [-4, -4]])
Matrix Product (A * B):
Matrix([[19, 22], [43, 50]])
Transpose of Matrix A:
Matrix([[1, 3], [2, 4]])
Determinant of Matrix A: -2
Inverse of Matrix A:
Matrix([[-2, 1], [3/2, -1/2]])
Eigenvalues of Matrix A: {5/2 - sqrt(33)/2: 1, 5/2 + sqrt(33)/2: 1}
Eigenvectors of Matrix A: [(5/2 - sqrt(33)/2, 1, [Matrix([
[-sqrt(33)/6 - 1/2],
[ 1]])]), (5/2 + sqrt(33)/2, 1, [Matrix([
[-1/2 + sqrt(33)/6],
[ 1]])])]
五、集合计算
import sympy as sp# 1. 创建集合
A = sp.FiniteSet(1, 2, 3) # 创建集合 A
B = sp.FiniteSet(3, 4, 5) # 创建集合 Bprint(f"Set A: {A}") # 输出集合 A
print(f"Set B: {B}") # 输出集合 B# 2. 集合运算
# 2.1 并集
union_set = A | B
print(f"Union of A and B: {union_set}") # 输出: Union of A and B: {1, 2, 3, 4, 5}# 2.2 交集
intersection_set = A & B
print(f"Intersection of A and B: {intersection_set}") # 输出: Intersection of A and B: {3}# 2.3 差集
difference_set = A - B # 使用 - 运算符计算差集
print(f"Difference of A and B (A - B): {difference_set}") # 输出: Difference of A and B (A - B): {1, 2}# 2.4 对称差集
symmetric_difference_set = A.symmetric_difference(B)
print(f"Symmetric Difference of A and B: {symmetric_difference_set}") # 输出: Symmetric Difference of A and B: {1, 2, 4, 5}# 3. 集合的基本操作
# 3.1 检查元素
print(f"Is 2 in A? {'Yes' if 2 in A else 'No'}") # 输出: Is 2 in A? Yes
print(f"Is 4 in A? {'Yes' if 4 in A else 'No'}") # 输出: Is 4 in A? No# 3.2 获取集合的大小
print(f"Size of A: {len(A)}") # 输出: Size of A: 3
print(f"Size of B: {len(B)}") # 输出: Size of B: 3Set A: {1, 2, 3}
Set B: {3, 4, 5}
Union of A and B: {1, 2, 3, 4, 5}
Intersection of A and B: {3}
Difference of A and B (A - B): {1, 2}
Symmetric Difference of A and B: {1, 2, 4, 5}
Is 2 in A? Yes
Is 4 in A? No
Size of A: 3
Size of B: 3
六、绘图
import sympy as sp
import matplotlib.pyplot as plt
import numpy as np# 1. 定义符号变量和函数
x = sp.symbols('x')# 定义多个函数
f1 = sp.sin(x) # 正弦函数
f2 = sp.cos(x) # 余弦函数
f3 = sp.tan(x) # 正切函数# 2. 生成数据点
x_vals = np.linspace(-2 * np.pi, 2 * np.pi, 400) # 从 -2π 到 2π 生成 400 个点
y_vals_f1 = [f1.evalf(subs={x: val}) for val in x_vals] # 计算正弦函数的 y 值
y_vals_f2 = [f2.evalf(subs={x: val}) for val in x_vals] # 计算余弦函数的 y 值
y_vals_f3 = [f3.evalf(subs={x: val}) for val in x_vals] # 计算正切函数的 y 值# 过滤掉无效值
y_vals_f1 = np.array(y_vals_f1, dtype=float) # 转换为 numpy 数组
y_vals_f2 = np.array(y_vals_f2, dtype=float) # 转换为 numpy 数组
y_vals_f3 = np.array(y_vals_f3, dtype=float) # 转换为 numpy 数组# 将无穷大值替换为 NaN
y_vals_f3[np.isinf(y_vals_f3)] = np.nan# 3. 绘制函数图像
plt.figure(figsize=(12, 8)) # 设置图像大小# 绘制正弦函数
plt.plot(x_vals, y_vals_f1, label='y = sin(x)', color='blue')
# 绘制余弦函数
plt.plot(x_vals, y_vals_f2, label='y = cos(x)', color='orange')
# 绘制正切函数
plt.plot(x_vals, y_vals_f3, label='y = tan(x)', color='red')# 4. 填充区域
plt.fill_between(x_vals, y_vals_f1, y_vals_f2, where=(y_vals_f1 > y_vals_f2),color='lightblue', alpha=0.5, label='Area between sin(x) and cos(x)')# 5. 添加标题和标签
plt.title('Plot of Trigonometric Functions', fontsize=16) # 图像标题
plt.xlabel('x', fontsize=14) # x 轴标签
plt.ylabel('y', fontsize=14) # y 轴标签# 6. 添加网格和图例
plt.axhline(0, color='black', linewidth=0.5, ls='--') # 添加 x 轴
plt.axvline(0, color='black', linewidth=0.5, ls='--') # 添加 y 轴
plt.grid(color='gray', linestyle='--', linewidth=0.5) # 添加网格
plt.legend() # 添加图例# 7. 设置 y 轴范围以避免正切函数的无限值
plt.ylim(-3, 3)# 8. 显示图像
plt.show()