Great Michael Trott & six-fold integral

• To: mathgroup at smc.vnet.net
• Subject: [mg30124] Great Michael Trott & six-fold integral
• From: seidovzf at yahoo.com (Zakir F. Seidov)
• Date: Sun, 29 Jul 2001 21:26:16 -0400 (EDT)
• Organization: The Math Forum
• References: <9jtk3e\$87k\$1@smc.vnet.net>
• Sender: owner-wri-mathgroup at wolfram.com

```Yes, Michael,

Your cunning usage of Mathematica is really amazing.

I'm not sure that ordinary user (as me, e.g.)
could understand your session, but I'd try...

I've just rerun your session and  got result in seconds!!

Unfortunately (but as usually) I've made a lot of misprints:
missing closing ")", sign and factor 2...

Thank you very much indeed for your great posting,

But I dare ask:
if you can, with your method, also solve
the general case of  arbitrary limits

Int_0 ^a Int_0 ^b Int_0 ^c Int_0^ A Int_0 ^B Int_0 ^C

which was solved in my preprint, please?!

Zakir F. Seidov

seidov at bgumail.bgu.ac.il

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Subject: [mg30124]      Re: six-fold integral
Author:       Michael Trott <mtrott at wolfram.com>
Organization: Steven M. Christensen and Associates, Inc and
MathTensor, Inc.

Your result

(2 3^(1/2)-2^(1/2)-1)/5 + pi/3 + ln[(2^(1/2)-1)(2-3^(1/2)]

misses a closing ')'. Assuming you mean

(2 3^(1/2)-2^(1/2)-1)/5+Pi/3 + Log[(2^(1/2)-1)(2-3^(1/2))] // N,

this is unfortunately negative (-0.941156).

Using the joint probability density for (x - y) it is straightforward
to calculate the following exact result in just a few lines:

In[1]:= logRewrite[expr_] := FixedPoint[PowerExpand[# //.
Log[a_] :> Log[Factor[Together[a]]]]&, TrigToExp[expr]]

In[2]:= indefiniteIntegral = Integrate[# 2(1 - z), z]& /@
Expand[logRewrite[
Integrate[Integrate[1/Sqrt[x^2 + y^2 + z^2] 2(1 - x), x] 2(1 - y),
y]]];

In[3]:= replaceLimits[expr_, xyz_] :=((# /. xyz -> 1) - (# /. xyz ->
0))&[
Expand[logRewrite[expr]] /. xyz^n_?Positive Log[xyz] -> 0]

In[4]:= Fold[replaceLimits, indefiniteIntegral, {z, y, x}] //
Together;

In[5]:= FullSimplify[%]

Out[5]= -((2*Pi)/3) + Log[3 + 2*Sqrt[2]] +
1/10*(4 + 4*Sqrt[2] - 8*Sqrt[3] - 25*Log[2] + 50*Log[1 + Sqrt[3]] -
5*Log[2 + Sqrt[3]])

In[6]:= N[%]

Out[6]= 1.88231

This agrees with a Monte-Carlo check:

In[7]:=
SeedRandom[888];
Compile[{{n, _Integer}},
Module[{sum = 0.},
Do[sum = sum + 1/Sqrt[(Random[] - Random[])^2 +
(Random[] - Random[])^2 +
(Random[] - Random[])^2],
{n}]; sum/n]][10^6]

Out[7]= 1.88295

```

• Prev by Date: Re: Graphics from Mathematica to MS Word
• Next by Date: Two factors of (10^71-1)/9 = R71
• Previous by thread: Re: Graphics from Mathematica to MS Word
• Next by thread: Two factors of (10^71-1)/9 = R71