数学建模1：解析方法与几何模型

这两天时间有限，简略看了第一章节内容，记录了其中坐标选择与线性代数相关的知识。之前在学习机器人与矩阵论的时候学习过这些内容，再次做了复习

坐标旋转

坐标旋转是指在二维或三维平面内，将一个点围绕原点旋转一定的角度后得到的新坐标。

二维旋转矩阵

二维平面内的旋转矩阵：$\left[\begin{array}{cc}\cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{array}\right]$

其中，$\theta$ 是旋转的角度。

设原坐标为 $(x,y)$，旋转后的新坐标为 $(x^{\prime },y^{\prime })$，则有：

import numpy as np

def rotate_point_2d(x, y, theta):
     """
    将点 (x, y) 绕原点旋转 theta 角度。

    参数：
    x (float): 原始x坐标
    y (float): 原始y坐标
    theta (float): 旋转角度，单位为弧度

    返回：
    (float, float): 旋转后的新坐标 (x', y')
    """
     # 计算旋转矩阵
    rotation_matrix = np.array([
        [np.cos(theta), -np.sin(theta)],
        [np.sin(theta), np.cos(theta)]
    ])
    
    # 原始坐标
    original_point = np.array([x, y])
    
    # 旋转后的新坐标
    new_point = np.dot(rotation_matrix, original_point)
    
    return new_point[0], new_point[1]

# 示例用法
x, y = 1.0, 0.0
theta = np.pi / 4  # 45度
x_new, y_new = rotate_point_2d(x, y, theta)
print(f"原始坐标: ({x}, {y})")
print(f"旋转 {theta} 弧度后的新坐标: ({x_new}, {y_new})")

三维旋转矩阵

三维空间内的旋转可以绕X轴、Y轴或Z轴进行，每个轴的旋转矩阵如下：

• 绕X轴旋转：$R_{\mathrm{x}}(\gamma )=\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \gamma & -\sin \gamma \\ 0& \sin \gamma & \cos \gamma \end{array}\right]$

• 绕Y轴旋转：$R_{\mathrm{y}}(\beta )=\left[\begin{array}{ccc}\cos \beta & 0& \sin \beta \\ 0& 1& 0\\ -\sin \beta & 0& \sin \beta \end{array}\right]$

• 绕Z轴旋转：$R_{\mathrm{z}}(\alpha )=\left[\begin{array}{ccc}\cos \alpha & -\sin \alpha & 0\\ sin\alpha & \cos \alpha & 0\\ 0& 0& 1\end{array}\right]$

• 旋转矩阵：$R=R_{\mathrm{z}}(\alpha )R_{\mathrm{y}}(\beta )R_{\mathrm{x}}(\gamma )$

设原坐标为 $(x,y,z)$，分别绕X轴、Y轴、Z轴旋转后的新坐标为 $(x^{\prime },y^{\prime },z^{\prime })$，公式如下：

• 绕X轴旋转：$$\{\begin{array}{l}{x}^{\prime }=x\\ y^{\prime }=y\cos \theta -z\sin \theta \\ z^{\prime }=y\sin \theta +z\cos \theta \end{array}$$

• 绕Y轴旋转：$$\{\begin{array}{l}{x}^{\prime }=x\cos \theta +z\sin \theta \\ y^{\prime }=y\\ z^{\prime }=-x\sin \theta +z\cos \theta \end{array}$$

• 绕Z轴旋转：$$\{\begin{array}{l}{x}^{\prime }=x\cos \theta -y\sin \theta \\ y^{\prime }=x\sin \theta +y\cos \theta \\ z^{\prime }=z\end{array}$$

import numpy as np

def rotate_point_3d(x, y, z, theta_x, theta_y, theta_z):
	"""
    将点 (x, y,z) 绕x轴旋转 theta_x 角度，绕y轴旋转 theta_y 角度，绕z轴旋转 theta_z 角度。

    参数：
    x (float): 原始x坐标
    y (float): 原始y坐标
    z (float): 原始z坐标
    theta_x (float): 绕x轴旋转角度，单位为弧度
	theta_y (float): 绕y轴旋转角度，单位为弧度
	theta_z (float): 绕z轴旋转角度，单位为弧度
	
    返回：
    (float, float, float): 旋转后的新坐标 (x', y', z')
    """
    # 绕X轴旋转矩阵
    rotation_matrix_x = np.array([
        [1, 0, 0],
        [0, np.cos(theta_x), -np.sin(theta_x)],
        [0, np.sin(theta_x), np.cos(theta_x)]
    ])
    # 绕Y轴旋转矩阵
    rotation_matrix_y = np.array([
        [np.cos(theta_y), 0, np.sin(theta_y)],
        [0, 1, 0],
        [-np.sin(theta_y), 0, np.cos(theta_y)]
    ])
    # 绕Z轴旋转矩阵
    rotation_matrix_z = np.array([
        [np.cos(theta_z), -np.sin(theta_z), 0],
        [np.sin(theta_z), np.cos(theta_z), 0],
        [0, 0, 1]
    ])
    
    # 先绕X轴旋转
    point = np.dot(rotation_matrix_x, np.array([x, y, z]))
    # 再绕Y轴旋转
    point = np.dot(rotation_matrix_y, point)
    # 最后绕Z轴旋转
    point = np.dot(rotation_matrix_z, point)
    
    return point[0], point[1], point[2]

# 示例用法
x, y, z = 1.0, 0.0, 0.0
theta_x = np.pi / 4  # 绕X轴旋转45度
theta_y = np.pi / 4  # 绕Y轴旋转45度
theta_z = np.pi / 4  # 绕Z轴旋转45度
x_new, y_new, z_new = rotate_point_3d(x, y, z, theta_x, theta_y, theta_z)
print(f"原始坐标: ({x}, {y}, {z})")
print(f"分别绕X、Y、Z轴旋转 {theta_x}, {theta_y}, {theta_z} 弧度后的新坐标: ({x_new}, {y_new}, {z_new})")

线性代数

numpy.linalg 模块提供了线性代数的相关函数，包括矩阵的分解、求解线性方程组、计算行列式等。

1. 矩阵乘法

矩阵乘法是线性代数中最基本的操作之一。numpy.dot 或 @ 运算符可以用于矩阵乘法

import numpy as np

# 定义两个矩阵
A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])

# 矩阵乘法
C = np.dot(A, B)
# 或者使用 @ 运算符
C = A @ B

print("矩阵 A:")
print(A)
print("矩阵 B:")
print(B)
print("矩阵 A * B:")
print(C)

2. 矩阵的逆

使用 numpy.linalg.inv 计算矩阵的逆。如果矩阵不可逆，将抛出 LinAlgError

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4]])

# 计算矩阵的逆
try:
    A_inv = np.linalg.inv(A)
    print("矩阵 A 的逆:")
    print(A_inv)
except np.linalg.LinAlgError:
    print("矩阵 A 不可逆")

3. 线性方程组求解

使用 numpy.linalg.solve 求解线性方程组 Ax = b

import numpy as np

# 定义矩阵 A 和向量 b
A = np.array([[1, 2], [3, 4]])
b = np.array([5, 6])

# 求解线性方程组 Ax = b
x = np.linalg.solve(A, b)
print("线性方程组的解 x:")
print(x)

4. 矩阵的特征值和特征向量

使用 numpy.linalg.eig 计算矩阵的特征值和特征向量

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4]])

# 计算特征值和特征向量
eigenvalues, eigenvectors = np.linalg.eig(A)

print("矩阵 A 的特征值:")
print(eigenvalues)
print("矩阵 A 的特征向量:")
print(eigenvectors)

5. 矩阵的行列式

使用 numpy.linalg.det 计算矩阵的行列式

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4]])

# 计算行列式
det_A = np.linalg.det(A)
print("矩阵 A 的行列式:")
print(det_A)

6. 奇异值分解 (SVD)

使用 numpy.linalg.svd 进行奇异值分解

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4], [5, 6]])

# 进行奇异值分解
U, S, V = np.linalg.svd(A)

print("矩阵 A 的左奇异向量矩阵 U:")
print(U)
print("矩阵 A 的奇异值 S:")
print(S)
print("矩阵 A 的右奇异向量矩阵 V:")
print(V)

7. 二范数

使用 numpy.linalg.norm 计算向量或矩阵的二范数

import numpy as np

# 定义一个向量和矩阵
v = np.array([1, 2, 3])
A = np.array([[1, 2], [3, 4]])

# 计算向量的二范数
v_norm = np.linalg.norm(v)
print("向量 v 的二范数:")
print(v_norm)

# 计算矩阵的二范数
A_norm = np.linalg.norm(A)
print("矩阵 A 的二范数:")
print(A_norm)

8. QR分解

使用 numpy.linalg.qr 进行QR分解，将矩阵分解为正交矩阵Q和上三角矩阵R

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4]])

# 进行QR分解
Q, R = np.linalg.qr(A)

print("矩阵 A 的正交矩阵 Q:")
print(Q)
print("矩阵 A 的上三角矩阵 R:")
print(R)

9. Cholesky分解

使用 numpy.linalg.cholesky 进行Cholesky分解，将正定矩阵分解为下三角矩阵和其转置的乘积

import numpy as np

# 定义一个正定矩阵
A = np.array([[4, 12], [12, 37]])

# 进行Cholesky分解
L = np.linalg.cholesky(A)

print("矩阵 A 的Cholesky分解矩阵 L:")
print(L)
print("L 的转置矩阵 L.T:")
print(L.T)

10. Moore-Penrose 伪逆

使用 numpy.linalg.pinv 计算矩阵的Moore-Penrose伪逆

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4], [5, 6]])

# 计算矩阵的Moore-Penrose伪逆
A_pinv = np.linalg.pinv(A)

print("矩阵 A 的Moore-Penrose伪逆:")
print(A_pinv)

11. 矩阵的条件数

使用 numpy.linalg.cond 计算矩阵的条件数

import numpy as np

# 定义一个矩阵
A = np.array([[1, 2], [3, 4]])

# 计算矩阵的条件数
cond_A = np.linalg.cond(A)
print("矩阵 A 的条件数:")
print(cond_A)

12. 矩阵的秩

使用 numpy.linalg.matrix_rank 计算矩阵的秩

import numpy as np

# 生成数据
np.random.seed(0)
X = np.random.rand(100, 2)
y = 3 * X[:, 0] + 5 * X[:, 1] + np.random.randn(100)

# 添加截距项
X = np.hstack([np.ones((X.shape[0], 1)), X])

# 求解最小二乘解
beta, _, _, _ = np.linalg.lstsq(X, y, rcond=None)
print("线性回归的系数:")
print(beta)