Bug introduced in 9.0 or earlier and persisting through 11.0.1 or later
I have looked at previous questions and I'm aware that this seems to be a known bug: Mathematica giving inconsistent results for symbolic integrals done in different ways. The origins for the bugs seem to be varied.
I'm having a hard time tracking down where Mathematica is going wrong with this particular one, especially since I am explicitly assuming that all my variables ($L$ and $\rho$) are positive. I've played around with the assumptions without any luck. Here's my code:
$Assumptions = {ρ > 0, L > 0};
limits = Sequence[{ux, 0, L}, {vx, 0, L}, {uy, 0, ux}, {vy, 0, vx}]
N1 = Integrate[a1 = Exp[-ρ (ux - uy) (vx - vy)], limits];
N2 = Integrate[a2 = Exp[-ρ (ux - uy) (vx - vy)] (ρ (ux - uy) (vx - vy)), limits];
N3 = Integrate[a3 = Exp[-ρ (ux - uy) (vx - vy)] (ρ (ux - uy) (vx - vy))^2, limits];
N4 = Integrate[a4 = Exp[-ρ (ux - uy) (vx - vy)] (ρ (ux - uy) (vx - vy))^3, limits];
int1 = Integrate[4 ρ (a1 - 7/2 a2 + 2 a3 - 1/4 a4), limits]
int2 = Simplify[4 ρ (N1 - 7/2 N2 + 2 N3 - 1/4 N4)] // Expand
Apologies for the large expressions but the particular simplifcations only happen when all four integrands are included with the correct factors. By comparing integrands it's easy to see that int1
and int2
are the same integrals. But for int1
I get
(1 - E^(-L^2 ρ)) L^2
whereas for int2
I get
L^2 - E^(-L^2 ρ) L^2 - (7 L^4 ρ)/4 + 7/2 L^4 ρ Log[L]
I know that the outcome for int1
is correct whereas that for int2
is incorrect.
Any thoughts would be very welcome.
Answer
It's a bug, of course. Mathematica gets dizzy by the "hanging" integration order (I believe).
As simple as it is, it was tough to find out but it gives the right result if you just change the order of limits
:
$Assumptions = {ρ > 0, L > 0};
limits = Sequence[{ux, 0, L}, {uy, 0, ux}, {vx, 0, L}, {vy, 0, vx}];
(* instead of
limits = Sequence[{ux, 0, L}, {vx, 0, L}, {uy, 0, ux}, {vy, 0, vx}]
*)
N1 = Integrate[a1 = Exp[-ρ (ux - uy) (vx - vy)], limits];
N2 = Integrate[a2 = Exp[-ρ (ux - uy) (vx - vy)] (ρ (ux - uy) (vx - vy)), limits];
N3 = Integrate[a3 = Exp[-ρ (ux - uy) (vx - vy)] (ρ (ux - uy) (vx - vy))^2, limits];
N4 = Integrate[a4 = Exp[-ρ (ux - uy) (vx - vy)] (ρ (ux - uy) (vx - vy))^3, limits];
int1 = Integrate[4 ρ (a1 - 7/2 a2 + 2 a3 - 1/4 a4), limits]
int2 = FullSimplify[4 ρ (N1 - 7/2 N2 + 2 N3 - 1/4 N4)]
(*
(1 - E^(-L^2 ρ)) L^2
(1 - E^(-L^2 ρ)) L^2
*)
Comments
Post a Comment