Weird results from sympy.integrate and Heaviside

171 views Asked by At

I'm getting weird results from integrate in sympy using Heaviside. Here's my code:

import sympy as sp
x = sp.symbols('x')
L=1
def q(x): return (x-1/L)*sp.Heaviside(x-1/L)
def V(foo): return sp.integrate(q(x),(x,foo))
display(V(x))

If I run this, it works. But if I change L to, say, 2, it crashes. Why would this matter?

1

There are 1 answers

0
Oscar Benjamin On BEST ANSWER

This is a bug but there is a way to work around it which is to explicitly rewrite as Piecewise:

In [39]: expr
Out[39]: (x - 0.5)⋅θ(x - 0.5)

In [40]: expr.rewrite(Piecewise)
Out[40]: 
          ⎛⎧ 0    for x - 0.5 < 0⎞
          ⎜⎪                     ⎟
(x - 0.5)⋅⎜⎨θ(0)  for x - 0.5 = 0⎟
          ⎜⎪                     ⎟
          ⎝⎩ 1    for x - 0.5 > 0⎠

In [41]: expr.rewrite(Piecewise).integrate(x)
Out[41]: 
⎧          0             for x < 0.5
⎪                                   
⎨     2                             
⎪0.5⋅x  - 0.5⋅x + 0.125   otherwise 
⎩  

Disabling the meijerg algorithm also gives a correct result:

In [44]: integrate(expr, meijerg=False)
Out[44]: 
⎛ 2                ⎞           
⎜x                 ⎟           
⎜── - 0.5⋅x + 0.125⎟⋅θ(x - 0.5)
⎝2                 ⎠