Repository: sympy/sympy
Pull request: Integrate gives a wrong result
<!-- Your title above should be a short description of what
was changed. Do not include the issue number in the title. -->

#### References to other Issues or PRs
<!-- If this pull request fixes an issue, write "Fixes #NNNN" in that exact
format, e.g. "Fixes #1234". See
https://github.com/blog/1506-closing-issues-via-pull-requests .-->

Fixes #14782 

#### Brief description of what is fixed or changed

`integrate(sqrt(-x**2 + 1)*(-x**2 + x), [x, -1, 1])` gives `0` instead `-pi/8`.
The antiderivative works correctly, but it returns a function containing
a `Piecewise()` expression, defined in the open interval `(-1, 1)`.
The function `_eval_interval()` does not work correctly, since the antiderivative
is not a `Piecewise()` expression.
Hence, separating the terms of type `Piecewise()` and using the `_eval_interval()` for
`Piecewise()` solves the issue.

#### Other comments

Added a test.

Unfortunately, a test in `test_ode.py` fails due to the fact that a minus sign compares
inside `Piecewise()` instead that outside. So I modified such a test, without changing the
final result.



#### Release Notes

<!-- Write the release notes for this release below. See
https://github.com/sympy/sympy/wiki/Writing-Release-Notes for more information
on how to write release notes. If there is no release notes entry for this PR,
write "NO ENTRY". The bot will check your release notes automatically to see
if they are formatted correctly. -->

<!-- BEGIN RELEASE NOTES -->
* integrals
    * fixed a issue related to a wrong result in integrate
<!-- END RELEASE NOTES -->


Problem statement:
integrate(sqrt(-x**2 + 1)*(-x**2 + x), [x, -1, 1]) is incorrect
```
>>> from sympy import *
>>> x = Symbol('x')
>>> integrate(sqrt(-x**2 + 1)*(-x**2 + x), [x, -1, 1])
0
>>> 
```
The answer should be `-pi/8`.

The integral from `0` to `1` is also incorrect and returns a weird empty `Piecewise` object:

```
>>> integrate(sqrt(-x**2 + 1)*(-x**2 + x), [x, 0, 1])
-Piecewise() + 1/3
```

Hints from issue discussion:

Commit message: Added piecewise_fold after Add(*piecewises)

Patch to review:
diff --git a/sympy/integrals/integrals.py b/sympy/integrals/integrals.py
index 3119d84522ba..f380b5473851 100644
--- a/sympy/integrals/integrals.py
+++ b/sympy/integrals/integrals.py
@@ -616,18 +616,22 @@ def eval_factored(f, x, a, b):
                                     args.append(g)
                             return Mul(*args)
 
-                        integrals, others = [], []
+                        integrals, others, piecewises = [], [], []
                         for f in Add.make_args(antideriv):
                             if any(is_indef_int(g, x)
                                    for g in Mul.make_args(f)):
                                 integrals.append(f)
+                            elif any(isinstance(g, Piecewise)
+                                     for g in Mul.make_args(f)):
+                                piecewises.append(piecewise_fold(f))
                             else:
                                 others.append(f)
                         uneval = Add(*[eval_factored(f, x, a, b)
                                        for f in integrals])
                         try:
                             evalued = Add(*others)._eval_interval(x, a, b)
-                            function = uneval + evalued
+                            evalued_pw = piecewise_fold(Add(*piecewises))._eval_interval(x, a, b)
+                            function = uneval + evalued + evalued_pw
                         except NotImplementedError:
                             # This can happen if _eval_interval depends in a
                             # complicated way on limits that cannot be computed
diff --git a/sympy/integrals/tests/test_integrals.py b/sympy/integrals/tests/test_integrals.py
index 4832988be5b3..4da5558bb09f 100644
--- a/sympy/integrals/tests/test_integrals.py
+++ b/sympy/integrals/tests/test_integrals.py
@@ -1353,3 +1353,8 @@ def test_issue_14877():
     f = exp(1 - exp(x**2)*x + 2*x**2)*(2*x**3 + x)/(1 - exp(x**2)*x)**2
     assert integrate(f, x) == \
         -exp(2*x**2 - x*exp(x**2) + 1)/(x*exp(3*x**2) - exp(2*x**2))
+
+def test_issue_14782():
+    f = sqrt(-x**2 + 1)*(-x**2 + x)
+    assert integrate(f, [x, -1, 1]) == - pi / 8
+    assert integrate(f, [x, 0, 1]) == S(1) / 3 - pi / 16
diff --git a/sympy/solvers/tests/test_ode.py b/sympy/solvers/tests/test_ode.py
index 4783e7eaf976..cee7c0e0fa73 100644
--- a/sympy/solvers/tests/test_ode.py
+++ b/sympy/solvers/tests/test_ode.py
@@ -1584,7 +1584,7 @@ def test_1st_homogeneous_coeff_ode3():
     sol = Eq(log(f(x)), C1 - Piecewise(
             (-acosh(f(x)/x), abs(f(x)**2)/x**2 > 1),
             (I*asin(f(x)/x), True)))
-    assert dsolve(eq, hint='1st_homogeneous_coeff_subs_indep_div_dep') == sol
+    assert simplify(dsolve(eq, hint='1st_homogeneous_coeff_subs_indep_div_dep') - sol) == 0
 
 
 def test_1st_homogeneous_coeff_corner_case():


Code review comments:
diff_hunk: @@ -1584,7 +1584,7 @@ def test_1st_homogeneous_coeff_ode3():
     sol = Eq(log(f(x)), C1 - Piecewise(
             (-acosh(f(x)/x), abs(f(x)**2)/x**2 > 1),
             (I*asin(f(x)/x), True)))
-    assert dsolve(eq, hint='1st_homogeneous_coeff_subs_indep_div_dep') == sol
+    assert simplify(dsolve(eq, hint='1st_homogeneous_coeff_subs_indep_div_dep') - sol) == 0
path: sympy/solvers/tests/test_ode.py
text: @user1:
The use of `simplify` is discouraged as it is an expensive operation. In this case, it is better to edit `sol`.