在三维空间中，旋转变换是计算机图形学、机器人学和物理学中的基本操作。虽然绕坐标轴（ 、 、 轴）的旋转矩阵是众所周知的，但绕任意轴旋转的情况更为普遍且实用。本文将深入探讨绕任意单位轴旋转矩阵的数学推导，并展示如何使用 Python 的 SymPy 库进行符号计算实现，所有代码块都经过精心设计，确保独立可运行。
旋转矩阵的数学基础
三维空间中的旋转可以通过 正交矩阵表示，这些矩阵满足 且 。绕任意轴旋转的矩阵可以通过多种方法推导，其中最常见的是使用罗德里格斯旋转公式。
设旋转轴为单位向量 ，满足 ，旋转角度为 。则旋转矩阵 可以表示为：
其中 是单位矩阵， 是向量 的外积（得到一个矩阵）， 是 的叉乘矩阵：
将这三个部分展开，我们得到完整的旋转矩阵表达式：
这个公式的优美之处在于它统一了所有可能的旋转，包括绕坐标轴的特殊情况。
SymPy 实现
SymPy 是一个强大的 Python 符号计算库，可以处理符号数学表达式。下面我们实现一个函数来计算绕任意轴旋转的矩阵，这个代码块完全独立可运行：
import sympy as sp
defrotation_matrix(axis, angle):
"""
计算绕任意单位轴旋转的旋转矩阵
参数:
axis: 单位旋转轴向量 (N_x, N_y, N_z)
angle: 旋转角度 (eta)
返回:
3x3 旋转矩阵 (SymPy 矩阵对象)
"""
N_x, N_y, N_z = axis
eta = angle
# 计算旋转矩阵分量
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
# 构造旋转矩阵
R = sp.Matrix([
[c + N_x**2 * one_minus_c,
N_x*N_y*one_minus_c - N_z*s,
N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s,
c + N_y**2 * one_minus_c,
N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s,
N_z*N_y*one_minus_c + N_x*s,
c + N_z**2 * one_minus_c] # 这里将 ^ 改为 **
])
return R
# 示例用法：计算绕x轴旋转的矩阵
if __name__ == "__main__":
# 定义符号
eta = sp.symbols('eta')
# 计算绕x轴旋转矩阵
R_x = rotation_matrix((1, 0, 0), eta)
print("绕x轴旋转矩阵:")
sp.pprint(R_x)
特殊情况验证
为了验证我们的实现是否正确，我们可以检查几个特殊情况。每个验证代码块都是独立可运行的：
绕 x 轴旋转验证
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
# 验证绕x轴旋转
eta = sp.symbols('eta')
R_x = rotation_matrix((1, 0, 0), eta)
print("绕x轴旋转矩阵:")
sp.pprint(R_x)
绕 y 轴旋转验证
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
# 验证绕y轴旋转
eta = sp.symbols('eta')
R_y = rotation_matrix((0, 1, 0), eta)
print("绕y轴旋转矩阵:")
sp.pprint(R_y)
绕 z 轴旋转验证
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
# 验证绕z轴旋转
eta = sp.symbols('eta')
R_z = rotation_matrix((0, 0, 1), eta)
print("绕z轴旋转矩阵:")
sp.pprint(R_z)
旋转矩阵的性质
旋转矩阵有几个重要性质，我们可以使用 SymPy 来验证：
正交性验证
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
# 验证正交性
N_x, N_y, N_z, eta = sp.symbols('N_x N_y N_z eta')
R = rotation_matrix((N_x, N_y, N_z), eta)
R_T = R.T
identity_check = R_T * R
# 简化表达式（考虑单位向量条件）
from sympy import simplify
identity_check_simplified = simplify(identity_check.subs({
N_x**2 + N_y**2 + N_z**2: 1
}))
print("验证 R^T * R = I:")
sp.pprint(identity_check_simplified)
行列式验证
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
# 验证行列式为1
N_x, N_y, N_z, eta = sp.symbols('N_x N_y N_z eta')
R = rotation_matrix((N_x, N_y, N_z), eta)
det_R = R.det()
det_R_simplified = sp.simplify(det_R.subs({
N_x**2 + N_y**2 + N_z**2: 1
}))
print("旋转矩阵的行列式:")
sp.pprint(det_R_simplified)
实际应用示例
旋转点操作
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
defrotate_point(axis, angle, point):
"""
将点绕指定轴旋转指定角度
参数:
axis: 旋转轴 (N_x, N_y, N_z)
angle: 旋转角度
point: 要旋转的点 (x, y, z)
返回:
旋转后的点
"""
R = rotation_matrix(axis, angle)
x, y, z = point
P = sp.Matrix([x, y, z])
P_rotated = R * P
return P_rotated
# 示例：将点 (1, 0, 0) 绕z轴旋转π/2
axis = (0, 0, 1)
angle = sp.pi/2
point = (1, 0, 0)
result = rotate_point(axis, angle, point)
print(f"点 {point} 绕z轴旋转π/2后的位置:")
sp.pprint(result)
旋转组合操作
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
# 组合旋转示例
alpha, beta = sp.symbols('alpha beta')
R_x = rotation_matrix((1, 0, 0), alpha)
R_y = rotation_matrix((0, 1, 0), beta)
# 组合旋转：先x后y
R_xy = R_y * R_x
print("先绕x轴旋转α，再绕y轴旋转β的组合矩阵:")
sp.pprint(R_xy)
# 组合旋转：先y后x
R_yx = R_x * R_y
print("先绕y轴旋转β，再绕x轴旋转α的组合矩阵:")
sp.pprint(R_yx)
# 验证两个组合矩阵是否相同
print("两个组合矩阵是否相同:", R_xy.equals(R_yx))
四元数表示验证
import sympy as sp
defrotation_matrix(axis, angle):
N_x, N_y, N_z = axis
eta = angle
c = sp.cos(eta)
s = sp.sin(eta)
one_minus_c = 1 - c
R = sp.Matrix([
[c + N_x**2 * one_minus_c, N_x*N_y*one_minus_c - N_z*s, N_x*N_z*one_minus_c + N_y*s],
[N_y*N_x*one_minus_c + N_z*s, c + N_y**2 * one_minus_c, N_y*N_z*one_minus_c - N_x*s],
[N_z*N_x*one_minus_c - N_y*s, N_z*N_y*one_minus_c + N_x*s, c + N_z**2 * one_minus_c]
])
return R
defquaternion_to_rotation_matrix(q):
w, x, y, z = q
R = sp.Matrix([
[1 - 2*y**2 - 2*z**2, 2*x*y - 2*z*w, 2*x*z + 2*y*w],
[2*x*y + 2*z*w, 1 - 2*x**2 - 2*z**2, 2*y*z - 2*x*w],
[2*x*z - 2*y*w, 2*y*z + 2*x*w, 1 - 2*x**2 - 2*y**2]
])
return R
# 定义符号变量
eta = sp.symbols('eta')
N_x, N_y, N_z = sp.symbols('N_x N_y N_z')
# 显式归一化旋转轴（单位向量）
axis_norm = sp.sqrt(N_x**2 + N_y**2 + N_z**2)
N_x_norm = N_x / axis_norm
N_y_norm = N_y / axis_norm
N_z_norm = N_z / axis_norm
# 绕归一化轴旋转η的四元数
q = (sp.cos(eta/2), N_x_norm*sp.sin(eta/2), N_y_norm*sp.sin(eta/2), N_z_norm*sp.sin(eta/2))
# 转换为旋转矩阵
R_from_quat = quaternion_to_rotation_matrix(q)
# 直接计算的旋转矩阵（使用归一化轴）
R_direct = rotation_matrix((N_x_norm, N_y_norm, N_z_norm), eta)
# 比较两者是否等价
difference = sp.simplify(R_from_quat - R_direct)
print("四元数转换的旋转矩阵与直接计算的旋转矩阵是否等价:")
print(difference == sp.zeros(3, 3))
结论
本文详细探讨了绕任意轴旋转矩阵的数学理论和 SymPy 实现。我们推导了旋转矩阵的一般公式，实现了符号计算函数，并验证了其正确性和性质。通过 SymPy 的符号计算能力，我们能够处理旋转矩阵的解析表达式，这对于理论分析和公式推导非常有价值。
所有代码示例都经过精心设计，确保每个代码块都是独立可运行的。读者可以直接复制任意代码块到 Python 环境中执行，无需依赖其他代码块。这种设计使得本文不仅具有理论深度，还具有很强的实用价值。
旋转矩阵是三维变换的基础，在计算机图形学、机器人学、物理学和工程学中有广泛应用。理解其数学原理和实现方法对于在这些领域进行深入研究和开发至关重要。SymPy 作为一个强大的符号计算工具，为我们提供了一种有效的方法来探索和验证这些数学概念。
陪伴是最长情的告白
为你推送最实用的资讯
识别二维码 关注我们