My guess is this homework was chosen to show you how inefficient Monte Carlo can be when used in a high number of dimensions. THe curse of dimensionality strikes here.
The known volume of a unit 10-sphere is simple enough to compute. Google is your friend here, but I'll just use a tool I've posted on the FEX.
format long g
spheresegmentvolume([],10)
ans =
2.55016403987735
Use of Monte Carlo to compute volume is a dartboard approach. Thus throw a series darts at a hypercube of edge length 2 that just contains the unit 10-sphere. Count the number of times the dart falls inside the sphere, and divide by the total number of darts thrown. Multiply by the known 10-cube volume, and you get an approximation of the 10-sphere volume.
The trick is to do this efficiently. Do NOT use a loop.
% the dimension of the problem
dim = 10
% the volume of a 10-cube of edge length 2
cubevol = 2^dim;
% sample size
N = 10e6;
% draw the samples
X = 2*rand(N,dim) - 1;
% Distance from the origin for each sample
% Note that I could have skipped the square root, since this is a unit sphere.
R = sqrt(sum(X.^2,2));
% count the number of hits
C = sum(R <= 1)
C =
24925
% compute the estimated volume
V = C/N*cubevol
V =
2.55232
As you can see, the computed volume was not terribly inaccurate compared to the known volume of 2.5502… 3 significant digits is as much as we can expect.
To be honest, if I were your instructor, I'd have specified a 20 or even 50 dimensional sphere, to make that hit rate vanishingly small. For the 10-sphere, the fraction of darts that actually hit the sphere is only about 0.25%. But that is a large enough fraction that we got almost 3 digits of accuracy from our pretty large sample size.
hitrate = spheresegmentvolume([],10)/cubevol
hitrate =
0.00249039457019272
A nice thing about Monte Carlo is is allows you to compute a measure of your uncertainty in the final estimate.
A binomial random variable with p = 0.00249039… has variance n*p*(1-p). So if we wish to put +/- 2*sigma confidence limits around the number of hits we should have expected:
C + [-2 2]*sqrt(N*hitrate*(1-hitrate))
ans =
24609.7735731207 25240.2264268793
And of course, the final volume will have confidence bounds of:
V + cubevol/N*[-2 2]*sqrt(N*hitrate*(1-hitrate))
ans =
2.52004081388756 2.58459918611244
Again, learn how to avoid loops when you can do so.
