Search code examples
python-2.xsympysubstitutionsymbolic-math

sympy: Incorrect substitution results


I am unable to evaluate and plot a simple (known) function in sympy within Ipython notebook:

y(x) = ( F / 6EI )( x^3 - 3Lx^2 ) where:

  • F = 10^6
  • E = 200E9
  • I = (1/12)(0.5*1^3)
  • L = 3

I define the expression using symbols() and substitute known values (using subs()).

import sympy
from sympy import symbols

sympy.init_printing()
from IPython.display import Math, Image
from IPython.display import display


F, E, L, I, x = symbols('F, E, L, I, x')    #Definition of sympy symbols

y = F/(6*E*I) * (x**3 - 3*L*x**2)    #Define beam deflection equation
display(Math("y=" + latex(y)))

y = y.subs({F:10**6, E:200E9, L:3, I:(1/12)*(0.5)*(1**3)})    #Substitute known values to define specific deflection curve
display(Math("y=" + latex(y)))

However, it results in:

y(x) = (infinity)x^3 -(infinity)x^2

This is clearly incorrect, as the coefficients should be rational numbers - this is the well understood cantilevered beam deflection equation. It should evaluate to:

y(x) = 2E-5*x^3 - 1.8E-4*x^2 where:

  • y(x=0) = 0
  • y(x=L) = FL^3 / 3EI.

Why does sympy produce this result and how can I modify my solution to attain the correct solution?


As noted in the comments, the above code works correctly for Python 3 but fails in Python 2.7 (the version I am running).


Solution

  • If you are using Python 2, 1/12 results in 0, due to integer division. You can use from __future__ import division, or use Python 3 (I recommend using Python 3), and it will return a float.

    Better, though, with SymPy, is to use a rational. If you use Rational(1, 12) you will get an exact result. In general, your answers with SymPy will be more accurate if you use exact numbers (rational numbers), and then use evalf() to convert the expression to a float at the end.