当前位置： 文江博客话题详情

等待一个小时后，方程式无法解决

发布于 2025-02-11 18:25:14 字数 773 浏览 2 评论 0原文

我试图求解一个由两个有两个未知数的方程组成的非线性系统，我尝试使用Sympy模块来完成任务。 在尝试运行下面的代码（On Atom和Google Consoil Notebook）之后，什么也不会发生，并且该过程永远不会终止（等待15分钟）。 有什么问题？我的系统要复杂吗？ 我尝试使用Scipy的FSOLVE，并且能够立即获得根源

这是下面的代码

import sympy as sym

sym.init_printing()
x,y = sym.symbols('x,y')
dmin = 70 #minimum of the actuator
dmax = 140 #maximum of the actuator

l3 = 10
l4 = 10

Amin = 0.523599 #min angle in radians
Amax = 2.0944 #max angle in radians

cos1 = sym.cos(Amin - sym.atan(l4/x) - sym.atan(l3/y))
cos2 = sym.cos(Amax - sym.atan(l4/x) - sym.atan(l3/y))

a2 = x**2 + l4**2
b2 = y**2 + l3**2

f=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos1, dmin**2)
g=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos2, dmax**2)

print(sym.solve([f,g],(x,y)))
原文

I was trying to solve a non-linear system that consists of two equations with 2 unknowns, I tried using Sympy module to complete the task.
After trying to run the code below (on Atom and Google Collab notebook), nothing is happening and the process is never terminated (waited for over that 15min).
What could be the problem? is my system to complicated?
I tried to use fsolve from scipy and i was able to get one of the roots in no time

This is the code below

import sympy as sym

sym.init_printing()
x,y = sym.symbols('x,y')
dmin = 70 #minimum of the actuator
dmax = 140 #maximum of the actuator

l3 = 10
l4 = 10

Amin = 0.523599 #min angle in radians
Amax = 2.0944 #max angle in radians

cos1 = sym.cos(Amin - sym.atan(l4/x) - sym.atan(l3/y))
cos2 = sym.cos(Amax - sym.atan(l4/x) - sym.atan(l3/y))

a2 = x**2 + l4**2
b2 = y**2 + l3**2

f=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos1, dmin**2)
g=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos2, dmax**2)

print(sym.solve([f,g],(x,y)))

如果你对这篇内容有疑问，欢迎到本站社区发帖提问 参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论， 你可以免费 注册 一个本站的账号。

评论（2）

呢古 2025-02-18 18:25:15

您的方程式类似：

In [160]: f
Out[160]: 
               __________    __________                                                 
 2    2       ╱  2          ╱  2           ⎛    ⎛10⎞       ⎛10⎞           ⎞             
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅cos⎜atan⎜──⎟ + atan⎜──⎟ - 0.523599⎟ + 200 = 4900
                                           ⎝    ⎝x ⎠       ⎝y ⎠           ⎠             

In [161]: g
Out[161]: 
               __________    __________                                                
 2    2       ╱  2          ╱  2           ⎛    ⎛10⎞       ⎛10⎞         ⎞              
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅cos⎜atan⎜──⎟ + atan⎜──⎟ - 2.0944⎟ + 200 = 19600
                                           ⎝    ⎝x ⎠       ⎝y ⎠         ⎠

通常是代数和超验功能的混合物，没有任何隔离的解决方案，但是在这种情况下，cos可以使用atan取消，以使纯代代代数纯化方程式：

In [162]: expand_trig(f)
Out[162]: 
               __________    __________                                                                                                                                                     
 2    2       ╱  2          ╱  2        ⎛      0.866025291583566                 5.00000194337561                  5.00000194337561                   86.6025291583566        ⎞             
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅⎜───────────────────────────── + ─────────────────────────────── + ─────────────────────────────── - ─────────────────────────────────⎟ + 200 = 4900
                                        ⎜     _________      _________          _________      _________          _________      _________            _________      _________⎟             
                                        ⎜    ╱     100      ╱     100          ╱     100      ╱     100          ╱     100      ╱     100            ╱     100      ╱     100 ⎟             
                                        ⎜   ╱  1 + ─── ⋅   ╱  1 + ───    y⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅y⋅   ╱  1 + ─── ⋅   ╱  1 + ─── ⎟             
                                        ⎜  ╱         2    ╱         2        ╱         2    ╱         2        ╱         2    ╱         2          ╱         2    ╱         2 ⎟             
                                        ⎝╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y         ╲╱         x   ╲╱         y  ⎠             

In [163]: expand_trig(g)
Out[163]: 
               __________    __________                                                                                                                                                        
 2    2       ╱  2          ╱  2        ⎛        0.500004241445914                 8.6602295497065                   8.6602295497065                    50.0004241445914        ⎞              
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅⎜- ───────────────────────────── + ─────────────────────────────── + ─────────────────────────────── + ─────────────────────────────────⎟ + 200 = 19600
                                        ⎜       _________      _________          _________      _________          _________      _________            _________      _________⎟              
                                        ⎜      ╱     100      ╱     100          ╱     100      ╱     100          ╱     100      ╱     100            ╱     100      ╱     100 ⎟              
                                        ⎜     ╱  1 + ─── ⋅   ╱  1 + ───    y⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅y⋅   ╱  1 + ─── ⋅   ╱  1 + ─── ⎟              
                                        ⎜    ╱         2    ╱         2        ╱         2    ╱         2        ╱         2    ╱         2          ╱         2    ╱         2 ⎟              
                                        ⎝  ╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y         ╲╱         x   ╲╱         y  ⎠

使用此功能，我们可以将系统按摩到一对多项式中，然后用groebner bases求解：

In [181]: from sympy.solvers.solvers import unrad

In [174]: eq1 = trigsimp(unrad(expand_trig(nsimplify(f)).rewrite(Add))[0]/(x*y)).factor().args[-1]

In [175]: eq2 = trigsimp(unrad(expand_trig(nsimplify(g)).rewrite(Add))[0]/(x*y)).factor().args[-1]

In [176]: gb = groebner([nsimplify(eq1.n()), nsimplify(eq2.n())], [x, y])

In [177]: for p in gb: print(p.n(2))
x - 2.9e-27*y**15 - 4.4e-25*y**14 + 2.8e-22*y**13 + 4.2e-20*y**12 - 1.2e-17*y**11 - 1.6e-15*y**10 + 2.7e-13*y**9 + 3.2e-11*y**8 - 3.6e-9*y**7 - 3.6e-7*y**6 + 2.7e-5*y**5 + 0.0021*y**4 - 0.099*y**3 - 6.2*y**2 + 1.3e+2*y + 6.5e+3
y**16 + 1.4e+2*y**15 - 1.0e+5*y**14 - 1.3e+7*y**13 + 4.4e+9*y**12 + 4.7e+11*y**11 - 1.1e+14*y**10 - 9.0e+15*y**9 + 1.5e+18*y**8 + 9.3e+19*y**7 - 1.3e+22*y**6 - 5.1e+23*y**5 + 6.3e+25*y**4 + 1.3e+27*y**3 - 1.5e+29*y**2 - 1.2e+30*y + 1.3e+32

In [178]: ry = [r.n() for r in real_roots(gb[-1])]

In [179]: ry
Out[179]: 
[-112.210054906433, -94.6849042632114, -53.6834203472101, -49.4635109350851, 48.579017974275, 50.194835
2176397, 82.5120279645666, 103.812504924146, 118.173390238157, 142.710285016351]

In [180]: for yval in ry:
     ...:     print('y =', yval, ', x =', solve(gb[0].subs(y, yval), x)[0])
     ...: 
y = -112.210054906433 , x = 48.5790179742034
y = -94.6849042632114 , x = -53.6834203471953
y = -53.6834203472101 , x = -94.6849042632093
y = -49.4635109350851 , x = 103.812504924146
y = 48.5790179742750 , x = -112.210054906433
y = 50.1948352176397 , x = 118.173390238151
y = 82.5120279645666 , x = 142.710285016423
y = 103.812504924146 , x = -49.4635109348164
y = 118.173390238157 , x = 50.1948352170293
y = 142.710285016351 , x = 82.5120279639959

重写为多项式引入虚假解决方案，因此只有其中一些实际上是f和g g g的解决方案。

如果您对所需的特定解决方案有初始猜测，则更简单，更快地计算这些解决方案的方法只是使用nsolde来数值求解系统：

In [194]: nsolve([f, g], [x, y], [118, 50])
Out[194]: 
⎡118.173390238157⎤
⎢                ⎥
⎣50.1948352176396⎦

Your equations are like:

In [160]: f
Out[160]: 
               __________    __________                                                 
 2    2       ╱  2          ╱  2           ⎛    ⎛10⎞       ⎛10⎞           ⎞             
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅cos⎜atan⎜──⎟ + atan⎜──⎟ - 0.523599⎟ + 200 = 4900
                                           ⎝    ⎝x ⎠       ⎝y ⎠           ⎠             

In [161]: g
Out[161]: 
               __________    __________                                                
 2    2       ╱  2          ╱  2           ⎛    ⎛10⎞       ⎛10⎞         ⎞              
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅cos⎜atan⎜──⎟ + atan⎜──⎟ - 2.0944⎟ + 200 = 19600
                                           ⎝    ⎝x ⎠       ⎝y ⎠         ⎠

Usually a mix of algebraic and transcendental functions like this would not have any anlytic solution but in this case the cos can cancel with the atan to make purely algebraic equations:

In [162]: expand_trig(f)
Out[162]: 
               __________    __________                                                                                                                                                     
 2    2       ╱  2          ╱  2        ⎛      0.866025291583566                 5.00000194337561                  5.00000194337561                   86.6025291583566        ⎞             
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅⎜───────────────────────────── + ─────────────────────────────── + ─────────────────────────────── - ─────────────────────────────────⎟ + 200 = 4900
                                        ⎜     _________      _________          _________      _________          _________      _________            _________      _________⎟             
                                        ⎜    ╱     100      ╱     100          ╱     100      ╱     100          ╱     100      ╱     100            ╱     100      ╱     100 ⎟             
                                        ⎜   ╱  1 + ─── ⋅   ╱  1 + ───    y⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅y⋅   ╱  1 + ─── ⋅   ╱  1 + ─── ⎟             
                                        ⎜  ╱         2    ╱         2        ╱         2    ╱         2        ╱         2    ╱         2          ╱         2    ╱         2 ⎟             
                                        ⎝╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y         ╲╱         x   ╲╱         y  ⎠             

In [163]: expand_trig(g)
Out[163]: 
               __________    __________                                                                                                                                                        
 2    2       ╱  2          ╱  2        ⎛        0.500004241445914                 8.6602295497065                   8.6602295497065                    50.0004241445914        ⎞              
x  + y  - 2⋅╲╱  x  + 100 ⋅╲╱  y  + 100 ⋅⎜- ───────────────────────────── + ─────────────────────────────── + ─────────────────────────────── + ─────────────────────────────────⎟ + 200 = 19600
                                        ⎜       _________      _________          _________      _________          _________      _________            _________      _________⎟              
                                        ⎜      ╱     100      ╱     100          ╱     100      ╱     100          ╱     100      ╱     100            ╱     100      ╱     100 ⎟              
                                        ⎜     ╱  1 + ─── ⋅   ╱  1 + ───    y⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅   ╱  1 + ─── ⋅   ╱  1 + ───    x⋅y⋅   ╱  1 + ─── ⋅   ╱  1 + ─── ⎟              
                                        ⎜    ╱         2    ╱         2        ╱         2    ╱         2        ╱         2    ╱         2          ╱         2    ╱         2 ⎟              
                                        ⎝  ╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y       ╲╱         x   ╲╱         y         ╲╱         x   ╲╱         y  ⎠

Using this we can massage the system into a pair of polynomials and then solve it with Groebner bases:

In [181]: from sympy.solvers.solvers import unrad

In [174]: eq1 = trigsimp(unrad(expand_trig(nsimplify(f)).rewrite(Add))[0]/(x*y)).factor().args[-1]

In [175]: eq2 = trigsimp(unrad(expand_trig(nsimplify(g)).rewrite(Add))[0]/(x*y)).factor().args[-1]

In [176]: gb = groebner([nsimplify(eq1.n()), nsimplify(eq2.n())], [x, y])

In [177]: for p in gb: print(p.n(2))
x - 2.9e-27*y**15 - 4.4e-25*y**14 + 2.8e-22*y**13 + 4.2e-20*y**12 - 1.2e-17*y**11 - 1.6e-15*y**10 + 2.7e-13*y**9 + 3.2e-11*y**8 - 3.6e-9*y**7 - 3.6e-7*y**6 + 2.7e-5*y**5 + 0.0021*y**4 - 0.099*y**3 - 6.2*y**2 + 1.3e+2*y + 6.5e+3
y**16 + 1.4e+2*y**15 - 1.0e+5*y**14 - 1.3e+7*y**13 + 4.4e+9*y**12 + 4.7e+11*y**11 - 1.1e+14*y**10 - 9.0e+15*y**9 + 1.5e+18*y**8 + 9.3e+19*y**7 - 1.3e+22*y**6 - 5.1e+23*y**5 + 6.3e+25*y**4 + 1.3e+27*y**3 - 1.5e+29*y**2 - 1.2e+30*y + 1.3e+32

In [178]: ry = [r.n() for r in real_roots(gb[-1])]

In [179]: ry
Out[179]: 
[-112.210054906433, -94.6849042632114, -53.6834203472101, -49.4635109350851, 48.579017974275, 50.194835
2176397, 82.5120279645666, 103.812504924146, 118.173390238157, 142.710285016351]

In [180]: for yval in ry:
     ...:     print('y =', yval, ', x =', solve(gb[0].subs(y, yval), x)[0])
     ...: 
y = -112.210054906433 , x = 48.5790179742034
y = -94.6849042632114 , x = -53.6834203471953
y = -53.6834203472101 , x = -94.6849042632093
y = -49.4635109350851 , x = 103.812504924146
y = 48.5790179742750 , x = -112.210054906433
y = 50.1948352176397 , x = 118.173390238151
y = 82.5120279645666 , x = 142.710285016423
y = 103.812504924146 , x = -49.4635109348164
y = 118.173390238157 , x = 50.1948352170293
y = 142.710285016351 , x = 82.5120279639959

Rewriting as polynomial introduces spurious solutions so only some of those are actually solutions of f and g.

A simpler and faster way to compute these solutions if you have an initial guess for a particular solution that you want is just to use nsolve to solve the system numerically:

In [194]: nsolve([f, g], [x, y], [118, 50])
Out[194]: 
⎡118.173390238157⎤
⎢                ⎥
⎣50.1948352176396⎦
回复 原文
不忘初心 2025-02-18 18:25:15

这是延续方法的另一个好候选人。讨论是在链接上的，但我将使用c进行以下操作，以控制方程的更非线性术语：

import sympy as sym
x,y = sym.symbols('x,y')
dmin = 70 #minimum of the actuator
dmax = 140 #maximum of the actuator

l3 = 10
l4 = 10

Amin = 0.523599 #min angle in radians
Amax = 2.0944 #max angle in radians
for c in range(10):
    c=c/S(10)
    cos1 = sym.cos(Amin - sym.atan(l4/x)*c - sym.atan(l3/y)*c)
    cos2 = sym.cos(Amax - sym.atan(l4/x)*c - sym.atan(l3/y)*c)
    a2 = x**2 + l4**2
    b2 = y**2 + l3**2
    f=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos1, dmin**2)
    g=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos2, dmax**2)
    if not c:
        gs = solve((f,g),(x,y),dict=False)
    else:
        sol = []
        for gi in gs:
            sol.append(nsolve((f,g),(x,y),gi))

gs

[(-107.967591836877, -48.6044884662980), 
(-107.967591836877, 48.6044884662980),
(-48.6044884662980, -107.967591836877),
(-48.6044884662980, 107.967591836877),
(48.6044884662980, -107.967591836877),
(48.6044884662980, 107.967591836877),
(107.967591836877, -48.6044884662980),
(107.967591836877, 48.6044884662980)]

' c转到1成为最终解决方案

>>> for i in ordered(sol):tuple(i)
...
(-111.790176548536, 48.5400429865187)
(-96.2635727943765, -52.7612279529407)
(-52.7612279529407, -96.2635727943765)
(-49.3470005072566, 104.216202951134)
(48.5400429865187, -111.790176548536)
(49.8069789979713, 117.241028072945)
(104.216202951134, -49.3470005072566)
(117.241028072945, 49.8069789979713)

This is another good candidate for the continuation method. The discussion is at the link but I would apply it as follows using c to control the more non-linear term of your equations:

import sympy as sym
x,y = sym.symbols('x,y')
dmin = 70 #minimum of the actuator
dmax = 140 #maximum of the actuator

l3 = 10
l4 = 10

Amin = 0.523599 #min angle in radians
Amax = 2.0944 #max angle in radians
for c in range(10):
    c=c/S(10)
    cos1 = sym.cos(Amin - sym.atan(l4/x)*c - sym.atan(l3/y)*c)
    cos2 = sym.cos(Amax - sym.atan(l4/x)*c - sym.atan(l3/y)*c)
    a2 = x**2 + l4**2
    b2 = y**2 + l3**2
    f=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos1, dmin**2)
    g=sym.Eq(a2 + b2 - 2*sym.sqrt(a2)*sym.sqrt(b2)*cos2, dmax**2)
    if not c:
        gs = solve((f,g),(x,y),dict=False)
    else:
        sol = []
        for gi in gs:
            sol.append(nsolve((f,g),(x,y),gi))

gs starts out as

[(-107.967591836877, -48.6044884662980), 
(-107.967591836877, 48.6044884662980),
(-48.6044884662980, -107.967591836877),
(-48.6044884662980, 107.967591836877),
(48.6044884662980, -107.967591836877),
(48.6044884662980, 107.967591836877),
(107.967591836877, -48.6044884662980),
(107.967591836877, 48.6044884662980)]

and is gradually modified as c goes to 1 to become the final solution

>>> for i in ordered(sol):tuple(i)
...
(-111.790176548536, 48.5400429865187)
(-96.2635727943765, -52.7612279529407)
(-52.7612279529407, -96.2635727943765)
(-49.3470005072566, 104.216202951134)
(48.5400429865187, -111.790176548536)
(49.8069789979713, 117.241028072945)
(104.216202951134, -49.3470005072566)
(117.241028072945, 49.8069789979713)
回复 原文
~没有更多了~

关于作者

温暖的光

暂无简介

文章
评论
27 人气
发私信

相关话题

更多

热门标签

操作系统 程序设计 IT运维 Linux系统管理 JavaScript 服务器应用 solaris C/C++ PHP Shell BSD Vue.js aix Oracle Python HTML 系统管理 HTML5 CSS 前端
更多

推荐作者

櫻之舞

文章 0 评论 0

弥枳

文章 0 评论 0

m2429

文章 0 评论 0

寻找一个思念的角度

文章 0 评论 0

野却迷人

文章 0 评论 0

我怀念的。

文章 0 评论 0

更多

友情链接

文江博客
    我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站，即表示您同意使用 Cookies 和您的相关数据。
    原文