# Compute volume of ellipsoid using Monte Carlo methods

I have following code in Mathematica: https://pastebin.com/zRB9vEx8

and I need to compute volume of ellipsoid with given parameters: $$\frac{x^2}{9} + \frac{y^2}{4} + \frac{z^2}{4} = 1$$

Can anyone give me any advice how to do it using code above?

Here a solution independent of your code which calculates the volume with MonteCarlo:

First define the region of the ellipsoid

elli = ImplicitRegion[x^2/9 + y^2/4 + z^2/4 <= 1, {x, y, z}]


Mathematica is able to calculate the volume of this region directly

Volume[elli]
(* 16 Pi*)


or alterntively with Nintegrate

NIntegrate[1, Element[{x, y, z}, elli], Method -> "MonteCarlo"]
(*50.660588822*)

• Or Integrate[1, Element[{x, y, z}, elli]] – Bob Hanlon Apr 18 '20 at 18:48
• Thank you for advice, but it doesn't solve my problem at all :/ Computing volume of ellipsoid with MC is piece of cake. The problem rise when you need to use simulations. Could you give me a hand with applying it into my code where are three common types of MC method? – ajire-bardyn Apr 19 '20 at 8:19
• Sorry , I answered your question concerning volume of ellipsoid. Inspecting your code I doesn't find a part Integrate. It's not so easy to guess ("problem rise when you need to use simulations") what you realy ask for ! – Ulrich Neumann Apr 19 '20 at 10:40

Make $$n$$ random points in the cube of side $$6$$. Then count how many are in the ellipse. Divide that by $$n$$ and multiply by the volume of the cube.

n = 10000;
mylist = RandomVariate[UniformDistribution[{{-3, 3}, {-3, 3}, {-3, 3}}],
n];
6^3 (Count[mylist, u_ /; u[[1]]^2/9 + u[[2]]^2/4 + u[[3]]^2/4 < 1]/n)


(* 63153/1250 *)

• You can make things computationally easier by only considering one octant, and multiplying by 8 at the end. – J. M.'s ennui Apr 20 '20 at 4:11
• But then to get the same accuracy you need $8$ times the number of points. No difference. – David G. Stork Apr 20 '20 at 4:57

In a comment, David claims that "to get the same accuracy you need 8 times the number of points."

Here's one way to prove or disprove the claim. Compare the two Monte Carlo routines that follow:

mcEllipsoid[{a_, b_, c_}, n_Integer?Positive] := Module[{s = 0, r, tr},
tr = RescalingTransform[ConstantArray[{0, 1}, 3], {{0, a}, {0, b}, {0, c}}];
Do[r = tr[RandomReal[1, {3}]];
If[(#.#) &[r/{a, b, c}] <= 1, s++], {n}];
(8 a b c) s/n]

mcEllipsoid2[{a_, b_, c_}, n_Integer?Positive] := Module[{s = 0, r, tr},
tr = RescalingTransform[ConstantArray[{0, 1}, 3], {{-a, a}, {-b, b}, {-c, c}}];
Do[r = tr[RandomReal[1, {3}]];
If[(#.#) &[r/{a, b, c}] <= 1, s++], {n}];
(8 a b c) s/n]


where the first routine only considers one octant and multiplies by $$8$$ at the end, while the other considers the entire ellipsoid and the cuboid enclosing it.

Some limited testing seems to indicate that for the same value of n, the estimates obtained by mcEllipsoid[{3., 2, 2}, 1*^7] are closer to the true value of the volume than the corresponding results of mcEllipsoid2[{3., 2, 2}, 1*^7].

But I'll let someone else with a better computer do more stringent tests. ;)