% Use Monte Carlo method to approximate triple integral \int\int\int
% sqrt(x^2 + z^2) dV, where integration domain is bounded by paraboloid
%  y = x^2 + z^2 and the plane y = 4
% The analytic solution is 128*pi/15

N = input('Enter number of sample points: '); 
f=inline('sqrt(x*x+z*z)');    %integrand function f = sqrt(x^2 + z^2) 

volumeofBox=4*4*4;   % the box is [-2, 2]*[0,4]*z[-2,2];
vol =0; mass = 0;

for i =1:N,         % loop over all sample points
 x=-2+4*rand;       % Gererate a point in the sourrounding box
 z=-2+4*rand;
 y=0+4*rand;
 
 if x*x+z*z <=y,               
  vol = vol+1;
  mass = mass + f(x,z);
 end;
end;

volumeofobj = (vol/N)* volumeofBox % fraction of pts inside times vol of box
massofObj = (mass/N)*volumeofBox   % Average mass times volume of box