Integration Error

[David M. Rogers]

Hello,

  I was double-checking an integral with sympy, and noticed that the software comes up with the wrong answer.  The integrand has five terms, which can be collected into three groups.  Two of the groups -- Ec1 and Ec2 -- cancel after integrating over z.  Strangely, sympy gets the answer right if you ask it to integrate the groups separately.

x, z, R, a = symbols('x z R a')
r = sqrt(x**2 + z**2)
u = erf(a*r/sqrt(2))/r
Ec = diff(u, z, z).subs([(x, sqrt(R*R-z*z))])

# compare
simplify(integrate(Ec, (z, -R, R)))

# to the analytical answer
-2*sqrt(2)*R*a**3*exp(-R**2*a**2/2)/(3*sqrt(pi))

# in particular, Ec = Ec1 + Ec2 + Ec3, where
Ec1 = sympify('3*z**2*erf(sqrt(2)*a*sqrt(R**2)/2)/(R**2)**(5/2) - erf(sqrt(2)*a*sqrt(R**2)/2)/(R**2)**(3/2)')
Ec2 = sympify('+ sqrt(2)*a*exp(-R**2*a**2/2)/(sqrt(pi)*R**2) - 3*sqrt(2)*a*z**2*exp(-R**2*a**2/2)/(sqrt(pi)*R**4)')

Ec3 = Ec - Ec1 - Ec2 # -sqrt(2)*a**3*z**2*exp(-R**2*a**2/2)/(sqrt(pi)*R**2)

simplify(integrate(Ec1, (z, -R, R))) # -> 0
simplify(integrate(Ec2, (z, -R, R))) # -> 0

integrate(Ec3, (z, -R, R))


~ David M. Rogers
  University of South Florida

[Julien Rioux]

On Sunday, 9 March 2014 00:32:09 UTC+1, David M. Rogers wrote:

Hello,

  I was double-checking an integral with sympy, and noticed that the software comes up with the wrong answer.  The integrand has five terms, which can be collected into three groups.  Two of the groups -- Ec1 and Ec2 -- cancel after integrating over z.  Strangely, sympy gets the answer right if you ask it to integrate the groups separately.

Could you come up with simpler expressions reproducing the problem? That will help a lot in determining what goes wrong and where.

Thanks,
Cheers,
Julien