如何使用 SymPy 求解方程组以找到一些角度？
How to solve equation sets with SymPy to find some angles?
我不太明白如何在 Python 中用 SymPy 求解我的方程组。我试过使用 solve() 和 nsolve()。 solve() 没有给我答案，它只是站着走，nsolve() 给我这个错误。
ValueError: Could not find root within given tolerance. (7378.8100001153525709 > 2.16840434497100886801e-19) Try another starting point or tweak arguments.
这是我尝试使用solve()时写的：
from sympy import *
theta_1, theta_2, theta_3 = symbols('theta_1 theta_2 theta_3')
eq1 = -136.2*sin(theta_2)*sin(theta_3)*cos(theta_1) + 136.2*cos(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*cos(theta_1)*cos(theta_2)
eq2 = -136.2*sin(theta_1)*sin(theta_2)*cos(theta_3) + 136.2*sin(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*sin(theta_1)*cos(theta_2)
eq3 = 136.2*sin(theta_2)*cos(theta_3)+ 222.1*sin(theta_2) + 136.2*sin(theta_3)*cos(theta_2) + 100.9
eq1 = Eq(eq1, 0)
eq2 = Eq(eq2, (-323.90))
eq3 = Eq(eq3, 176.71)
ans = solve((eq1, eq2, eq3), (theta_1, theta_2, theta_3))
print(ans)
这是我尝试使用nsolve()时写的：
from sympy import *
theta_1, theta_2, theta_3 = symbols('theta_1 theta_2 theta_3')
eq1 = -136.2*sin(theta_2)*sin(theta_3)*cos(theta_1) + 136.2*cos(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*cos(theta_1)*cos(theta_2)
eq2 = -136.2*sin(theta_1)*sin(theta_2)*cos(theta_3) + 136.2*sin(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*sin(theta_1)*cos(theta_2)
eq3 = 136.2*sin(theta_2)*cos(theta_3)+ 222.1*sin(theta_2) + 136.2*sin(theta_3)*cos(theta_2) + 100.9
ans = nsolve((eq1, eq2, eq3), (theta_1, theta_2, theta_3), (0, -323.9, 176.71))
print(ans)
我想找到 theta1、theta2 和 theta3 的表达式。最好是 arctan2 的形式。有谁知道我该怎么做？或者为什么这不起作用？
我假设您的输入应该是以度为单位的角度。这里 nsolve 可以给你一个答案，如果你以弧度为单位输入：
In [48]: p = (pi/180).evalf()
In [49]: nsolve((eq1, eq2, eq3), (theta_1, theta_2, theta_3), (0, -323.9*p, 176.71*p))
Out[49]:
⎡-1.5707963267949 ⎤
⎢ ⎥
⎢-9.36640725122338⎥
⎢ ⎥
⎣-2.49008907692194⎦
我们可以将您的系统简化为具有如下有理系数的多项式系统：
In [60]: c1, c2, c3, s1, s2, s3 = symbols('c1:4 s1:4')
...: rep = {cos(theta_1): c1, cos(theta_2): c2, cos(theta_3): c3,
...: sin(theta_1): s1, sin(theta_2): s2, sin(theta_3): s3}
...:
...: eqs = [eq.subs(rep) for eq in [eq1, eq2, eq3]] + [
...: c1**2 + s1**2 - 1, c2**2 + s2**2 - 1, c3**2 + s3**2 - 1]
...: eqs = [nsimplify(eq) for eq in eqs]
In [61]: for eq in eqs:
...: print(eq)
...:
Eq(681*c1*c2*c3/5 - 681*c1*s2*s3/5 + 2221*c1/10, 0)
Eq(681*c2*c3*s1/5 - 681*c3*s1*s2/5 + 2221*s1/10, -3239/10)
Eq(681*c2*s3/5 + 681*c3*s2/5 + 2221*s2/10 + 1009/10, 17671/100)
c1**2 + s1**2 - 1
c2**2 + s2**2 - 1
c3**2 + s3**2 - 1
从那里您可以使用 SymPy 的 groebner 函数来计算多项式系统解的 Groebner 基。为了对称，我们将添加一个分隔元素 t:
In [62]: gb = groebner(eqs, order='grevlex').fglm('lex')
这就像一个新的分离方程列表，其解是您想要的不同 theta 的 cos 和 sin 值。
最终方程是 c3 的 20 次多项式，其解可以代入其他方程得到系统的完整解。多项式有 4 个实根，每个实根都是一个较小的 8 次不可约多项式的根。由于 Abel-Ruffini 根不能有根表达式，但可以用数值计算它们：
In [70]: c3num = [r.n() for r in real_roots(gb[-1])]
In [71]: c3num
Out[71]: [-0.795172954862362, -0.552080421425289, 0.702791508766261, 0.99923774387104]
其中的每一个都可以用来解决其余的 Groebner 基础：
In [73]: for c3n in c3num:
...: for eq in gb[:-1]:
...: print(solve([eq.subs(c3, c3n)], dict=True))
...: print(solve([c3 - c3n], dict=True))
...: print()
...:
[{s1: -1.00000000000000}]
[{s2: -0.0583375690202956}]
[{s3: -0.606382694211788}]
[{c1: 0}]
[{c1: 0.0}]
[{c2: -0.998296913792728}]
[{c3: -0.795172954862362}]
[{s1: -1.00000000000000}]
[{s2: 0.881316344165583}]
[{s3: 0.833790866034178}]
[{c1: 0}]
[{c1: 0.0}]
[{c2: -0.472526720395763}]
[{c3: -0.552080421425289}]
[{s1: -1.00000000000000}]
[{s2: -0.0656752871541357}]
[{s3: 0.711395878031908}]
[{c1: 0}]
[{c1: 0.0}]
[{c2: 0.997841047763359}]
[{c3: 0.702791508766261}]
[{s1: -1.00000000000000}]
[{s2: 0.226102985521720}]
[{s3: -0.0390375619754195}]
[{c1: -3.05175781250000e-5}, {c1: 3.05175781250000e-5}]
[{c1: 0.0}]
[{c2: 0.974103403161280}]
[{c3: 0.999237743871040}]
最后一种情况似乎在数值上不一致，但这可能是舍入误差，因为这在数值上非常不稳定，也许解应该是 c1=0。我认为发生这种情况是因为您的方程式稍微退化了。
上面给出了角度的 cos 和 sin 的所有解决方案，您可以使用 atan2 从那里得到角度本身。您可以使用 nsolve 根据需要精确地优化上述任何一项。 4 个案例中的第一个是上面 nsolve 找到的。
如果您想在此处使用显式解析表达式，那似乎不太可能。
我不太明白如何在 Python 中用 SymPy 求解我的方程组。我试过使用 solve() 和 nsolve()。 solve() 没有给我答案，它只是站着走，nsolve() 给我这个错误。
ValueError: Could not find root within given tolerance. (7378.8100001153525709 > 2.16840434497100886801e-19) Try another starting point or tweak arguments.
这是我尝试使用solve()时写的：
from sympy import *
theta_1, theta_2, theta_3 = symbols('theta_1 theta_2 theta_3')
eq1 = -136.2*sin(theta_2)*sin(theta_3)*cos(theta_1) + 136.2*cos(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*cos(theta_1)*cos(theta_2)
eq2 = -136.2*sin(theta_1)*sin(theta_2)*cos(theta_3) + 136.2*sin(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*sin(theta_1)*cos(theta_2)
eq3 = 136.2*sin(theta_2)*cos(theta_3)+ 222.1*sin(theta_2) + 136.2*sin(theta_3)*cos(theta_2) + 100.9
eq1 = Eq(eq1, 0)
eq2 = Eq(eq2, (-323.90))
eq3 = Eq(eq3, 176.71)
ans = solve((eq1, eq2, eq3), (theta_1, theta_2, theta_3))
print(ans)
这是我尝试使用nsolve()时写的：
from sympy import *
theta_1, theta_2, theta_3 = symbols('theta_1 theta_2 theta_3')
eq1 = -136.2*sin(theta_2)*sin(theta_3)*cos(theta_1) + 136.2*cos(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*cos(theta_1)*cos(theta_2)
eq2 = -136.2*sin(theta_1)*sin(theta_2)*cos(theta_3) + 136.2*sin(theta_1)*cos(theta_2)*cos(theta_3) + 222.1*sin(theta_1)*cos(theta_2)
eq3 = 136.2*sin(theta_2)*cos(theta_3)+ 222.1*sin(theta_2) + 136.2*sin(theta_3)*cos(theta_2) + 100.9
ans = nsolve((eq1, eq2, eq3), (theta_1, theta_2, theta_3), (0, -323.9, 176.71))
print(ans)
我想找到 theta1、theta2 和 theta3 的表达式。最好是 arctan2 的形式。有谁知道我该怎么做？或者为什么这不起作用？
我假设您的输入应该是以度为单位的角度。这里 nsolve 可以给你一个答案，如果你以弧度为单位输入：
In [48]: p = (pi/180).evalf()
In [49]: nsolve((eq1, eq2, eq3), (theta_1, theta_2, theta_3), (0, -323.9*p, 176.71*p))
Out[49]:
⎡-1.5707963267949 ⎤
⎢ ⎥
⎢-9.36640725122338⎥
⎢ ⎥
⎣-2.49008907692194⎦
我们可以将您的系统简化为具有如下有理系数的多项式系统：
In [60]: c1, c2, c3, s1, s2, s3 = symbols('c1:4 s1:4')
...: rep = {cos(theta_1): c1, cos(theta_2): c2, cos(theta_3): c3,
...: sin(theta_1): s1, sin(theta_2): s2, sin(theta_3): s3}
...:
...: eqs = [eq.subs(rep) for eq in [eq1, eq2, eq3]] + [
...: c1**2 + s1**2 - 1, c2**2 + s2**2 - 1, c3**2 + s3**2 - 1]
...: eqs = [nsimplify(eq) for eq in eqs]
In [61]: for eq in eqs:
...: print(eq)
...:
Eq(681*c1*c2*c3/5 - 681*c1*s2*s3/5 + 2221*c1/10, 0)
Eq(681*c2*c3*s1/5 - 681*c3*s1*s2/5 + 2221*s1/10, -3239/10)
Eq(681*c2*s3/5 + 681*c3*s2/5 + 2221*s2/10 + 1009/10, 17671/100)
c1**2 + s1**2 - 1
c2**2 + s2**2 - 1
c3**2 + s3**2 - 1
从那里您可以使用 SymPy 的 groebner 函数来计算多项式系统解的 Groebner 基。为了对称，我们将添加一个分隔元素 t:
In [62]: gb = groebner(eqs, order='grevlex').fglm('lex')
这就像一个新的分离方程列表，其解是您想要的不同 theta 的 cos 和 sin 值。
最终方程是 c3 的 20 次多项式，其解可以代入其他方程得到系统的完整解。多项式有 4 个实根，每个实根都是一个较小的 8 次不可约多项式的根。由于 Abel-Ruffini 根不能有根表达式，但可以用数值计算它们：
In [70]: c3num = [r.n() for r in real_roots(gb[-1])]
In [71]: c3num
Out[71]: [-0.795172954862362, -0.552080421425289, 0.702791508766261, 0.99923774387104]
其中的每一个都可以用来解决其余的 Groebner 基础：
In [73]: for c3n in c3num:
...: for eq in gb[:-1]:
...: print(solve([eq.subs(c3, c3n)], dict=True))
...: print(solve([c3 - c3n], dict=True))
...: print()
...:
[{s1: -1.00000000000000}]
[{s2: -0.0583375690202956}]
[{s3: -0.606382694211788}]
[{c1: 0}]
[{c1: 0.0}]
[{c2: -0.998296913792728}]
[{c3: -0.795172954862362}]
[{s1: -1.00000000000000}]
[{s2: 0.881316344165583}]
[{s3: 0.833790866034178}]
[{c1: 0}]
[{c1: 0.0}]
[{c2: -0.472526720395763}]
[{c3: -0.552080421425289}]
[{s1: -1.00000000000000}]
[{s2: -0.0656752871541357}]
[{s3: 0.711395878031908}]
[{c1: 0}]
[{c1: 0.0}]
[{c2: 0.997841047763359}]
[{c3: 0.702791508766261}]
[{s1: -1.00000000000000}]
[{s2: 0.226102985521720}]
[{s3: -0.0390375619754195}]
[{c1: -3.05175781250000e-5}, {c1: 3.05175781250000e-5}]
[{c1: 0.0}]
[{c2: 0.974103403161280}]
[{c3: 0.999237743871040}]
最后一种情况似乎在数值上不一致，但这可能是舍入误差，因为这在数值上非常不稳定，也许解应该是 c1=0。我认为发生这种情况是因为您的方程式稍微退化了。
上面给出了角度的 cos 和 sin 的所有解决方案，您可以使用 atan2 从那里得到角度本身。您可以使用 nsolve 根据需要精确地优化上述任何一项。 4 个案例中的第一个是上面 nsolve 找到的。
如果您想在此处使用显式解析表达式，那似乎不太可能。