How to substitute a subexpression with a function?

36 views
Skip to first unread message

Audrius-St

unread,
Jan 24, 2022, 4:28:50 PM1/24/22
to sympy
Hello, In the following test code I am attempting to simplify the resulting expression
by replacing a subexpression with a function

import symengine as se
import sympy as sp

# Potential V
def V(x, y):
    v = (x*x + y*y)/2 + y*(x*x - (y*y)/3)
    return v

def main():

    # Symbols
    x, y, Ho = sp.symbols('x, y, Ho')

    expr = sp.Pow(Ho - V(x, y), -1)
    dexpr_dx = se.Derivative(expr, x, 1).xreplace({y*(x**2 + (-1/3)*y**2) + (-1/2)*x**2 + (-1/2)*y**2: V(x, y)})

    print("dexpr_dx =", dexpr_dx)

if __name__ == "__main__":
    main()

Taking in derivative w.r.t. x resulting in V(x, y) being replace the expression
dexpr_dx = -(-x - 2*x*y)/(Ho - y*(x**2 + (-1/3)*y**2) + (-1/2)*x**2 + (-1/2)*y**2)**2

which I would like to rewrite as -(-x - 2*x*y)/(Ho - V(x, y))**2 as I will be repeating the differentiation in the actual code.
In other words, I would like to reduce the output code bloat.

I've tried subs(), replace(), and xreplace() with no success. Any insight in how to do this would be appreciated.







Chris Smith

unread,
Jan 25, 2022, 3:53:27 PM1/25/22
to sympy
There are a few reasons this is not working.

1) if you have floats like (1/3) those will not match the rationals like Rational(1, 3)
2) xreplace will target full nodes, not sub-nodes. So when trying to match an Add that is a subexpression, it will fail.
3) you have a sign error in the key of the replacement dictionary: the first term should be `-y*...`
4) the replacement value, `V(x, y)` recreates the thing you are replacing.

Try the following:

>>> d = Derivative(expr, x, 1)
>>> d.subs(V(x,y), Symbol('V(x,y)'))
Derivative(1/(Ho - V(x,y)), x)


/c

Hanspeter Schmid

unread,
Jan 26, 2022, 5:52:21 AM1/26/22
to sympy
Regarding 1), if you write 1/3, then Python converts this into a float before it is handed over to sympy. To prevent it, you could import S (Singleton) from sympy and write S(1)/3.

dexpr_dx = se.Derivative(expr, x, 1).xreplace({y*(x**2 + (-S(1)/3)*y**2) + (-S(1)/2)*x**2 + (-S(1)/2)*y**2: V(x, y)})

Audrius-St

unread,
Jan 26, 2022, 1:20:39 PM1/26/22
to sympy

After a bit of experimentation, the following code fragment worked:

    dexpr_dx = se.Derivative(expr, x, 1).subs(
                    (Ho - y*(x**2 + sp.Rational(-1, 3)*y**2) +
                     sp.Rational(-1, 2)*x**2 +
                     sp.Rational(-1, 2)*y**2), sp.Symbol('(Ho - V(x, y))'))

yielding

dexpr_dx -(-x - 2*x*y)/(Ho - V(x, y))**2
Thank you both for your advice.
Reply all
Reply to author
Forward
0 new messages