【Numpy核心编程攻略：Python数据处理、分析详解与科学计算】2.22 多项式运算：从求根到拟合的数值方法

2.22 多项式运算：从求根到拟合的数值方法

目录

2.22.1 多项式表示形式

2.22.1.1 多项式的定义

多项式是一种数学表达式，由一个或多个项组成，每个项都是常数与变量的乘积。多项式可以表示为：

[ P(x) = a_n x^n + a_{n-1} x^{n-1} + \cdots + a_1 x + a_0 ]

其中，( a_n, a_{n-1}, \ldots, a_0 ) 是多项式的系数，( x ) 是变量，( n ) 是多项式的次数。

2.22.1.2 NumPy 中的多项式表示

NumPy 提供了 numpy.polynomial 模块，用于处理多项式。多项式可以使用 numpy.polynomial.Polynomial 类来表示。

import numpy as np
from numpy.polynomial import Polynomial

# 创建一个多项式 P(x) = 2x^2 + 3x + 1
p = Polynomial([1, 3, 2])  # [1, 3, 2] 表示常数项、一次项和二次项的系数
print(p)  # 输出: 1.0 + 3.0·x + 2.0·x²

2.22.1.3 多项式的基本操作

2.22.1.3.1 多项式求值

# 求多项式在 x = 2 处的值
value = p(2)  # 求 P(2)
print(value)  # 输出: 15.0

2.22.1.3.2 多项式加法

# 创建另一个多项式 Q(x) = x^2 + 2
q = Polynomial([2, 0, 1])  # [2, 0, 1] 表示常数项、一次项和二次项的系数

# 多项式加法
r = p + q  # P(x) + Q(x)
print(r)  # 输出: 3.0 + 3.0·x + 3.0·x²

2.22.1.3.3 多项式乘法

# 多项式乘法
s = p * q  # P(x) * Q(x)
print(s)  # 输出: 2.0 + 6.0·x + 5.0·x² + 3.0·x³ + 2.0·x⁴

2.22.2 特征值求根法

2.22.2.1 特征值求根法的基本原理

特征值求根法是通过将多项式问题转化为矩阵的特征值问题来求解多项式的根。具体步骤如下：

1. 构造多项式的 companion 矩阵。
2. 求解 companion 矩阵的特征值。
3. 特征值即为多项式的根。

2.22.2.2 companion 矩阵

对于多项式 ( P(x) = a_n x^n + a_{n-1} x^{n-1} + \cdots + a_1 x + a_0 )，其 companion 矩阵 ( C ) 定义为：

[ C = \begin{bmatrix}
0 & 0 & \cdots & 0 & -\frac{a_0}{a_n} \
1 & 0 & \cdots & 0 & -\frac{a_1}{a_n} \
0 & 1 & \cdots & 0 & -\frac{a_2}{a_n} \
\vdots & \vdots & \ddots & \vdots & \vdots \
0 & 0 & \cdots & 1 & -\frac{a_{n-1}}{a_n}
\end{bmatrix} ]

2.22.2.3 代码实现

import numpy as np
from numpy.linalg import eig

# 定义多项式的系数
coefficients = [1, -6, 11, -6]  # 表示 x^3 - 6x^2 + 11x - 6

# 构造 companion 矩阵
n = len(coefficients) - 1  # 多项式的次数
companion_matrix = np.zeros((n, n))

# 填充 companion 矩阵
for i in range(n-1):
    companion_matrix[i, i+1] = 1
companion_matrix[-1, :] = -np.array(coefficients[:-1]) / coefficients[-1]

# 求解特征值
eigenvalues, _ = eig(companion_matrix)
roots = eigenvalues
print(f"多项式的根: {roots}")  # 输出: 多项式的根

2.22.3 最小二乘拟合优化

2.22.3.1 最小二乘法的基本原理

最小二乘法是一种常用的拟合方法，通过最小化误差平方和来求解拟合参数。具体步骤如下：

1. 定义拟合多项式。
2. 计算拟合多项式与实际数据之间的误差平方和。
3. 通过求解最小化问题来确定拟合参数。

2.22.3.2 代码实现

2.22.3.2.1 生成数据

import numpy as np

# 生成样本数据
x = np.linspace(0, 10, 100)  # 生成 100 个在 [0, 10] 之间的等间距点
y = 2 * x**2 + 3 * x + 1 + np.random.normal(0, 10, 100)  # 生成带有噪声的多项式数据

# 绘制数据
plt.figure(figsize=(10, 6))
plt.scatter(x, y, label='Sample Data')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()

2.22.3.2.2 多项式拟合

# 使用 np.polyfit 进行多项式拟合
degrees = [1, 2, 3]  # 拟合的多项式次数
fits = {}

for degree in degrees:
    # 拟合多项式
    coefficients = np.polyfit(x, y, degree)  # 拟合多项式
    fit_polynomial = np.poly1d(coefficients)  # 创建多项式对象
    fits[degree] = fit_polynomial

    # 计算拟合值
    y_fit = fit_polynomial(x)  # 计算拟合值

    # 绘制拟合结果
    plt.figure(figsize=(10, 6))
    plt.scatter(x, y, label='Sample Data')
    plt.plot(x, y_fit, label=f'Fit {degree} Degree')
    plt.xlabel('x')
    plt.ylabel('y')
    plt.legend()
    plt.title(f'Polynomial Fit of Degree {degree}')
    plt.show()

2.22.3.2.3 拟合效果评估

# 计算拟合误差
for degree, fit_polynomial in fits.items():
    y_fit = fit_polynomial(x)
    error = np.sum((y - y_fit)**2)  # 计算误差平方和
    print(f"拟合 {degree} 次多项式的误差: {error:.2f}")  # 输出: 拟合误差

2.22.4 控制系统建模案例

2.22.4.1 控制系统简介

控制系统用于管理、指挥、指导和调节一个系统的运行，使其按期望的方式工作。多项式在控制系统中常用于描述系统的传递函数和模型。

2.22.4.2 多项式在控制系统的应用

2.22.4.2.1 传递函数

传递函数是一种描述线性时不变系统（LTI）在频域中的数学表达式。形式为：

[ H(s) = \frac{B(s)}{A(s)} ]

其中，( A(s) ) 和 ( B(s) ) 是多项式，( s ) 是复变量。

2.22.4.2.2 代码实现

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import lti, step

# 定义传递函数的系数
numerator = [1, 1]  # 分子多项式的系数
denominator = [1, 5, 6]  # 分母多项式的系数

# 创建传递函数对象
system = lti(numerator, denominator)

# 计算阶跃响应
t, y = step(system)

# 绘制阶跃响应
plt.figure(figsize=(10, 6))
plt.plot(t, y)
plt.title('Step Response of the Control System')
plt.xlabel('Time [s]')
plt.ylabel('Response')
plt.grid(True)
plt.show()

2.22.5 与SymPy对比

2.22.5.1 SymPy 简介

SymPy 是一个 Python 的符号计算库，可以用于处理符号数学问题，包括多项式的符号求解、微积分、代数方程等。

2.22.5.2 多项式的符号求解

2.22.5.2.1 SymPy 中的多项式求根

import sympy as sp

# 定义变量
x = sp.symbols('x')

# 定义多项式
polynomial = x**3 - 6*x**2 + 11*x - 6

# 求解多项式的根
roots = sp.solve(polynomial, x)
print(f"多项式的根: {roots}")  # 输出: 多项式的根

2.22.5.2.2 SymPy 与 NumPy 的对比

• NumPy：适用于数值计算，速度快，适用于大规模数据处理。
• SymPy：适用于符号计算，可以处理复杂的数学问题，但速度较慢。

2.22.5.3 多项式拟合的符号求解

# 生成带噪声的数据
x_data = np.linspace(0, 10, 100)
y_data = 2 * x_data**2 + 3 * x_data + 1 + np.random.normal(0, 10, 100)

# 转换为 SymPy 的符号变量
x_sym = sp.symbols('x')
y_sym = 2 * x_sym**2 + 3 * x_sym + 1

# 求解拟合参数
A = np.vander(x_data, 3)  # 生成 Vandermonde 矩阵
b = y_data  # 实际数据
coefficients = np.linalg.lstsq(A, b, rcond=None)[0]  # 求解最小二乘问题
fit_polynomial = np.poly1d(coefficients)  # 创建多项式对象

# 绘制拟合结果
plt.figure(figsize=(10, 6))
plt.scatter(x_data, y_data, label='Sample Data')
plt.plot(x_data, fit_polynomial(x_data), label='Fit Polynomial', color='red')
plt.title('Polynomial Fit Using SymPy and NumPy')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()

2.22.6 总结与参考文献

2.22.6.1 总结

本文详细介绍了 NumPy 中多项式的基本操作、特征值求根法、最小二乘拟合优化、控制系统建模案例以及与 SymPy 的对比。通过这些内容，读者可以深入理解如何在实际应用中高效地处理多项式运算。

2.22.6.2 参考文献

资料名称	链接
NumPy 官方文档	https://numpy.org/doc/
SciPy 官方文档	https://docs.scipy.org/doc/scipy/reference/
SymPy 官方文档	https://docs.sympy.org/latest/index.html
最小二乘法原理	https://en.wikipedia.org/wiki/Least_squares
多项式拟合原理	https://www.mathsisfun.com/data/least-squares-regression.html
控制系统建模与分析	https://www.mathworks.com/help/control/ug/modeling-control-systems.html
companion 矩阵	https://en.wikipedia.org/wiki/Companion_matrix
特征值求解	https://www.cs.cornell.edu/~bindel/class/cs4220-s12/lec16.pdf
随机数生成与拟合	https://www.stat.cmu.edu/~cshalizi/402/lectures/06-least-squares/lecture-06.pdf
Python 数据科学手册	https://www.data-science-handbook.com/
随机数生成器的高级应用	https://www.high-performance-computing.com/
多项式的数值计算	https://www.cs.princeton.edu/courses/archive/fall18/cos323/notes/COS323-polynomials.pdf
控制系统的数学原理	https://www.control.lth.se/
Python 科学计算库对比	https://www.pydata.org/

---

希望本文对您理解 NumPy 中多项式的数值方法及其应用有所帮助。这篇文章包含了详细的原理介绍、代码示例、源码注释以及案例等。希望这对您有帮助。如果有任何问题请随私信或评论告诉我。