Search code examples
pythonfloating-pointsympynumerical-methodslambdify

lambdified sympy expression returns incorrect result

I encountered following problem, first I will explain briefly what is going on:

f(x)
g(x, y) = f(x) - y

From there we expect

g(x, f(x)) = f(x) - f(x) = 0

lambdified g(x,y) returns something very close to zero instead of zero. Here's a code that reproduces the problem. It arrives only when I put sufficient amount of log evaluations in f(x)

gist: https://gist.github.com/marekyggdrasil/39a24213ebaba6293464d116821cc334

source:

from sympy import Symbol, pprint, log, lambdify

# setting symbols
g1 = Symbol("gamma1")
g2 = Symbol("gamma2")
g3 = Symbol("gamma3")
g4 = Symbol("gamma4")
rt = Symbol("rt")

# setting expressions
criteria  = (g1 * log(g1, 2.0))/2.0
criteria += (g2 * log(g2, 2.0))/2.0
criteria += (g3 * log(g3, 2.0))/2.0
criteria += (g4 * log(g4, 2.0))/2.0
rooteq = criteria - rt

print "\ncriteria function: "
pprint(criteria)

print "\ncriteria function - rt: "
pprint(rooteq)

# lambdifying expressions to callable functions
tsymbols = [g1, g2, g3, g4, rt]
lambfun_criteria = lambdify(tsymbols, criteria)
lambfun_rooteq = lambdify(tsymbols, rooteq)

# example point x
x = [0.25006462253641376, 2.2501938662000542, 2.2501938662000542, 2.2501938662000542, 0.0]

# evaluating of criteria on x
rootval = lambfun_criteria(*x)

# setting rt to this evaluation
x[4] = rootval

print "\nactual evaluation of rooteq: " + str(lambfun_rooteq(*x))
print "\nexpected evaluation of rooteq: " + str(- x[4] + lambfun_criteria(*x))

output

$ python lambdifytest.py 

criteria function: 
0.721347520444482⋅γ₁⋅log(γ₁) + 0.721347520444482⋅γ₂⋅log(γ₂) + 0.721347520444482⋅γ₃⋅log(γ₃) + 0.721347520444482⋅γ₄⋅log(γ₄)

criteria function - rt: 
0.721347520444482⋅γ₁⋅log(γ₁) + 0.721347520444482⋅γ₂⋅log(γ₂) + 0.721347520444482⋅γ₃⋅log(γ₃) + 0.721347520444482⋅γ₄⋅log(γ₄) - rt

actual evaluation of rooteq: 4.4408920985e-16

expected evaluation of rooteq: 0.0
Solution

  • I found that a custom math module can be specified when SymPy expression is lambdified. In my case mpmath helped.

    Import mpmath and set the precision to something you find sufficient

    from mpmath import mp
mp.dps = 15

    Tell SymPy you want your lambda function to use mpmath for handling floating point arithmetic

    lambfun_criteria = lambdify(tsymbols, criteria, "mpmath")
lambfun_rooteq = lambdify(tsymbols, rooteq, "mpmath")

    Then finally I get the output

    actual evaluation of rooteq: 0.0
expected evaluation of rooteq: 0.0

    The cause of the issue was the following, the default floating point arithmetic used in created lambda function was not same as one used in the script that compares the results.

    I updated the gist with full solution file if someone needs copy-paste

    https://gist.github.com/marekyggdrasil/39a24213ebaba6293464d116821cc334#file-lambdify_test_solution-py

    Edit: forgot to give docs reference

    http://docs.sympy.org/dev/modules/utilities/lambdify.html