0
$\begingroup$

I want to integrate a function of two variables that is defined via a piecewise function. I can compute the integral analytically as follows:

Integrate[Exp[-y x^α1 ], {x, 0, x0}, Assumptions -> {1 < α1, x0 > 0, y > y0, y0 > 0}] + Integrate[Exp[-x0^(α1 - α2) y x^α2], {x, x0, Infinity}, Assumptions -> {1 < α1, 0 < α2 < 1, x0 > 0, y > y0, y0 > 0, x0 > 0}]

that gives the solution

(x0 ExpIntegralE[(-1 + α2)/α2, x0^α1 y])/α2 +(y^(-1/α1) (Gamma[1/α1] -Gamma[1/α1, x0^α1 y]))/α1

Then, I integrate over y

Integrate[(x0 ExpIntegralE[(-1 + α2)/α2, x0^α1 y])/α2 + (y^(-1/α1) (Gamma[1/α1] - Gamma[1/α1, x0^α1 y]))/α1, {y, y0, Infinity}, Assumptions -> {α1 > 1, 1 > α2 > 0, x0 > 1, y0 > 0}]

and obtain

(y0^(1 - 1/α1) Gamma[1/α1])/(1 - α1) + (y0^((-1 + α1)/α1) Gamma[1/α1, x0^α1 y0])/(-1 + α1) + ( x0 ((E^(-x0^α1 y0)x0^-α1 (α1 - α2))/(-1 + α1) - x0^(-(α1/α2)) y0^((-1 + α2)/α2) Gamma[1/α2, x0^α1 y0]))/(-1 + α2)

I can evaluate it with the parameters that interest me and obtain a result:

With[{α1 = 1.4, α2 = 0.8, x0 = 8, y0 = 1}, ( y0^(1 - 1/α1) Gamma[1/α1])/(1 - α1) + ( y0^((-1 + α1)/α1)   Gamma[1/α1, x0^α1 y0])/(-1 + α1) + (  x0 ((E^(-x0^α1 y0)       x0^-α1 (α1 - α2))/(-1 + α1) -      x0^(-(α1/α2)) y0^((-1 + α2)/α2)      Gamma[1/α2, x0^α1 y0]))/(-1 + α2)]

which evaluates as -3.18998. This can't be true, because the original function was positive in all its domain.

If I try to integrate the function numerically, it warns me about possible singularities.

With[{α1 = 1.4, α2 = 0.8, x0 = 8, R = 10, y0 = 1},  NIntegrate[  Piecewise[{{Exp[-y x^α1 ],      x <= x0}, {Exp[-x0^(α1 - α2) y x^α2],      x > x0}}], {x, 0, Infinity}, {y, y0, Infinity}]]

And gives me a completely different result: 3.21302*10^10

My questions are: - Why does Mathematica produce an analytical solution that appears to be wrong? - How can I obtain a good estimate of the numerical integral? I've tried all Methods and played around with the options, but couldn't get an estimate that I could trust.

Thank you

$\endgroup$
1
  • $\begingroup$ Can you also state the original problem in LaTeX with details on assumptions. $\endgroup$
    – m0nhawk
    Commented Mar 30, 2015 at 9:44

1 Answer 1

2
$\begingroup$

In Mathematica 10 I can evaluate this directly, without any intermediate steps.

With[{α1 = 1.4, α2 = 0.8, x0 = 8, y0 = 1}, 
 Integrate[
  Piecewise[{{Exp[-y x^α1], 
     x <= x0}, {Exp[-x0^(α1 - α2) y x^α2], 
     x > x0}}], {x, 0, Infinity}, {y, y0, Infinity}, 
  PrincipalValue -> True]]
(* ∞ *)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.