I need to write MATLAB code that will integrate over a R^5 hypercube using Monte Carlo. I have a basic algorithm that works when I have a generic function. But the function I need to integrate is:
∫dA
A is an element of R^5.
If I had ∫f(x)dA then I think my algorithm would work.
Here is the algorithm:
% Writen by Jerome W Lindsey III
clear;
n = 10000;
% Make a matrix of the same dimension
% as the problem. Each row is a dimension
A = rand(5,n);
% Vector to contain the solution
B = zeros(1,n);
for k = 1:n
% insert the integrand here
% I don't know how to enter a function {f(1,n), f(2,n), … f(5n)} that
% will give me the proper solution
% I threw in a function that will spit out 5!
% because that is the correct solution.
B(k) = 1 / (2 * 3 * 4 * 5);
end
mean(B)
%Question 2, problem set 1
% Writen by Jerome W Lindsey III
clear;
n = 10000;
% Make a matrix of the same dimension
% as the problem. Each row is a dimension
A = rand(5,n);
% Vector to contain the solution
B = zeros(1,n);
for k = 1:n
% insert the integrand here
% this bit of code works as the integrand
if sum(A(:,k)) < 1
B(k) = 1;
end
end
clear k;
clear A;
% Begin error estimation calculations
std_mc = std(B);
clear n;
clear B;
% using the error I calculate a new random
% vector of corect length
N_new = round(std_mc ^ 2 * 3.291 ^ 2 * 1000000);
A_new = rand(5, N_new);
B_new = zeros(1,N_new);
clear std_mc;
for k = 1:N_new
if sum(A_new(:,k)) < 1
B_new(k) = 1;
end
end
clear k;
clear A_new;
% collect descriptive statisitics
M_new = mean(B_new);
std_new = std(B_new);
MC_new_error_999 = std_new * 3.921 / sqrt(N_new);
clear N_new;
clear B_new;
clear std_new;
% Display Results
disp('Integral in question #2 is');
disp(M_new);
disp(' ');
disp('Monte Carlo Error');
disp(MC_new_error_999);