# Monte Carlo integration of functions to the 3rd power

37 views (last 30 days)
Denise Veldhuis on 2 Dec 2020
Edited: Denise Veldhuis on 6 Dec 2020
I have to create a monte carlo integration of the function above and this is what i've tried so far;
f = @(t) (-0.1).*(t.^3) + (0.58).*(t.^2).*(cosh((-1)./(t+1))) + exp((t)./(2.7)) + 9.6;
n = rand(10,1); %deze waardes moeten in de domain van f passen
upper_lim = 12;
lower_lim = 0;
y_max = 200;
N = 1000; %total number of points
N0 = 66; %number of points underneath the graph
y = (10/n).*(y_max.*(upper_lim - lower_lim))
monte_carlo = sum(y);
plot(n,y); hold on;
end
This gives me 2 straight lines but I am struggeling to translate this into the correct monte carlo formation. I have succeeded in creating a random graph of poitns which one would need but can't seem to get it coded into the function above, namely;
f = @(t)(-0.1).*(t.^3) + (0.58).*(t.^2).*(cosh((-1)./(t+1))) + exp((t)./(2.7)) + 9.6;
%yt = sort(rand(2, N)), 1, 'descend');
%t = 12*yt(1.:);
n = rand(10,1); %deze waardes moeten in de domain van f passen
N = 1000; %total number of points
t = 12*rand(1, N); %domain x = [0,12]
y = 200*rand(1, N); %y = [0,200]
plot(t,y,'.')
y_max = 200;
%N0 = ; %number of pointa underneed the graph
y1 = (n_max/n).*(y_max.*(max(t) - min(t)))
monte_carlo = sum(y1);
end
If someone can help or point me in the right direction it would be much appreciated!

Alan Stevens on 2 Dec 2020
Edited: Alan Stevens on 2 Dec 2020
Here's a simple MC calculation of the integration of your function
f = @(t) -0.1*t.^3 + 0.58*t.^2.*cosh(-1./(t+1)) + exp(t/2.7) + 9.6;
t = 0:0.1:12;
% Always a good idea to plot the function first
plot(t,f(t)),grid
%The function values lie below 20, so scatter random points
% with t values between 0 and 12, and y values between 0 and 20.
% The total rectangular area is
Arect = 12*20;
% The approximate area under the curve is given by
N = 10^4; % Number of MC trials
t = 12*rand(N,1); % Random values of t
y = 20*rand(N,1); % Random values of y
% Now find values of y that are less than or equal to the function.
ylo = y(y<=f(t));
% Approximate area
A = Arect*numel(ylo)/N;
disp(A)
(This could be made more efficient by noticing that the curve doesn't drop below 5 over the domain of interest.)

#### 1 Comment

Denise Veldhuis on 6 Dec 2020
Thank you so much! I changed some parts of it but I see what I was messing up, thank you!

R2020a

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by