3
$\begingroup$

I want to Plot a definite integral

myIntegral[x_] := NIntegrate[
  Sqrt[(a - b)^2 + (c - d)^2]/(a - b)
  , {a, 0, x}
  , {b, x, 1}
  , {c, 0, 1}
  , {d, 0, 1}
  ]

Plot[myIntegral[x], {x, 0, 1}]

I am using Mathematica 11.3. When I evaluate the code, Mathematica immediately gives some warnings about numerical integration converging too slowly, then it seems stuck in computation. I have tried some options of NIntegrate, for example as suggested here, but have had no success.

Is there a way to obtain the result of Plot in an acceptable time (e.g. some minutes)?

$\endgroup$
3
  • 1
    $\begingroup$ Method -> "MonteCarlo" allows you to plot the function in less than a minute, but the plot is not perfectly smooth: i.imgur.com/sEuL6oN.png increasing the number of integration points makes the plot slower but nicer. $\endgroup$ Jun 19 '18 at 14:56
  • 1
    $\begingroup$ Thank you for the hint. It seems that Method -> "QuasiMonteCarlo" can give a smooth plot in this case. In addition, it is slightly faster (few seconds less). $\endgroup$
    – renato
    Jun 19 '18 at 16:15
  • 1
    $\begingroup$ @renato Self-answers are encouraged on StackExchange, so you would be welcome to post your comment as an answer! $\endgroup$
    – MarcoB
    Jun 19 '18 at 17:46
4
$\begingroup$

The integral is within reach of Integrate

Assuming[0 < x < 1, 
  Simplify[
    Integrate[Sqrt[(a-b)^2 + (c-d)^2]/(a-b), {a,0,x}, {b,x,1}, {c,0,1}, {d,0,1}]]]

which gives you

(-11*Sqrt[2] + 11*Sqrt[2+(-2+x)*x] + x*(-8-7*Sqrt[2+(-2+x)*x] + 13*Sqrt[1+x^2] +
 2*x*(6-3*Sqrt[2+(-2+x)*x] + x*(-4+2*x+Sqrt[2+(-2+x)*x] - Sqrt[1+x^2]))) +
 18*x^2*ArcCoth[Sqrt[2+(-2+x)*x]] + 9*ArcSinh[1] - 3*ArcSinh[1-x] -
 3*ArcSinh[x] + 6*Log[1-x] + 9*x*Log[-1+Sqrt[2+(-2+x)*x]] + 3*((-2+x)*Log[1 +
 Sqrt[2 + (-2 + x)*x]] - 4*x*Log[-(((-1 + x)*(1 + Sqrt[1 + x^2]))/x)] +
 x^3*(Log[2+x^2+2*Sqrt[1+x^2]] + 2*Log[(1-x)/(x+x*Sqrt[2-2*x+x^2])])))/36

and that will be far faster than doing thousands of NIntegrate inside Plot

$\endgroup$
3
$\begingroup$

You can first calculate the integral over c and d to get:

res = Integrate[{Sqrt[t + (c - d)^2]}, {c, 0, 1}, {d, 0, 1}];
out=FullSimplify[res, {a >= 0, b >= 0, t > 0}]

{1/3 (2 t^(3/2) + Sqrt[1 + t] - 2 t Sqrt[1 + t] + 3 t ArcCoth[Sqrt[1 + t]])}

Now put (a-b)^2 for t and dividing by (a-b).

out2=FullSimplify[(out /. t -> (a - b)^2)/(a - b), {a >= 0,b >= 0, t >= 0}];

So you can then write a simplified version of you integral as:

myIntegral[x_] := NIntegrate[out2, {a, 0, x}, {b, x, 1}]

Plot[myIntegral[x], {x, 0, 1}]

enter image description here

However, to get the first integral i put in some numbers for t=10,100,1000... and then guessed the solution. That was a bit faster, but mathematica is also able to find this solution if you wait for a while.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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