使用 sympy 获取椭球体与平面相交产生的椭圆的特征向量

Question

给定一个具有特征向量、特征值、质心和平面的椭球，我试图找到由平面和椭球相交产生的椭圆的特征向量（如果该椭圆存在）。我现在有点被 Sympy 困住了。 采取以下示例数据：

import numpy as np
from sympy import symbols
from sympy.solvers import solve

centroid = np.array([313.81153387, 252.73655237,  78.81324428])
eigenvectors = np.array([[-0.17245704,  0.75883261,  0.62803792],
       [-0.82066049, -0.46331271,  0.33445133],
       [ 0.54477053, -0.45772742,  0.70264548]]) # where each row is one eigenvector
eigenvalues = np.array([4.03632232, 3.80325721, 7.21909427])
# Plane equation: z = 78

coords = symbols('x_prime y_prime z_prime', positive = True)
z_coor, y_coor, x_coor = np.dot(eigenvectors, coords)

# Ellipsoids have a standard form of: 
# (x - h)^2 / a^2 + (y - i)^2 / b^2 + (z - j)^2 / c^2 = 1

ellipsoid_eq_rotated = ((x_coor - centroid[0])**2 / eigenvalues[0]) + \
                       ((y_coor - centroid[1])**2 / eigenvalues[1]) + \
                       ((z_coor - centroid[2])**2 / eigenvalues[2]) - 1

intersection_eq = ellipsoid_eq_rotated.subs(z_coor, 78)
intersection_points = solve(intersection_eq, (x_coor, y_coor))

但是，这给了我一个空数组（即没有解决方案）。 当我在椭球的质心处相交时，这是没有意义的。 我哪里做错了？

Answer 1

您正在传递

solve

三个变量中的单个方程，并请求两个 表达式 的解。如有疑问，请将其打印出来。您也没有使用

intersection_eq

这个表达方式。

我还建议仅使用 SymPy，如

from sympy import symbols
from sympy.solvers import solve

centroid = [313.81153387, 252.73655237,  78.81324428]
eigenvectors = [[-0.17245704,  0.75883261,  0.62803792],
       [-0.82066049, -0.46331271,  0.33445133],
       [ 0.54477053, -0.45772742,  0.70264548]] # where each row is one eigenvector
eigenvalues = [4.03632232, 3.80325721, 7.21909427]
# Plane equation: z = 78

coords = symbols('x_prime y_prime z_prime', positive = True)
z_coor, y_coor, x_coor = Matrix(eigenvectors)*Matrix(coords)

# Ellipsoids have a standard form of: 
# (x - h)^2 / a^2 + (y - i)^2 / b^2 + (z - j)^2 / c^2 = 1

ellipsoid_eq_rotated = ((x_coor - centroid[0])**2 / eigenvalues[0]) + \
                       ((y_coor - centroid[1])**2 / eigenvalues[1]) + \
                       ((z_coor - centroid[2])**2 / eigenvalues[2]) - 1

intersection_eq = ellipsoid_eq_rotated.subs(z_coor, 78)
intersection_points = solve(ellipsoid_eq_rotated, (x_coor, y_coor))

（上面和以前有同样的问题，但现在只使用SymPy。）