我不能把一个线性方程象征性地解成一个变量，即使手工操作非常简单

Question

我现在很沮丧，因为我有一个简单的方程，上面有符号，不能解成变量。即使手动执行，在替换时也会出现错误（数值顺序为10^{-17}）

我无法发布mwe示例，因此我将发布问题的最简化形式

import sympy as sy
from sympy import sqrt, Derivative
r          = sy.Symbol('r')
uf         = sy.Function('uf')
Expression =-4.5*r*sqrt((r**2 - 2*r + 0.81)/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) - (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62))*(r**2 - 2*r + 0.81)*uf(r)/((-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*sqrt(r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r - (r**2 + 3.03701355989026e-33)*uf(r)**2 + 1.0*(r**2 + 0.81)**2 - 0.6561)) - 1.0*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*(r**2 - 2*r + 0.81)*(1.0*r*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*(r**2 - 2*r + 0.81)*uf(r)/((r**2 + 3.03701355989026e-33)**(3/2)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)) + 0.5*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62))*(r**2 - 2*r + 0.81)*uf(r)/((r**2 + 3.03701355989026e-33)**(3/2)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)**2) + 0.5*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*(-2*r*(r**2 - 2*r + 0.81)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*uf(r) + 4*(r - 1)*(r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*uf(r) + (r**2 + 3.03701355989026e-33)*(2*(1 - r)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 - 2*r + 0.81)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62))*uf(r) + 2*(r**2 + 3.03701355989026e-33)*(r**2 - 2*r + 0.81)*(-4*r*(r**2 + 0.81) + 1.62*r - 1.62)*uf(r) + 2*(r**2 + 3.03701355989026e-33)*(r**2 - 2*r + 0.81)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*Derivative(uf(r), r) + 2*(r*(r**2 - 2*r + 0.81) + (1 - r)*(r**2 + 3.03701355989026e-33))*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*uf(r))/((r**2 + 3.03701355989026e-33)**(3/2)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)**2))*uf(r)/(sqrt(r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)) + 0.5*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*(r**2 - 2*r + 0.81)*(-1.8*r*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*sqrt((r**2 - 2*r + 0.81)/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 + 3.03701355989026e-33)**(7/2)*(2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) - (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62))*(r**2 - 2*r + 0.81)*sqrt(1/(r**2 - 2*r + 0.81))*uf(r) - 2.0*r*(-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))**2*(r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)**2 + 1.0*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*sqrt((r**2 - 2*r + 0.81)/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 + 3.03701355989026e-33)**(7/2)*(1.8*r**2 - 5.46662440780247e-33)*(r**2 - 2*r + 0.81)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*sqrt(1/(r**2 - 2*r + 0.81))*uf(r) + 1.0*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*sqrt((r**2 - 2*r + 0.81)/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 + 3.03701355989026e-33)**(7/2)*(1.8*r*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62)) + (1.8*r**2 - 5.46662440780247e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 - 2*r + 0.81)*sqrt(1/(r**2 - 2*r + 0.81))*uf(r) + 2.0*sqrt((-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)/(r**2 - 2*r + 0.81))*(r**2 + 3.03701355989026e-33)**4*(-0.9*r*(1.8*r*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62)) + (1.8*r**2 - 5.46662440780247e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)) - (0.9*r*(1.8*r**2 - 5.46662440780247e-33) - (r**2 - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 - 2*r + 0.81)**2*sqrt(1/(r**2 - 2*r + 0.81)))/((-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))**2*(r**2 + 3.03701355989026e-33)**5*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)**2*sqrt(1/(r**2 - 2*r + 0.81))) + 2.5*sqrt((r**2 - 2*r + 0.81)/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(1.8*r**2 - 5.46662440780247e-33)*(r**2 - 2*r + 0.81)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*uf(r)/((-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*sqrt(r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r - (r**2 + 3.03701355989026e-33)*uf(r)**2 + 1.0*(r**2 + 0.81)**2 - 0.6561)) + 2.5*sqrt((r**2 - 2*r + 0.81)/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(1.8*r*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62)) + (1.8*r**2 - 5.46662440780247e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 - 2*r + 0.81)*uf(r)/((-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*sqrt(r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r - (r**2 + 3.03701355989026e-33)*uf(r)**2 + 1.0*(r**2 + 0.81)**2 - 0.6561)) - 2.5*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62))*(r**2 - 2*r + 0.81)*uf(r)**2/((r**2 + 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)*(-0.81*r**2 + 1.62*r - (r**2 + 3.03701355989026e-33)*uf(r)**2 + 1.0*(r**2 + 0.81)**2 - 0.6561)) + (-0.9*r*(r**2 + 3.03701355989026e-33)*(1.8*r*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62)) + (1.8*r**2 - 5.46662440780247e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 - 2*r + 0.81)/((-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))**2*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)) - (r**2 + 3.03701355989026e-33)*(0.9*r*(1.8*r**2 - 5.46662440780247e-33) - (r**2 - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(r**2 - 2*r + 0.81)/(-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))**2)*(r**2 - 2*r + 0.81)/(r**2 + 3.03701355989026e-33) + (-0.9*r*(1.8*r*(-2*r*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + (r**2 + 3.03701355989026e-33)*(4*r*(r**2 + 0.81) - 1.62*r + 1.62)) + (1.8*r**2 - 5.46662440780247e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561)) - (0.9*r*(1.8*r**2 - 5.46662440780247e-33) - (r**2 - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))*(1.0*(r**2 + 3.03701355989026e-33)*uf(r)**2/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + 1)*(r**2 - 2*r + 0.81)**2/((-3.24*r**2 + (-r**2 + 2*r - 3.03701355989026e-33)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))**2*(-1.0*(r**2 + 3.03701355989026e-33)*uf(r)**2/(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561) + 1.0)*(-0.81*r**2 + 1.62*r + (r**2 + 0.81)**2 - 0.6561))
sy.solve(Expression,sy.diff(uf(r),r)) 

您可以验证关于sy.diff(uf(r),r)的方程非常容易求解，因为sy.diff(Expression,sy.diff(uf(r),r),1)是非零的sy.diff(Expression,sy.diff(uf(r),r),2)是0

上面的解算器在使用jupyter时会压碎我的内核，我必须手动解算 -solution = (Expression-sy.diff(Expression,sy.diff(uf(r),r),1)*sy.diff(uf(r),r))/(-sy.diff(Expression,sy.diff(uf(r),r),1))-仍然不完全准确

Answer 1

即使这是一个笨拙的表达式，也可以使用Symphy工具解决这个问题。目标是减少必须跟踪和操作的符号数量：

0）将导数替换为z，以便于跟踪

>>> from sympy import diff
>>> from sympy.abc import z
>>> y = Expression.subs(diff(uf(r), r), z)

1）提取公共子表达式

>>> from sympy import cse
>>> reps, e = cse(y)

2）e是一个列表，获取表达式

>>> e = e[0]

3）扩展e

>>> e = e.expand(); type(e)
Add

4）隔离z

>>> i, d = e.as_independent(z); type(d)
Mul

5）求解z

>>> sol = i/d.coeff(z)

6）恢复提取的模式（按相反顺序）

>>> sol = sol.subs(reps[::-1]); sol.count_ops()
1584

综上所述，我们有一个退路：

def unwieldy_linear_solve(eq, x):
    from sympy import cse
    if not x.is_Symbol:
        d = Dummy()
        return unwieldy_linear_solve(eq.subs(x, d), d)
    r, e = cse(eq)
    i, d = e[0].expand().as_independent(x)
    assert d == x or d.is_Mul and x in d.args
    c, z = d.as_independent(x)
    assert z == x
    sol = i/c
    for o,n in r[::-1]:
        sol = sol.xreplace({o:n})
    return sol

>>> count_ops(unwieldy_linear_solve(Expression, diff(uf(r),r)))
1584