Search code examples
matlabnormal-distributionmontecarlonumerical-integration

Monte Carlo integration of exp(-x^2/2) from x=-infinity to x=+infinity


I want to integrate

f(x) = exp(-x^2/2)

from x=-infinity to x=+infinity

by using the Monte Carlo method. I use the function randn() to generate all x_i for the function f(x_i) = exp(-x_i^2/2) I want to integrate to calculate afterwards the mean value of f([x_1,..x_n]). My problem is, that the result depends on what values I choose for my borders x1 and x2 (see below). My result is going far away from the real value by increasing the value of x1 and x2. Actually the result should be better and better by increasing x1 and x2.

Does someone see my mistake?

Here is my Matlab code

clear all;
b=10;                 % border
x1 = -b;              % left border
x2 = b;               % right border
n = 10^6;             % number of random numbers
x = randn(n,1);
f = ones(n,1);
g = exp(-(x.^2)/2); 
F = ((x2-x1)/n)*f'*g;

The right value should be ~2.5066.

Thanks


Solution

  • Try this:

    clear all;
    b=10;                 % border
    x1 = -b;              % left border
    x2 = b;               % right border
    n = 10^6;             % number of random numbers
    
    x = sort(abs(x1 - x2) * rand(n,1) + x1);
    f = exp(-x.^2/2);
    F = trapz(x,f)
    
    F =
    
        2.5066