Write a Matlab script that simulates the throwing of a large

Write a Matlab script that simulates the throwing of a large number of darts that land anywhere on a square with equal probability. The square has a

side 2 and a center (0,0). Use your simulation to produce an estimate of .

also using 500 Monte Carlo simulations, generate a sequence of approximations of and display the result in a graph.

Solution

Matlab script that simulates the throwing of a large number of darts that land anywhere on a square with equal probability.

Computing Area of a Shape by Random Sampling

The dart has been thrown repeatedly N times, randomly at a circular dart board (radius r=1) in a square cabinet of side 2r.

If a dart lands outside the cabinet we throw again without counting.

The key here is that the darts are thrown with uniform probability of landing at any point in the square cabinet.

N large enough we expect that the ratio of board hits to total throws will tend towards:

                      Area(board)/Area(cabinet)

Recalling formulae for the area of a circle and the area of a square this ratio is

which is about 0.79.

We would make N throws and count the M successful hits. Then we would use that the ratio of the areas tends to the ratio of hits to total throws:

Area(board)/Area(cabinet) ~= M/N

Area(board) ~= M*Area(cabinet)/N = M*4*(r^2)/N

Code:

>> N = 10; // N random points in the x-y plane

>> x = rand(N,1);

>> y = rand(N,1);

>> radii = sqrt(x.^2+y.^2); // we used .^ to compute the entry-wise square of each entry of x and y

>>

>> hits = find(radii<1); // We use find to locate the random nodes that have radius less than one

>>

>> M = length(hits);

>>

>> approxarea = (M/N)*4 // we compute the ratio of the number of nodes with radius less than one and the total number of random nodes.

approxarea =

   2.800000000000000

>>

We can also compute using for loop.

>> for N=[10 100 1000 1e4 1e5 1e6 1e7]

        x = rand(N,1);

        y = rand(N,1);

        radii = sqrt(x.^2+y.^2);

        hits = find(radii<1);

        M = length(hits);

        approxarea = (M/N)*4

    end

approxarea =

   2.800000000000000

approxarea =

   3.360000000000000

approxarea =

   3.160000000000000

approxarea =

   3.127200000000000

approxarea =

   3.147920000000000

approxarea =

   3.141000000000000

approxarea =

   3.142194400000000

>>