Thanks Aaron. The main assumption that we can make for our package is that these variables are real numbers, so that will still leave us with divide by zero issues. Maybe once I have a list of expressions that are common to our problems, we can figure out how to deal with that subset.
Jason moorepants.info +01 530-601-9791 On Wed, Sep 11, 2013 at 1:45 PM, Aaron Meurer <[email protected]> wrote: > With your example, it matters what order you substitute the variables. > If you substitute vx=0 first, you get 0/sqrt(vy**2 + vz**2) == 0. If > you substitute vy=vz=0, you get vx/sqrt(vx**2), which should simplify > to sqrt(vx) for vx >= 0 (which it is, since vx=0). Can you assume > nonnegative for these variables. Of course, if you substitute all > three at once, you get 0/0. > > There's no general solution to this problem, but you aren't working in > a completely general situation. I think the solution you should adopt > will depend on just how much structure you know about your > expressions. > > Sometimes you can rewrite expressions so that they don't have 0/0 > (like you learn in calculus). I don't see how to do it for this > expression, though. > > Aaron Meurer > > On Tue, Sep 10, 2013 at 2:43 PM, Jason Moore <[email protected]> wrote: > > Also, why does .evalf give zero in my previous example? Shouldn't that > also > > raise a DivideByZero error? > > > > > > Jason > > moorepants.info > > +01 530-601-9791 > > > > > > On Tue, Sep 10, 2013 at 4:41 PM, Jason Moore <[email protected]> > wrote: > >> > >> I'll don't have a concrete example to show at the moment, but I feel > like > >> we get these singularity issues when we do use simplification and don't > when > >> we don't or when we expand everything out. I'll try to put together an > >> example where this is the case. > >> > >> > >> Jason > >> moorepants.info > >> +01 530-601-9791 > >> > >> > >> On Tue, Sep 10, 2013 at 4:14 PM, Stefan Krastanov > >> <[email protected]> wrote: > >>> > >>> > In this case, real and positive assumptions don't help for the > lambdify > >>> > call. I guess making lambdify assumption aware could be an option. > >>> > >>> I was actually expecting that one would always use `simplify` before > >>> `lambdify`. This would make a big difference, because with the > >>> appropriate assumptions `simplify` would remove all such removable > >>> singularities[1], which is the technical term for these point. > >>> Unfortunately a quick google search for "numerical removable > >>> singularity" does not show anything promising. > >>> > >>> [1]: https://en.wikipedia.org/wiki/Removable_singularity > >>> > >>> On 10 September 2013 14:19, Jason Moore <[email protected]> wrote: > >>> > In this case, real and positive assumptions don't help for the > lambdify > >>> > call. I guess making lambdify assumption aware could be an option. > >>> > > >>> > I think we should define dynamicsymbols in mechanics to have a > default > >>> > real > >>> > assumption, this could help in simplification. > >>> > > >>> > The different basis is an interesting idea. Do you have an example of > >>> > that? > >>> > > >>> > > >>> > Jason > >>> > moorepants.info > >>> > +01 530-601-9791 > >>> > > >>> > > >>> > On Tue, Sep 10, 2013 at 2:01 PM, Stefan Krastanov > >>> > <[email protected]> wrote: > >>> >> > >>> >> The poorman's solution that I used when I had this problem was to do > >>> >> the calculations in a different basis. I do not know how > >>> >> general/automatic this solution would be. > >>> >> > >>> >> However, a common reason for bad/no simplification is lack of > >>> >> appropriate assumptions. Define all you symbols as real/positive and > >>> >> half of the issues will go away. > >>> >> > >>> >> On 10 September 2013 12:17, Jason Moore <[email protected]> > wrote: > >>> >> > In the mechanics package we sometimes end up with symbolic > >>> >> > expressions > >>> >> > that > >>> >> > can give a divide by zero error if we use the expression in a > >>> >> > numerical > >>> >> > function. Below is a simple example of an expression that could be > >>> >> > generated > >>> >> > by the mechanics package. It should evaluate to zero (i.e the > limit > >>> >> > as > >>> >> > vx -> > >>> >> > zero. evalf says it's zero but if I use lambdify (or probably any > >>> >> > code > >>> >> > generator) to transform the expression to something that is purely > >>> >> > numeric > >>> >> > we get divide by zeros. I feel like this must be a common issue > when > >>> >> > transforming symbolic expressions to numerical functions. Does > >>> >> > anyone > >>> >> > know > >>> >> > how or if this is dealt with? It be nice if we could generate > >>> >> > equations > >>> >> > of > >>> >> > motion in the mechanics package and generate code from them that > >>> >> > wouldn't > >>> >> > have divide by zero errors, requiring some fix on the numerical > >>> >> > side. > >>> >> > > >>> >> > In [1]: from sympy import symbols, sqrt, lambdify > >>> >> > > >>> >> > In [2]: k, vx, vy, vz = symbols('k vx vy vz') > >>> >> > > >>> >> > In [3]: f = k * sqrt(vx ** 2 + vy ** 2 + vz ** 2) > >>> >> > > >>> >> > In [4]: dfdvx = f.diff(vx) > >>> >> > > >>> >> > In [5]: dfdvx > >>> >> > Out[5]: k*vx/sqrt(vx**2 + vy**2 + vz**2) > >>> >> > > >>> >> > In [6]: dfdvx.evalf(subs={k: 1.0, vx: 0.0, vy: 0.0, vz: 0.0}) > >>> >> > Out[6]: 0 > >>> >> > > >>> >> > In [7]: numeric_dfdvx = lambdify((k, vx, vy, vz), dfdvx) > >>> >> > > >>> >> > In [8]: numeric_dfdvx(1.0, 0.0, 0.0, 0.0) > >>> >> > > >>> >> > > >>> >> > > --------------------------------------------------------------------------- > >>> >> > ZeroDivisionError Traceback (most recent > >>> >> > call > >>> >> > last) > >>> >> > <ipython-input-8-55e44f8ab1bf> in <module>() > >>> >> > ----> 1 numeric_dfdvx(1.0, 0.0, 0.0, 0.0) > >>> >> > > >>> >> > /usr/local/lib/python2.7/dist-packages/numpy/__init__.pyc in > >>> >> > <lambda>(_Dummy_26, _Dummy_27, _Dummy_28, _Dummy_29) > >>> >> > > >>> >> > ZeroDivisionError: float division by zero > >>> >> > > >>> >> > In [9]: numeric_dfdvx = lambdify((k, vx, vy, vz), dfdvx, > >>> >> > modules='numpy') > >>> >> > > >>> >> > In [10]: numeric_dfdvx(1.0, 0.0, 0.0, 0.0) > >>> >> > /usr/local/lib/python2.7/dist-packages/numpy/__init__.py:1: > >>> >> > RuntimeWarning: > >>> >> > invalid value encountered in double_scalars > >>> >> > """ > >>> >> > Out[10]: nan > >>> >> > > >>> >> > In [11]: numeric_dfdvx = lambdify((k, vx, vy, vz), dfdvx, > >>> >> > modules='mpmath') > >>> >> > > >>> >> > In [12]: numeric_dfdvx(1.0, 0.0, 0.0, 0.0) > >>> >> > > >>> >> > > >>> >> > > --------------------------------------------------------------------------- > >>> >> > ZeroDivisionError Traceback (most recent > >>> >> > call > >>> >> > last) > >>> >> > <ipython-input-12-55e44f8ab1bf> in <module>() > >>> >> > ----> 1 numeric_dfdvx(1.0, 0.0, 0.0, 0.0) > >>> >> > > >>> >> > <string> in <lambda>(k, vx, vy, vz) > >>> >> > > >>> >> > > >>> >> > > /usr/local/lib/python2.7/dist-packages/sympy/mpmath/ctx_mp_python.pyc in > >>> >> > __rdiv__(s, t) > >>> >> > 206 if t is NotImplemented: > >>> >> > 207 return t > >>> >> > --> 208 return t / s > >>> >> > 209 > >>> >> > 210 def __rpow__(s, t): > >>> >> > > >>> >> > > >>> >> > > /usr/local/lib/python2.7/dist-packages/sympy/mpmath/ctx_mp_python.pyc in > >>> >> > __div__(self, other) > >>> >> > > >>> >> > > /usr/local/lib/python2.7/dist-packages/sympy/mpmath/libmp/libmpf.pyc > >>> >> > in > >>> >> > mpf_div(s, t, prec, rnd) > >>> >> > 928 if not sman or not tman: > >>> >> > 929 if s == fzero: > >>> >> > --> 930 if t == fzero: raise ZeroDivisionError > >>> >> > 931 if t == fnan: return fnan > >>> >> > 932 return fzero > >>> >> > > >>> >> > ZeroDivisionError: > >>> >> > > >>> >> > In [19]: dfdvx.limit(vx, 0.0) > >>> >> > Out[19]: 0 > >>> >> > > >>> >> > > >>> >> > Jason > >>> >> > moorepants.info > >>> >> > +01 530-601-9791 > >>> >> > > >>> >> > -- > >>> >> > You received this message because you are subscribed to the Google > >>> >> > Groups > >>> >> > "sympy" group. > >>> >> > To unsubscribe from this group and stop receiving emails from it, > >>> >> > send > >>> >> > an > >>> >> > email to [email protected]. > >>> >> > To post to this group, send email to [email protected]. > >>> >> > Visit this group at http://groups.google.com/group/sympy. > >>> >> > For more options, visit https://groups.google.com/groups/opt_out. > >>> >> > >>> >> -- > >>> >> You received this message because you are subscribed to the Google > >>> >> Groups > >>> >> "sympy" group. > >>> >> To unsubscribe from this group and stop receiving emails from it, > send > >>> >> an > >>> >> email to [email protected]. > >>> >> To post to this group, send email to [email protected]. > >>> >> Visit this group at http://groups.google.com/group/sympy. > >>> >> For more options, visit https://groups.google.com/groups/opt_out. > >>> > > >>> > > >>> > -- > >>> > You received this message because you are subscribed to the Google > >>> > Groups > >>> > "sympy" group. > >>> > To unsubscribe from this group and stop receiving emails from it, > send > >>> > an > >>> > email to [email protected]. > >>> > To post to this group, send email to [email protected]. > >>> > Visit this group at http://groups.google.com/group/sympy. > >>> > For more options, visit https://groups.google.com/groups/opt_out. > >>> > >>> -- > >>> You received this message because you are subscribed to the Google > Groups > >>> "sympy" group. > >>> To unsubscribe from this group and stop receiving emails from it, send > an > >>> email to [email protected]. > >>> To post to this group, send email to [email protected]. > >>> Visit this group at http://groups.google.com/group/sympy. > >>> For more options, visit https://groups.google.com/groups/opt_out. > >> > >> > > > > -- > > You received this message because you are subscribed to the Google Groups > > "sympy" group. > > To unsubscribe from this group and stop receiving emails from it, send an > > email to [email protected]. > > To post to this group, send email to [email protected]. > > Visit this group at http://groups.google.com/group/sympy. > > For more options, visit https://groups.google.com/groups/opt_out. > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To unsubscribe from this group and stop receiving emails from it, send an > email to [email protected]. > To post to this group, send email to [email protected]. > Visit this group at http://groups.google.com/group/sympy. > For more options, visit https://groups.google.com/groups/opt_out. > -- You received this message because you are subscribed to the Google Groups "sympy" group. To unsubscribe from this group and stop receiving emails from it, send an email to [email protected]. To post to this group, send email to [email protected]. Visit this group at http://groups.google.com/group/sympy. For more options, visit https://groups.google.com/groups/opt_out.
