10
$\begingroup$
Integrate[(4 a b/Pi) (a^2 + b^2 - 2 a b Cos[c])^(1/2), {a, 0, 1}, {b, 
  0, 1}, {c, 0, Pi}]

enter image description here

I'm using basic plan, it gives me the result like that. The approximate is 0.9054... but I want exact value with closed form.

$\endgroup$
2
  • 1
    $\begingroup$ I'm afraid it's impossible $\endgroup$
    – yode
    May 17 at 8:35
  • 1
    $\begingroup$ Actually it is possible, see my answer below, and the result is exactly $128/(45\pi)$. $\endgroup$
    – Hans Olo
    May 17 at 9:14

3 Answers 3

23
$\begingroup$

Doing the Integrals separately for a, b and c gives the answer $\frac{128}{45 \pi}\sim 0.905415$, which agrees with the numerical estimate from NIntegrate:

NIntegrate[(4 a b/Pi) (a^2 + b^2 - 2 a b Cos[c])^(1/2), {a, 0, 1}, {b, 0, 1}, {c, 0, Pi}]

0.905415

In a nutshell, the trick is to do the integrals separately as indefinite and then take the limits properly. This is a somewhat straight-forward integral and the trick works, in more general cases, caution is advised!

In more detail: First do the integral over $a$:

inta = Integrate[(4 a b/Pi) (a^2 + b^2 - 2 a b Cos[c])^(1/2), a]

(1/(3 [Pi]))b (Sqrt[a^2 + b^2-2 a b Cos[c]] (4 a^2 + b^2 - 2 a b Cos[c] - 3 b^2 Cos[2 c]) + 6 b^3 Cos[c] Log[a - b Cos[c] + Sqrt[a^2 + b^2 - 2 a b Cos[c]]] Sin[c]^2)

and then the one over b:

intb = Integrate[(inta /. a -> 1) - (inta /. a -> 0), b] // PowerExpand

(1/(15 [Pi]))(b^5 (-1 + 3 Cos[2 c]) + Sqrt[1 + b^2 - 2 b Cos[c]] (1 + 8 b^2 + b^4 - 2 (b + b^3) Cos[c] - 3 (1 + b^4) Cos[2 c]) - 6 b^5 Cos[c] Log[b - b Cos[c]] Sin[c]^2 + 6 Cos[c] Log[b - Cos[c] + Sqrt[1 + b^2 - 2 b Cos[c]]] Sin[c]^2 + 6 b^5 Cos[c] Log[1 - b Cos[c] + Sqrt[1 + b^2 - 2 b Cos[c]]] Sin[c]^2)

Then the one over c:

intc = Limit[intb, b -> 1] - (Series[intb, {b, 0, 0}] // FullSimplify // Normal) //  FullSimplify

-(1/(15 [Pi])) 2 (1 - 5 Sqrt[2 - 2 Cos[c]] + 3 (-1 + Sqrt[2 - 2 Cos[c]]) Cos[2 c] + 2 Cos[c] (Sqrt[2 - 2 Cos[c]] + 3 (Log[1 - Cos[c]] - Log[1 + Sqrt[2 - 2 Cos[c]] - Cos[c]]) Sin[c]^2))

Integrate[intc, {c, 0, Pi}]
% // N

128/(45 [Pi])

0.905415

Note: In some cases one has to use Limit or even a series expansion to extract the correct value of the definite integral.

$\endgroup$
3
  • 1
    $\begingroup$ I'd like to inform that Integrate[(4 a b/Pi) (a^2 + b^2 - 2 a b Cos[c])^(1/2), {a, 0, 1}, {b, 0, 1}, Assumptions -> c > -Infinity] produces a correct huge conditional expression in several hours as NIntegrate[%[[1]], {c, 0, Pi}] confirms with the result 0.905415.Unfortunately, Integrate[%%[[1]], {c, 0, Pi}] returns the input (I think because of Arg[-E^(-I c) (-1 + E^(I c))^2] and similar expressions.). $\endgroup$
    – user64494
    May 17 at 17:30
  • 3
    $\begingroup$ For those who wonder why "caution is advised" when using this method, read this blog post from a Mathematica developer (which I never miss an opportunity to link to.) $\endgroup$ May 17 at 18:40
  • $\begingroup$ @MichaelSeifert Great link thanks! I was playing fast and loose with the calculus in the post, hence just the warning ;-) $\endgroup$
    – Hans Olo
    May 17 at 19:12
10
$\begingroup$

The unevaluated integral that Mathematica produced in the OP can be massaged a bit. First, note that due to symmetry, we can instead consider the integral

(16/π) Integrate[a b (a + b) EllipticE[(4 a b)/(a + b)^2], {a, 0, 1}, {b, 0, a}]

and using the conventional definition for EllipticE, we have

(16/π) Integrate[a b (a + b) Sqrt[1 - (4 a b)/(a + b)^2 Sin[t]^2],
                 {t, 0, π/2}, {a, 0, 1}, {b, 0, a}]

Due again to symmetry, this is equivalent to the integral

(8/π) Integrate[a b (a + b) Sqrt[1 - (4 a b)/(a + b)^2 Sin[t]^2],
                {t, -π/2, π/2}, {a, 0, 1}, {b, 0, a}]

and letting w == Sin[t], we have

(8/π) Integrate[a b (a + b) Sqrt[(1 - (4 a b w^2)/(a + b)^2)/(1 - w^2)],
                {w, -1, 1}, {a, 0, 1}, {b, 0, a}]

which happens to evaluate pretty quickly to 128/(45 π).


It should be noted that the equivalent integrals

16/π Integrate[a b (a + b) Sqrt[(1 - (4 a b)/(a + b)^2 w^2)/(1 - w^2)],
               {w, 0, 1}, {a, 0, 1}, {b, 0, a}]

and

4/π Integrate[a b (a + b) Sqrt[(1 - (4 a b)/(a + b)^2 w^2)/(1 - w^2)],
              {w, -1, 1}, {a, 0, 1}, {b, 0, 1}]

do not seem to finish as quickly (I waited ~ 2 minutes for each before giving up.)

$\endgroup$
8
$\begingroup$

If you regard symmetry of integrand (two times integral from 0<b<a) and rearange to do integration over c at last, you get the result in a non splitted integration in 10 seconds.

Original integrations does not work may be because integrate can not handle the case a == b.

f[a_, b_, c_] = (4 a b/Pi) (a^2 + b^2 - 2 a b Cos[c])^(1/2);

f[a, b, c] == f[b, a, c]     (*   True   *)

{(int = 2*Integrate[f[a, b, c], {c, 0, Pi}, {a, 0, 1}, {b, 0, a}]), 
  int // N} // Timing

(*   {10.141, {128/(45 \[Pi]), 0.905415}}   *)
$\endgroup$
1
  • $\begingroup$ +1. Simply and strongly. On my comp {138.766, {128/(45 \[Pi]), 0.905415}}. $\endgroup$
    – user64494
    May 18 at 8:03

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.