How do you build an array to plot it?

I'm writing code using monte carlo integration to approximate pi. I need to made a plot of the convergence to pi as I increase the number of iterations. I'm trying to build an array but I can't figure out how to initialize and plot it all without Matlab saying "Subscript indices must either be real positive integers or logicals." I don't know what initializations I should be using for the array, I've never used an array in Matlab before.

Here is what I have so far

i = 0;
p = 0;
k = 1;
g = 1;
n = max(size(k));
p(i) = (1:k);
Ndom(i) =(1:k);
Nint = 0;
Ndom = 0;
Adom = d*d';

tic
while k < 1000000

        x = (d*rand)-r;         %Run random number sequence for x
        y = (d*rand)-r;         %Run random number sequence for y

         if (x.^2)+(y.^2)<= r.^2

              Nint = Nint + 1;   %Increase number of points within circle 

          else

         end

         if g*10.^g == k

                  i = i + 1;

                  Ndom(i) = k;
                  Aint = (Adom*Nint)./Ndom(i);  
                  p(i) = Aint./(r.^2);

                  %for every 10th iteration calc p, plot it using semilogx or log(x) vs. y
                  %use array

                  semilogx(Ndom(i), p(i));
                  xlabel('Number of points');
                  ylabel('Pi calculation');
                  title('Calculation of pi using Monte Carlo integration');

                  n = n+1;

         else 
              end

           k = k+1;

end

Ndom = k;
Aint = (Adom*Nint)./Ndom;  
p = Aint./ (r.^2);

toc

end

I built the code before trying to figure out how to plot it, just to make sure that the math was correct within an accepted error tolerance, and it works fine. I just can't figure out how to plot it.

[Kirsten Rawls]

I should note, I caught a couple mistakes in my code not related to the question, replaced n = n + 1 with g = g + 1.

[Walter Roberson]

I recommend

   semilogx(Ndom(1:i), p(1:i));

as your current code only plots one point per graph.

Are you sure you only want to plot at iterations 10, 200, 3000, 40000, 500000 ?? Because that's what the test if g*10.^g == k is going to select for.

Answer 1

[Star Strider]

I didn’t run it, but after editing it (to include all your code as [{} code]) and looking at it, observed that the hold function could be your friend here.

[Kirsten Rawls]

Holy crap that is much better!! I don't know why my teacher said to put it into an array and then plot it. This works, except that it's plotting the individual points, but I need then connected. How can I use this, and have my points connected by a line?

[Star Strider]

Thank you!

I would put the plot outside (after) the loop, for example:

figure(1)
semilogx(Ndom, p);
xlabel('Number of points');
ylabel('Pi calculation');
title('Calculation of pi using Monte Carlo integration');integration');

I am not certain exactly where you want to put it (that is to say the inner loop or the outer loop), but that should do what you want, in terms of plotting a continuous (rather than point-wise) curve.

I would also consider replacing:

while k < 1000000

with some sort of convergence criterion (for example the difference between consecutive calculations, perhaps 1E-10) instead of testing the counter itself, that is the equivalent of a for loop.

Just a thought...

[Walter Roberson]

To get them to join,

semilogx(Ndom(1:i), p(1:i));

with which "hold on" would probably not be needed.

Oh yes, put in

drawnow()

after the title() call.

[Star Strider]

Thank you, Walter!