0

I want to implement the monte carlo algorithm in matlab in a slightly different way than I've seen anywhere else. I have a working example where I generate a vector of 1000 points and calculate pi from that. Now I want a model where I can add extra points to the same example.

This is my code so far:

a=rand(1000,1); %initial random vectors
b=rand(1000,1);
n=1000;
m=0;   %Number of points inside circle

while true
  x1=a-0.5;
  y1=b-0.5; %cirle has centre at (0.5,0.5)
  r=x1.^2+y1.^2;

  for i=1:n
    if r(i)<=0.25
       m=m+1;
   end
 end
pi=m/(0.25*n);
a=rand(1000,1);
b=rand(1000,1);
n=n+1000;
pause(1);
end

But this doesn't work because of the for-loop, where you check the variables... the r(i) value shoud shift to the next 1000 values I create at the bottom of the while-loop... Somebody knows the solution for this?

wietjes
  • 65
  • 6

1 Answers1

3

The issue is in the range of your for loop. You are going from 1 to n every time but r only has 1000 entries each time through the loop (the size of a and b). This will not be an issue the first time through the loop since n == 1000 and numel(r) == 1000 but in successive times through the while loop, n increases to 2000, 3000, etc. and this exceeds the size of r.

You will want to change your for loop so that you only go from 1 to numel(r).

for k = 1:numel(r)
    if r(k) <= 0.25;
        m = m + 1;
    end
end

Or you can remove the for loop altogether.

m = m + sum(r <= 0.25);
Suever
  • 64,497
  • 14
  • 82
  • 101