function [vol,err] = volumeSphereND(nd,np,ns)
% Finds the volume of a sphere of unit radius in nd dimensions by MC method
% Input:
%   nd   : number of spacial dimensions
%   np   : number of random points per sample
%   ns   : number of samples used to calculate err
% Output:
%   vol  : sphere volume
%   err  : statistical error in volume
% J.M.Soler. Mar.2021

% Loop on samples
v = zeros(1,ns);
for is = 1:ns
    p = rand(nd,np);              % point coordinates
    p2 = sum(p.^2,1);             % distances to origin squared
    inside = sum(p2<1);           % points inside sphere of unit radius
    volCube = 2^nd;               % volume of hypercube of side [-1:+1]
    v(is) = volCube * inside/np;  % volume of sphere of unit radius
end

% Set output for unit radius
vol = sum(v)/ns;                        % average of volume
err = sqrt(sum((v-vol).^2)/ns/(ns-1));  % statistical error

end % function




