Monte Carlo methods: Difference between revisions
(added ocaml) |
(Ada solution added) |
||
Line 5: | Line 5: | ||
Write a function to run a simulation like this with a variable number of random points to select. Also, show the results of a few different sample sizes. |
Write a function to run a simulation like this with a variable number of random points to select. Also, show the results of a few different sample sizes. |
||
=={{header|Ada}}== |
|||
<ada> |
|||
with Ada.Text_IO; use Ada.Text_IO; |
|||
with Ada.Numerics.Float_Random; use Ada.Numerics.Float_Random; |
|||
procedure Test_Monte_Carlo is |
|||
Dice : Generator; |
|||
function Pi (Throws : Positive) return Float is |
|||
Inside : Natural := 0; |
|||
begin |
|||
for Throw in 1..Throws loop |
|||
if Random (Dice) ** 2 + Random (Dice) ** 2 <= 1.0 then |
|||
Inside := Inside + 1; |
|||
end if; |
|||
end loop; |
|||
return 4.0 * Float (Inside) / Float (Throws); |
|||
end Pi; |
|||
begin |
|||
Put_Line (" 10_000:" & Float'Image (Pi ( 10_000))); |
|||
Put_Line (" 100_000:" & Float'Image (Pi ( 100_000))); |
|||
Put_Line (" 1_000_000:" & Float'Image (Pi ( 1_000_000))); |
|||
Put_Line (" 10_000_000:" & Float'Image (Pi ( 10_000_000))); |
|||
Put_Line ("100_000_000:" & Float'Image (Pi (100_000_000))); |
|||
end Test_Monte_Carlo; |
|||
</ada> |
|||
The implementation uses built-in uniformly distributed on [0,1] random numbers. Note that the accuracy of the result depends on the quality of the pseudo random generator: its circle length and correlation to the function being simulated. Sample output: |
|||
<pre> |
|||
10_000: 3.13920E+00 |
|||
100_000: 3.14684E+00 |
|||
1_000_000: 3.14197E+00 |
|||
10_000_000: 3.14215E+00 |
|||
100_000_000: 3.14151E+00 |
|||
</pre> |
|||
=={{header|Java}}== |
=={{header|Java}}== |
||
<java>public class MC { |
<java>public class MC { |
Revision as of 10:02, 28 September 2008
![Task](http://static.miraheze.org/rosettacodewiki/thumb/b/ba/Rcode-button-task-crushed.png/64px-Rcode-button-task-crushed.png)
You are encouraged to solve this task according to the task description, using any language you may know.
A Monte Carlo Simulation is a way of approximating the value of a function where calculating the actual value is difficult or impossible. It uses random sampling to define constraints on the value and then makes a sort of "best guess."
A simple Monte Carlo Simulation can be used to calculate the value for pi. If you had a circle and a square where the length of a side of the square was the same as the diameter of the circle, the ratio of the area of the circle to the area of the square would be pi/4. So, if you put this circle inside the square and select many random points inside the square, the number of points inside the circle divided by the number of points inside the square and the circle would be approximately pi/4.
Write a function to run a simulation like this with a variable number of random points to select. Also, show the results of a few different sample sizes.
Ada
<ada> with Ada.Text_IO; use Ada.Text_IO; with Ada.Numerics.Float_Random; use Ada.Numerics.Float_Random;
procedure Test_Monte_Carlo is
Dice : Generator; function Pi (Throws : Positive) return Float is Inside : Natural := 0; begin for Throw in 1..Throws loop if Random (Dice) ** 2 + Random (Dice) ** 2 <= 1.0 then Inside := Inside + 1; end if; end loop; return 4.0 * Float (Inside) / Float (Throws); end Pi;
begin
Put_Line (" 10_000:" & Float'Image (Pi ( 10_000))); Put_Line (" 100_000:" & Float'Image (Pi ( 100_000))); Put_Line (" 1_000_000:" & Float'Image (Pi ( 1_000_000))); Put_Line (" 10_000_000:" & Float'Image (Pi ( 10_000_000))); Put_Line ("100_000_000:" & Float'Image (Pi (100_000_000)));
end Test_Monte_Carlo; </ada> The implementation uses built-in uniformly distributed on [0,1] random numbers. Note that the accuracy of the result depends on the quality of the pseudo random generator: its circle length and correlation to the function being simulated. Sample output:
10_000: 3.13920E+00 100_000: 3.14684E+00 1_000_000: 3.14197E+00 10_000_000: 3.14215E+00 100_000_000: 3.14151E+00
Java
<java>public class MC { public static void main(String[] args) { System.out.println(getPi(10000)); System.out.println(getPi(100000)); System.out.println(getPi(1000000)); System.out.println(getPi(10000000)); System.out.println(getPi(100000000));
} public static double getPi(int numThrows){ int inCircle= 0; for(int i= 0;i < numThrows;i++){ //a square with a side of length 2 centered at 0 has //x and y range of -1 to 1 double randX= (Math.random() * 2) - 1;//range -1 to 1 double randY= (Math.random() * 2) - 1;//range -1 to 1 //distance from (0,0) = sqrt((x-0)^2+(y-0)^2) double dist= Math.sqrt(randX * randX + randY * randY); if(dist < 1){//circle with diameter of 2 has radius of 1 inCircle++; } } return 4.0 * inCircle / numThrows; } }</java> Output:
3.1396 3.14256 3.141516 3.1418692 3.14168604
OCaml
<ocaml>let get_pi throws =
let rec helper i count = if i = throws then count else let rand_x = Random.float 2.0 -. 1.0 and rand_y = Random.float 2.0 -. 1.0 in let dist = sqrt (rand_x *. rand_x +. rand_y *. rand_y) in if dist < 1.0 then helper (i+1) (count+1) else helper (i+1) count in float (4 * helper 0 0) /. float throws</ocaml>
Example:
# get_pi 10000;; - : float = 3.15 # get_pi 100000;; - : float = 3.13272 # get_pi 1000000;; - : float = 3.143808 # get_pi 10000000;; - : float = 3.1421704 # get_pi 100000000;; - : float = 3.14153872
Python
Python 2.6rc2 (r26rc2:66507, Sep 18 2008, 14:27:33) [MSC v.1500 32 bit (Intel)] on win32
IDLE 2.6rc2
<python>>>> import random, math >>> throws = 1000 >>> 4.0 * sum(1 for p in xrange(throws) if math.hypot(*[random.random()*2-1 for q in [0,1]]) < 1) / float(throws) 3.1520000000000001 >>> throws = 1000000 >>> 4.0 * sum(1 for p in xrange(throws) if math.hypot(*[random.random()*2-1 for q in [0,1]]) < 1) / float(throws) 3.1396359999999999 >>> throws = 100000000 >>> 4.0 * sum(1 for p in xrange(throws) if math.hypot(*[random.random()*2-1 for q in [0,1]]) < 1) / float(throws) 3.1415666400000002</python>