Monte Carlo methods: Difference between revisions

m
m (→‎{{header|Wren}}: Minor tidy)
 
(16 intermediate revisions by 11 users not shown)
Line 29:
 
=={{header|11l}}==
<langsyntaxhighlight lang="11l">F monte_carlo_pi(n)
V inside = 0
L 1..n
Line 38:
R 4.0 * inside / n
 
print(monte_carlo_pi(1000000))</langsyntaxhighlight>
 
{{out}}
Line 46:
 
=={{header|360 Assembly}}==
<langsyntaxhighlight lang="360asm">* Monte Carlo methods 08/03/2017
MONTECAR CSECT
USING MONTECAR,R13 base register
Line 127:
WP DS PL16 packed decimal 16
YREGS
END MONTECAR</langsyntaxhighlight>
{{out}}
<pre>
Line 139:
{{libheader|Action! Tool Kit}}
{{libheader|Action! Real Math}}
<langsyntaxhighlight Actionlang="action!">INCLUDE "H6:REALMATH.ACT"
 
DEFINE PTR="CARD"
Line 219:
Test(1000)
Test(10000)
RETURN</langsyntaxhighlight>
{{out}}
[https://gitlab.com/amarok8bit/action-rosetta-code/-/raw/master/images/Monte_Carlo_methods.png Screenshot from Atari 8-bit computer]
Line 231:
 
=={{header|Ada}}==
<langsyntaxhighlight lang="ada">with Ada.Text_IO; use Ada.Text_IO;
with Ada.Numerics.Float_Random; use Ada.Numerics.Float_Random;
 
Line 253:
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;</langsyntaxhighlight>
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.
Line 271:
 
{{works with|ELLA ALGOL 68|Any (with appropriate job cards) - tested with release 1.8.8d.fc9.i386}}
<langsyntaxhighlight lang="algol68">PROC pi = (INT throws)REAL:
BEGIN
INT inside := 0;
Line 286:
print ((" 1 000 000:",pi ( 1 000 000),new line));
print ((" 10 000 000:",pi ( 10 000 000),new line));
print (("100 000 000:",pi (100 000 000),new line))</langsyntaxhighlight>
{{out}}
<pre>
Line 298:
=={{header|Arturo}}==
 
<langsyntaxhighlight lang="rebol">Pi: function [throws][
inside: new 0.0
do.times: throws [
Line 307:
loop [100 1000 10000 100000 1000000] 'n ->
print [pad to :string n 8 "=>" Pi n]</langsyntaxhighlight>
 
{{out}}
Line 320:
{{AutoHotkey case}}
Source: [http://www.autohotkey.com/forum/topic44657.html AutoHotkey forum] by Laszlo
<langsyntaxhighlight lang="autohotkey">
MsgBox % MontePi(10000) ; 3.154400
MsgBox % MontePi(100000) ; 3.142040
Line 333:
Return 4*p/n
}
</syntaxhighlight>
</lang>
 
=={{header|AWK}}==
<syntaxhighlight lang="awk">
<lang AWK>
# --- with command line argument "throws" ---
 
Line 348:
Pi = 3.14333
 
</syntaxhighlight>
</lang>
 
=={{header|BASIC}}==
{{works with|QuickBasic|4.5}}
{{trans|Java}}
<langsyntaxhighlight lang="qbasic">DECLARE FUNCTION getPi! (throws!)
CLS
PRINT getPi(10000)
Line 374:
NEXT i
getPi = 4! * inCircle / throws
END FUNCTION</langsyntaxhighlight>
{{out}}
<pre>
Line 383:
</pre>
 
==={{header|BASIC256}}===
{{works with|basic256|1.1.4.0}}
<langsyntaxhighlight lang="basic">
# Monte Carlo Simulator
# Determine value of pi
Line 407:
next i
 
print float(4*in_c/tosses)</langsyntaxhighlight>
{{out}}
<pre>
Line 417:
</pre>
 
====Other solution:====
=={{header|BBC BASIC}}==
<syntaxhighlight lang="freebasic">print " Number of throws Ratio (Pi) Error"
<lang bbcbasic> PRINT FNmontecarlo(1000)
 
for pow = 2 to 8
n = 10 ^ pow
pi_ = getPi(n)
error_ = 3.141592653589793238462643383280 - pi_
print rjust(string(int(n)), 17); " "; ljust(string(pi_), 13); " "; ljust(string(error_), 13)
next
end
 
function getPi(n)
incircle = 0.0
for throws = 0 to n
incircle = incircle + (rand()^2 + rand()^2 < 1)
next
return 4.0 * incircle / throws
end function</syntaxhighlight>
{{out}}
<pre> Number of throws Ratio (Pi) Error
100 2.970297 0.17129562389
1000 3.14085914086 0.00073351273
10000 3.13208679132 0.00950586227
100000 3.14428855711 -0.00269590352
1000000 3.14041685958 0.00117579401
10000000 3.14094968591 0.000643
100000000 3.14153 0.00006264501</pre>
 
==={{header|BBC BASIC}}===
<syntaxhighlight lang="bbcbasic"> PRINT FNmontecarlo(1000)
PRINT FNmontecarlo(10000)
PRINT FNmontecarlo(100000)
Line 430 ⟶ 458:
IF RND(1)^2 + RND(1)^2 < 1 n% += 1
NEXT
= 4 * n% / t%</langsyntaxhighlight>
{{out}}
<pre>
Line 439 ⟶ 467:
3.1412816
</pre>
 
==={{header|FreeBASIC}}===
<syntaxhighlight lang="freebasic">' version 23-10-2016
' compile with: fbc -s console
 
Randomize Timer 'seed the random function
 
Dim As Double x, y, pi, error_
Dim As UInteger m = 10, n, n_start, n_stop = m, p
 
Print
Print " Mumber of throws Ratio (Pi) Error"
Print
 
Do
For n = n_start To n_stop -1
x = Rnd
y = Rnd
If (x * x + y * y) <= 1 Then p = p +1
Next
Print Using " ############, "; m ;
pi = p * 4 / m
error_ = 3.141592653589793238462643383280 - pi
Print RTrim(Str(pi),"0");Tab(35); Using "##.#############"; error_
m = m * 10
n_start = n_stop
n_stop = m
Loop Until m > 1000000000 ' 1,000,000,000
 
 
' empty keyboard buffer
While Inkey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End</syntaxhighlight>
{{out}}
<pre> Mumber of throws Ratio (Pi) Error
 
10 3.2 -0.0584073464102
100 3.16 -0.0184073464102
1,000 3.048 0.0935926535898
10,000 3.1272 0.0143926535898
100,000 3.13672 0.0048726535898
1,000,000 3.14148 0.0001126535898
10,000,000 3.1417668 -0.0001741464102
100,000,000 3.14141 0.0001826535898
1,000,000,000 3.14169192 -0.0000992664102</pre>
 
 
==={{header|Liberty BASIC}}===
<syntaxhighlight lang="lb">
for pow = 2 to 6
n = 10^pow
print n, getPi(n)
next
 
end
 
function getPi(n)
incircle = 0
for throws=0 to n
scan
incircle = incircle + (rnd(1)^2+rnd(1)^2 < 1)
next
getPi = 4*incircle/throws
end function
</syntaxhighlight>
 
{{out}}
<pre>
100 2.89108911
1000 3.12887113
10000 3.13928607
100000 3.13864861
1000000 3.13945686
</pre>
 
==={{header|Locomotive Basic}}===
 
<syntaxhighlight lang="locobasic">10 mode 1:randomize time:defint a-z
20 input "How many samples";n
30 u=n/100+1
40 r=100
50 for i=1 to n
60 if i mod u=0 then locate 1,3:print using "##% done"; i/n*100
70 x=rnd*2*r-r
80 y=rnd*2*r-r
90 if sqr(x*x+y*y)<r then m=m+1
100 next
110 pi2!=4*m/n
120 locate 1,3
130 print m;"points in circle"
140 print "Computed value of pi:"pi2!
150 print "Difference to real value of pi: ";
160 print using "+#.##%"; (pi2!-pi)/pi*100</syntaxhighlight>
 
[[File:Monte Carlo, 200 points, Locomotive BASIC.png]]
[[File:Monte Carlo, 5000 points, Locomotive BASIC.png]]
 
==={{header|Run BASIC}}===
{{trans|Liberty BASIC}}
<syntaxhighlight lang="freebasic">for pow = 2 to 6
n = 10 ^ pow
print n; chr$(9); getPi(n)
next
end
 
function getPi(n)
incircle = 0
for throws = 0 to n
incircle = incircle + (rnd(1)^2 + rnd(1)^2 < 1)
next
getPi = 4 * incircle / throws
end function</syntaxhighlight>
{{out}}
<pre>100 3.12
1000 3.108
10000 3.1652
100000 3.14248
1000000 3.1435</pre>
 
==={{header|True BASIC}}===
{{trans|BASIC}}
<syntaxhighlight lang="qbasic">FUNCTION getpi(throws)
LET incircle = 0
FOR i = 1 to throws
!a square with a side of length 2 centered at 0 has
!x and y range of -1 to 1
LET randx = (rnd*2)-1 !range -1 to 1
LET randy = (rnd*2)-1 !range -1 to 1
!distance from (0,0) = sqrt((x-0)^2+(y-0)^2)
LET dist = sqr(randx^2+randy^2)
IF dist < 1 then !circle with diameter of 2 has radius of 1
LET incircle = incircle+1
END IF
NEXT i
LET getpi = 4*incircle/throws
END FUNCTION
 
CLEAR
PRINT getpi(10000)
PRINT getpi(100000)
PRINT getpi(1000000)
PRINT getpi(10000000)
END</syntaxhighlight>
{{out}}
<pre>3.1304
3.14324
3.141716
3.1416452</pre>
 
==={{header|PureBasic}}===
<syntaxhighlight lang="purebasic">OpenConsole()
Procedure.d MonteCarloPi(throws.d)
inCircle.d = 0
For i = 1 To throws.d
randX.d = (Random(2147483647)/2147483647)*2-1
randY.d = (Random(2147483647)/2147483647)*2-1
dist.d = Sqr(randX.d*randX.d + randY.d*randY.d)
If dist.d < 1
inCircle = inCircle + 1
EndIf
Next i
pi.d = (4 * inCircle / throws.d)
ProcedureReturn pi.d
EndProcedure
 
PrintN ("'built-in' #Pi = " + StrD(#PI,20))
PrintN ("MonteCarloPi(10000) = " + StrD(MonteCarloPi(10000),20))
PrintN ("MonteCarloPi(100000) = " + StrD(MonteCarloPi(100000),20))
PrintN ("MonteCarloPi(1000000) = " + StrD(MonteCarloPi(1000000),20))
PrintN ("MonteCarloPi(10000000) = " + StrD(MonteCarloPi(10000000),20))
 
PrintN("Press any key"): Repeat: Until Inkey() <> ""
</syntaxhighlight>
{{out}}
<pre>'built-in' #PI = 3.14159265358979310000
MonteCarloPi(10000) = 3.17119999999999980000
MonteCarloPi(100000) = 3.14395999999999990000
MonteCarloPi(1000000) = 3.14349599999999980000
MonteCarloPi(10000000) = 3.14127720000000020000
Press any key</pre>
 
=={{header|C}}==
<langsyntaxhighlight Clang="c">#include <stdio.h>
#include <stdlib.h>
#include <math.h>
Line 476 ⟶ 689:
printf("Pi is %f\n", pi(3e-4)); /* set to 1e-4 for some fun */
return 0;
}</langsyntaxhighlight>
 
=={{header|C sharp|C#}}==
<langsyntaxhighlight lang="csharp">using System;
 
class Program {
Line 502 ⟶ 715:
}
}
}</langsyntaxhighlight>
 
{{out}}
Line 514 ⟶ 727:
 
=={{header|C++}}==
<langsyntaxhighlight lang="cpp">
#include<iostream>
#include<math.h>
Line 536 ⟶ 749:
cout<<""<<4*double(hit)/double(imax)<<endl; } // Print out Pi number
}
</syntaxhighlight>
</lang>
 
=={{header|Clojure}}==
<langsyntaxhighlight lang="lisp">(defn calc-pi [iterations]
(loop [x (rand) y (rand) in 0 total 1]
(if (< total iterations)
Line 545 ⟶ 758:
(double (* (/ in total) 4)))))
 
(doseq [x (take 5 (iterate #(* 10 %) 10))] (println (str (format "% 8d" x) ": " (calc-pi x))))</langsyntaxhighlight>
 
{{out}}
Line 556 ⟶ 769:
</pre>
 
<langsyntaxhighlight lang="lisp">(defn experiment
[]
(if (<= (+ (Math/pow (rand) 2) (Math/pow (rand) 2)) 1) 1 0))
Line 565 ⟶ 778:
 
(pi-estimate 10000)
</syntaxhighlight>
</lang>
 
{{out}}
Line 573 ⟶ 786:
 
=={{header|Common Lisp}}==
<langsyntaxhighlight lang="lisp">(defun approximate-pi (n)
(/ (loop repeat n count (<= (abs (complex (random 1.0) (random 1.0))) 1.0)) n 0.25))
 
(dolist (n (loop repeat 5 for n = 1000 then (* n 10) collect n))
(format t "~%~8d -> ~f" n (approximate-pi n)))</langsyntaxhighlight>
 
{{out}}
Line 590 ⟶ 803:
=={{header|Crystal}}==
{{trans|Ruby}}
<langsyntaxhighlight lang="ruby">def approx_pi(throws)
times_inside = throws.times.count {Math.hypot(rand, rand) <= 1.0}
4.0 * times_inside / throws
Line 597 ⟶ 810:
[1000, 10_000, 100_000, 1_000_000, 10_000_000].each do |n|
puts "%8d samples: PI = %s" % [n, approx_pi(n)]
end</langsyntaxhighlight>
{{out}}
<pre> 1000 samples: PI = 3.1
Line 607 ⟶ 820:
 
=={{header|D}}==
<langsyntaxhighlight lang="d">import std.stdio, std.random, std.math;
 
double pi(in uint nthrows) /*nothrow*/ @safe /*@nogc*/ {
Line 620 ⟶ 833:
foreach (immutable p; 1 .. 8)
writefln("%10s: %07f", 10 ^^ p, pi(10 ^^ p));
}</langsyntaxhighlight>
{{out}}
<pre> 10: 3.200000
Line 641 ⟶ 854:
 
===More Functional Style===
<langsyntaxhighlight lang="d">void main() {
import std.stdio, std.random, std.math, std.algorithm, std.range;
 
Line 649 ⟶ 862:
foreach (immutable p; 1 .. 8)
writefln("%10s: %07f", 10 ^^ p, pi(10 ^^ p));
}</langsyntaxhighlight>
{{out}}
<pre> 10: 3.200000
Line 663 ⟶ 876:
From example at [https://www.dartlang.org/ Dart Official Website]
 
<langsyntaxhighlight lang="dart">
import 'dart:async';
import 'dart:html';
Line 714 ⟶ 927:
bool get isInsideUnitCircle => x * x + y * y <= 1;
}
</syntaxhighlight>
</lang>
{{out}}
The script give in reality an output formatted in HTML
<pre>π ≅ 3.14139</pre>
 
=={{header|Delphi}}==
{{works with|Delphi|6.0}}
{{libheader|SysUtils,StdCtrls}}
 
 
<syntaxhighlight lang="Delphi">
 
function MonteCarloPi(N: cardinal): double;
{Approximate Pi by seeing if points fall inside circle}
var I,InsideCnt: integer;
var X,Y: double;
begin
InsideCnt:=0;
for I:=1 to N do
begin
{Random X,Y = 0..1}
X:=Random;
Y:=Random;
{See if it falls in Unit Circle}
if X*X + Y*Y <= 1 then Inc(InsideCnt);
end;
{Because X and Y are squared, they only fall with 1/4 of the circle}
Result:=4 * InsideCnt / N;
end;
 
 
procedure ShowOneSimulation(Memo: TMemo; N: cardinal);
var MyPi: double;
begin
MyPi:=MonteCarloPi(N);
Memo.Lines.Add(Format('Samples: %15.0n Pi= %2.15f',[N+0.0,MyPi]));
end;
 
 
procedure ShowMonteCarloPi(Memo: TMemo);
begin
ShowOneSimulation(Memo,1000);
ShowOneSimulation(Memo,10000);
ShowOneSimulation(Memo,100000);
ShowOneSimulation(Memo,1000000);
ShowOneSimulation(Memo,10000000);
ShowOneSimulation(Memo,100000000);
end;
 
 
</syntaxhighlight>
{{out}}
<pre>
Samples: 1,000 Pi= 3.156000000000000
Samples: 10,000 Pi= 3.152000000000000
Samples: 100,000 Pi= 3.142920000000000
Samples: 1,000,000 Pi= 3.140864000000000
Samples: 10,000,000 Pi= 3.141990800000000
Samples: 100,000,000 Pi= 3.141426720000000
</pre>
 
 
=={{header|E}}==
Line 723 ⟶ 993:
This computes a single quadrant of the described square and circle; the effect should be the same since the other three are symmetric.
 
<langsyntaxhighlight lang="e">def pi(n) {
var inside := 0
for _ ? (entropy.nextFloat() ** 2 + entropy.nextFloat() ** 2 < 1) in 1..n {
Line 729 ⟶ 999:
}
return inside * 4 / n
}</langsyntaxhighlight>
 
Some sample runs:
Line 752 ⟶ 1,022:
 
=={{header|EasyLang}}==
<syntaxhighlight lang="text">
<lang>func mc n . .
func mc n .
for i range n
for xi = randomf1 to n
y x = randomf
if x * x + y * y <= 1randomf
hitif x * x += y * y < 1
. hit += 1
.
.
print 4.0 * hit / n
return 4.0 * hit / n
.
numfmt 4 0
fscale 4
callprint mc 10000
callprint mc 100000
callprint mc 1000000
callprint mc 10000000</lang>
</syntaxhighlight>
Output:
3.1292
Line 772 ⟶ 1,044:
3.1407
3.1413
 
=={{header|EDSAC order code}}==
Because real numbers on EDSAC were restricted to the interval [-1,1), this solution estimates pi/10 instead of pi. With 100,000 trials the program would have taken about 3.5 hours on the original EDSAC.
<syntaxhighlight lang="edsac">
[Monte Carlo solution for Rosetta Code.]
[EDSAC program, Initial Orders 2.]
 
[Arrange the storage]
T45K P56F [H parameter: library s/r P1 to print real number]
T46K P78F [N parameter: library s/r P7 to print integer]
T47K P210F [M parameter: main routine]
T48K P114F [& (delta) parameter: library s/r C6 (division)]
T49K P150F [L parameter: library subroutine R4 to read data]
T51K P172F [G parameter: generator for pseudo-random numbers]
 
[Library subroutine M3, runs at load time and is then overwritten.
Prints header; here the header sets teleprinter to figures.]
PFGKIFAFRDLFUFOFE@A6FG@E8FEZPF
*!!!!TRIALS!!EST!PI#X10@&..
PK [after header, blank tape and PK (WWG, 1951, page 91)]
 
[================ Main routine ====================]
E25K TM GK
[Variables]
[0] PF PF [target count: print result when count = target]
[2] PF PF [count of points]
[4] PF PF [count of hits (point inside circle)]
[6] PF PF [x-coordinate - 1/2]
[Constants]
T8#Z PF T10#Z PF [clear sandwich bits in 35-bit constants]
T8Z [resume normal loading]
[8] PD PF [35-bit constant 1]
[10] L1229F Y819F [35-bit constant 2/5 (near enough)]
[12] IF [1/2]
[13] RF [1/4]
[14] #F [figures shift]
[15] MF [dot (decimal point) in figures mode]
[16] @F [carriage return]
[17] &F [line feed]
[18] !F [space]
 
[Enter with acc = 0]
[19] A19@ GL [read seed for LCG into 0D]
AD T4D [pass seed to LCG in 4D]
[23] A23@ GG [initialize LCG]
T2#@ T4#@ [zero trials and hits]
[Outer loop: round target counts]
[27] TF [clear acc]
[28] A28@ GL [read next target count into 0D]
SD [acc := -target]
E85@ [exit if target = 0]
T#@ [store negated target]
[Inner loop : round points in the square]
[33] TF T4D [pass LCG range = 0 to return random real in [0,1)]
[35] A35@ G1G [call LCG, 0D := random x]
AD S12@ T6#@ [store x - 1/2 over next call]
T4D
[41] A41@ G1G [call LCG, 0D := random y]
AD S12@ TD [store y - 1/2]
H6#@ V6#@ [acc := (x - 1/2)^2]
HD VD [acc := acc := (x - 1/2)^2 + (y - 1/2)^2]
S13@ [test for point inside circle, i.e. acc < 1/4]
E56@ [skip if not]
TF A4#@ A8#@ T4#@ [inc number of hits]
[56] TF A2#@ A8#@ U2#@ [inc number of trials]
A#@ [add negated target]
G33@ [if not reached target, loop back]
A2#@ TD [pass number of trials to print s/r]
[64] A64@ GN [print number of trials]
A4#@ TD A2#@ T4D [pass hits and trials to division s/r]
[70] A70@ G& [0D := hits/trials, estimated value of pi/4]
HD V10#@ TD [times 2/5; pass estimated pi/10 to print s/r]
O18@ O18@ O8@ O15@ [print ' 0.']
[79] A79@ GH P5F [print estimated pi/10 to 5 decimals]
O16@ O17@ [print CR, LF]
E27@ [loop back for new target]
[85] O14@ [exit: print dummy character to flush printer buffer]
ZF [halt program]
 
[==================== Generator for pseudo-random numbers ===========]
[Linear congruential generator, same algorithm as Delphi 7 LCG.
38 locations]
E25K TG
GK G10@ G15@ T2#Z PF T2Z I514D P257F T4#Z PF T4Z PD PF T6#Z PF T6Z PF RF A6#@ S4#@ T6#@ E25F E8Z PF T8Z PF PF A3F T14@ A4D T8#@ ZF A3F T37@ H2#@ V8#@ L512F L512F L1024F A4#@ T8#@ H6#@ C8#@ T8#@ S4D G32@ TD A8#@ E35@ H4D TD V8#@ L1F TD ZF
 
[==================== LIBRARY SUBROUTINES ============================]
[D6: Division, accurate, fast.
36 locations, workspace 6D and 8D.
0D := 0D/4D, where 4D <> 0, -1.]
E25K T& GK
GKA3FT34@S4DE13@T4DSDTDE2@T4DADLDTDA4DLDE8@RDU4DLDA35@
T6DE25@U8DN8DA6DT6DH6DS6DN4DA4DYFG21@SDVDTDEFW1526D
 
[R4: Input of one signed integer at runtime.
22 storage locations; working positions 4, 5, and 6.]
E25K TL
GKA3FT21@T4DH6@E11@P5DJFT6FVDL4FA4DTDI4FA4FS5@G7@S5@G20@SDTDT6FEF
 
[P1: Prints non-negative fraction in 0D, without '0.']
E25K TH
GKA18@U17@S20@T5@H19@PFT5@VDUFOFFFSFL4FTDA5@A2FG6@EFU3FJFM1F
 
[P7, prints long strictly positive integer;
10 characters, right justified, padded left with spaces.
Even address; 35 storage locations; working position 4D.]
E25K TN
GKA3FT26@H28#@NDYFLDT4DS27@TFH8@S8@T1FV4DAFG31@SFLDUFOFFFSF
L4FT4DA1FA27@G11@XFT28#ZPFT27ZP1024FP610D@524D!FO30@SFL8FE22@
 
[===================================================================]
[The following, without the comments and white space, might have
been input from a separate tape.]
E25K TM GK
E19Z [define entry point]
PF [acc = 0 on entry]
[Integers supplied by user: (1) seed for LCG; (2) list of numbers of trials
for which to print result; increasing order, terminated by 0.
To be read by library subroutine R4; sign comes after value.]
987654321+100+1000+10000+100000+0+
</syntaxhighlight>
{{out}}
<pre>
TRIALS EST PI/10
100 0.32400
1000 0.31319
10000 0.31371
100000 0.31410
</pre>
 
=={{header|Elixir}}==
<langsyntaxhighlight lang="elixir">defmodule MonteCarlo do
def pi(n) do
count = Enum.count(1..n, fn _ ->
Line 787 ⟶ 1,187:
Enum.each([1000, 10000, 100000, 1000000, 10000000], fn n ->
:io.format "~8w samples: PI = ~f~n", [n, MonteCarlo.pi(n)]
end)</langsyntaxhighlight>
 
{{out}}
Line 800 ⟶ 1,200:
=={{header|Erlang}}==
===With inline test===
<syntaxhighlight lang="erlang">
<lang ERLANG>
-module(monte).
-export([main/1]).
Line 818 ⟶ 1,218:
 
main(N) -> io:format("PI: ~w~n", [ monte(N) ]).
</syntaxhighlight>
</lang>
{{out}}
<pre>
Line 830 ⟶ 1,230:
</pre>
===With test in a function===
<syntaxhighlight lang="erlang">
<lang ERLANG>
-module(monte2).
-export([main/1]).
Line 851 ⟶ 1,251:
 
main(N) -> io:format("PI: ~w~n", [ monte(N) ]).
</syntaxhighlight>
</lang>
{{out}}
<pre>Xcoord
Line 865 ⟶ 1,265:
 
=={{header|ERRE}}==
<syntaxhighlight lang="erre">
<lang ERRE>
PROGRAM RANDOM_PI
 
Line 889 ⟶ 1,289:
MONTECARLO(1000000->RES) PRINT(RES)
MONTECARLO(10000000->RES) PRINT(RES)
END PROGRAM</langsyntaxhighlight>
{{out}}
<pre>
Line 900 ⟶ 1,300:
 
=={{header|Euler Math Toolbox}}==
<syntaxhighlight lang="euler math toolbox">
<lang Euler Math Toolbox>
>function map MonteCarloPI (n,plot=false) ...
$ X:=random(1,n);
Line 915 ⟶ 1,315:
3.14159265359
>MonteCarloPI(10000,true):
</syntaxhighlight>
</lang>
 
[[File:Test.png]]
Line 922 ⟶ 1,322:
There is some support and test expressions.
 
<langsyntaxhighlight lang="fsharp">
let print x = printfn "%A" x
 
Line 942 ⟶ 1,342:
MonteCarloPiGreco 10000 |> print
MonteCarloPiGreco 100000 |> print
</syntaxhighlight>
</lang>
{{out}}
<pre>
Line 953 ⟶ 1,353:
Since Factor lets the user choose the range of the random generator, we use 2^32.
 
<langsyntaxhighlight lang="factor">USING: kernel math math.functions random sequences ;
 
: limit ( -- n ) 2 32 ^ ; inline
: in-circle ( x y -- ? ) limit [ sq ] tri@ [ + ] [ <= ] bi* ;
: rand ( -- r ) limit random ;
: pi ( n -- pi ) [ [ drop rand rand in-circle ] count ] keep / 4 * >float ;</langsyntaxhighlight>
 
Example use:
 
<langsyntaxhighlight lang="factor">10000 pi .
3.1412</langsyntaxhighlight>
 
=={{header|Fantom}}==
 
<langsyntaxhighlight lang="fantom">
class MontyCarlo
{
Line 991 ⟶ 1,391:
}
}
</syntaxhighlight>
</lang>
 
{{out}}
Line 1,024 ⟶ 1,424:
=={{header|Fortran}}==
{{works with|Fortran|90 and later}}
<langsyntaxhighlight lang="fortran">MODULE Simulation
IMPLICIT NONE
Line 1,058 ⟶ 1,458:
END DO
END PROGRAM MONTE_CARLO</langsyntaxhighlight>
 
{{out}}
Line 1,069 ⟶ 1,469:
</pre>
{{works with|Fortran|2008 and later}}
<langsyntaxhighlight lang="fortran">
program mc
integer :: n,i
Line 1,086 ⟶ 1,486:
pi = 4.d0 * dble( count( hypot(x(1,:),x(2,:)) <= 1.d0 ) ) / n
end function
</syntaxhighlight>
</lang>
 
=={{header|FreeBASIC}}==
<lang freebasic>' version 23-10-2016
' compile with: fbc -s console
 
Randomize Timer 'seed the random function
 
Dim As Double x, y, pi, error_
Dim As UInteger m = 10, n, n_start, n_stop = m, p
 
Print
Print " Mumber of throws Ratio (Pi) Error"
Print
 
Do
For n = n_start To n_stop -1
x = Rnd
y = Rnd
If (x * x + y * y) <= 1 Then p = p +1
Next
Print Using " ############, "; m ;
pi = p * 4 / m
error_ = 3.141592653589793238462643383280 - pi
Print RTrim(Str(pi),"0");Tab(35); Using "##.#############"; error_
m = m * 10
n_start = n_stop
n_stop = m
Loop Until m > 1000000000 ' 1,000,000,000
 
 
' empty keyboard buffer
While Inkey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End</lang>
{{out}}
<pre> Mumber of throws Ratio (Pi) Error
 
10 3.2 -0.0584073464102
100 3.16 -0.0184073464102
1,000 3.048 0.0935926535898
10,000 3.1272 0.0143926535898
100,000 3.13672 0.0048726535898
1,000,000 3.14148 0.0001126535898
10,000,000 3.1417668 -0.0001741464102
100,000,000 3.14141 0.0001826535898
1,000,000,000 3.14169192 -0.0000992664102</pre>
 
=={{header|Futhark}}==
Line 1,139 ⟶ 1,492:
Since Futhark is a pure language, random numbers are implemented using Sobol sequences.
 
<syntaxhighlight lang="futhark">
<lang Futhark>
import "futlib/math"
 
Line 1,188 ⟶ 1,541:
let inside = reduce (+) 0 bs
in 4.0f32*f32(inside)/f32(n)
</syntaxhighlight>
</lang>
 
=={{header|Go}}==
'''Using standard library math/rand:'''
<langsyntaxhighlight lang="go">package main
 
import (
Line 1,224 ⟶ 1,577:
fmt.Println(getPi(10000000))
fmt.Println(getPi(100000000))
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,236 ⟶ 1,589:
 
For very careful Monte Carlo studies, you might consider the subrepository rand library. The random number generator there has some advantages such as better known statistical properties and better use of memory.
<langsyntaxhighlight lang="go">package main
 
import (
Line 1,269 ⟶ 1,622:
fmt.Println(getPi(10000000))
fmt.Println(getPi(100000000))
}</langsyntaxhighlight>
 
=={{header|Haskell}}==
<langsyntaxhighlight lang="haskell">import Control.Monad
import System.Random
 
Line 1,284 ⟶ 1,637:
let dist :: Double
dist = sqrt (rand_x * rand_x + rand_y * rand_y)
return (if dist < 1 then 1 else 0)</langsyntaxhighlight>
{{Out}}
<pre>Example:
Line 1,296 ⟶ 1,649:
Or, using foldM, and dropping sqrt:
 
<langsyntaxhighlight lang="haskell">import Control.Monad (foldM, (>=>))
import System.Random (randomRIO)
import Data.Functor ((<&>))
Line 1,318 ⟶ 1,671:
mapM_
(monteCarloPi >=> print)
[1000, 10000, 100000, 1000000]</langsyntaxhighlight>
{{Out}}
For example:
Line 1,327 ⟶ 1,680:
 
=={{header|HicEst}}==
<langsyntaxhighlight HicEstlang="hicest">FUNCTION Pi(samples)
inside = 0
DO i = 1, samples
Line 1,338 ⟶ 1,691:
WRITE(ClipBoard) Pi(1E5) ! 3.14204
WRITE(ClipBoard) Pi(1E6) ! 3.141672
WRITE(ClipBoard) Pi(1E7) ! 3.1412856</langsyntaxhighlight>
 
=={{header|Icon}} and {{header|Unicon}}==
<langsyntaxhighlight Iconlang="icon">procedure main()
every t := 10 ^ ( 5 to 9 ) do
printf("Rounds=%d Pi ~ %r\n",t,getPi(t))
Line 1,354 ⟶ 1,707:
incircle +:= 1
return 4 * incircle / rounds
end</langsyntaxhighlight>
 
{{libheader|Icon Programming Library}}
Line 1,368 ⟶ 1,721:
=={{header|J}}==
'''Explicit Solution:'''
<langsyntaxhighlight lang="j">piMC=: monad define "0
4* y%~ +/ 1>: %: +/ *: <: +: (2,y) ?@$ 0
)</langsyntaxhighlight>
 
'''Tacit Solution:'''
<langsyntaxhighlight lang="j">piMCt=: (0.25&* %~ +/@(1 >: [: +/&.:*: _1 2 p. 0 ?@$~ 2&,))"0</langsyntaxhighlight>
 
'''Examples:'''
<langsyntaxhighlight lang="j"> piMC 1e6
3.1426
piMC 10^i.7
4 2.8 3.24 3.168 3.1432 3.14256 3.14014</langsyntaxhighlight>
 
'''Alternative Tacit Solution:'''
<syntaxhighlight lang="j">pimct=. (4 * +/ % #)@:(1 >: |)@:(? j. ?)@:($&0)"0
(,. pimct) 10 ^ 3 + i.6
1000 3.168
10000 3.122
100000 3.13596
1e6 3.1428
1e7 3.14158
1e8 3.14154</syntaxhighlight>
 
=={{header|Java}}==
<langsyntaxhighlight lang="java">public class MC {
public static void main(String[] args) {
System.out.println(getPi(10000));
Line 1,407 ⟶ 1,771:
return 4.0 * inCircle / numThrows;
}
}</langsyntaxhighlight>
{{out}}
3.1396
Line 1,415 ⟶ 1,779:
3.14168604
{{works with|Java|8+}}
<langsyntaxhighlight lang="java">package montecarlo;
 
import java.util.stream.IntStream;
Line 1,458 ⟶ 1,822:
return (4.0 * inCircle) / numThrows;
}
}</langsyntaxhighlight>
{{out}}
3.1556
Line 1,468 ⟶ 1,832:
=={{header|JavaScript}}==
===ES5===
<langsyntaxhighlight JavaScriptlang="javascript">function mcpi(n) {
var x, y, m = 0;
 
Line 1,487 ⟶ 1,851:
console.log(mcpi(100000));
console.log(mcpi(1000000));
console.log(mcpi(10000000));</langsyntaxhighlight>
<pre>3.168
3.1396
Line 1,496 ⟶ 1,860:
 
===ES6===
<langsyntaxhighlight JavaScriptlang="javascript">(() => {
"use strict";
 
Line 1,506 ⟶ 1,870:
const [x, y] = [rnd(), rnd()];
 
return (x ** x2) + (y ** y2) < 1 ? (
1 + a
) : a;
Line 1,534 ⟶ 1,898:
);
});
})();</langsyntaxhighlight>
 
{{Out}} For example:
<pre>1000 samples: 3.064
Line 1,550 ⟶ 1,913:
jq does not have a built-in PRNG so we will use /dev/urandom
as a source of entropy by invoking jq as follows:
<langsyntaxhighlight lang="sh"># In case gojq is used, trim leading 0s:
function prng {
cat /dev/urandom | tr -cd '0-9' | fold -w 10 | sed 's/^0*\(.*\)*\(.\)*$/\1\2/'
}
 
prng | jq -nMr -f program.jq</langsyntaxhighlight>
 
'''program.jq'''
<langsyntaxhighlight lang="jq">def lpad($len): tostring | ($len - length) as $l | (" " * $l)[:$l] + .;
 
def percent: "\(100000 * . | round / 1000)%";
Line 1,581 ⟶ 1,944:
| ($p | mcPi) as $mcpi
| ((($pi - $mcpi)|length) / $pi) as $error
| "\($p|lpad(10)) \($mcpi|lpad(10)) \($error|percent|lpad(6))" )</langsyntaxhighlight>
{{out}}
<pre>
Line 1,596 ⟶ 1,959:
=={{header|Jsish}}==
From Javascript ES5 entry, with PRNG seeded during unit testing for reproducibility.
<langsyntaxhighlight lang="javascript">/* Monte Carlo methods, in Jsish */
function mcpi(n) {
var x, y, m = 0;
Line 1,627 ⟶ 1,990:
mcpi(1000000) ==> 3.142124
=!EXPECTEND!=
*/</langsyntaxhighlight>
 
{{out}}
Line 1,634 ⟶ 1,997:
 
=={{header|Julia}}==
<langsyntaxhighlight lang="julia">using Printf
 
function monteπ(n)
Line 1,643 ⟶ 2,006:
for n in 10 .^ (3:8)
p = monteπ(n)
println("$(lpad(n, 9)): π ≈ $(lpad(p, 10)), pct.err = ", @sprintf("%2.5f%%", 100 * abs(p - π) / π))
end</langsyntaxhighlight>
 
{{out}}
<pre> 1000: π ≈ 3.224, pct.err = 0.02623%
10000: π ≈ 3.1336, pct.err = 0.00254254%
100000: π ≈ 3.13468, pct.err = 0.00220220%
1000000: π ≈ 3.14156, pct.err = 0.00001001%
10000000: π ≈ 3.1412348, pct.err = 0.00011011%
100000000: π ≈ 3.14123216, pct.err = 0.00011011%</pre>
 
=={{header|K}}==
<langsyntaxhighlight Klang="k"> sim:{4*(+/{~1<+/(2_draw 0)^2}'!x)%x}
 
sim 10000
Line 1,661 ⟶ 2,024:
 
sim'10^!8
4 2.8 3.4 3.072 3.1212 3.14104 3.14366 3.1413</langsyntaxhighlight>
 
=={{header|Kotlin}}==
<langsyntaxhighlight lang="scala">// version 1.1.0
 
fun mcPi(n: Int): Double {
Line 1,686 ⟶ 2,049:
n *= 10
}
}</langsyntaxhighlight>
Sample output:
{{out}}
Line 1,699 ⟶ 2,062:
100000000 -> 3.14160244 -> 0.0003
</pre>
 
=={{header|Liberty BASIC}}==
<lang lb>
for pow = 2 to 6
n = 10^pow
print n, getPi(n)
next
 
end
 
function getPi(n)
incircle = 0
for throws=0 to n
scan
incircle = incircle + (rnd(1)^2+rnd(1)^2 < 1)
next
getPi = 4*incircle/throws
end function
</lang>
 
{{out}}
<pre>
100 2.89108911
1000 3.12887113
10000 3.13928607
100000 3.13864861
1000000 3.13945686
</pre>
 
=={{header|Locomotive Basic}}==
 
<lang locobasic>10 mode 1:randomize time:defint a-z
20 input "How many samples";n
30 u=n/100+1
40 r=100
50 for i=1 to n
60 if i mod u=0 then locate 1,3:print using "##% done"; i/n*100
70 x=rnd*2*r-r
80 y=rnd*2*r-r
90 if sqr(x*x+y*y)<r then m=m+1
100 next
110 pi2!=4*m/n
120 locate 1,3
130 print m;"points in circle"
140 print "Computed value of pi:"pi2!
150 print "Difference to real value of pi: ";
160 print using "+#.##%"; (pi2!-pi)/pi*100</lang>
 
[[File:Monte Carlo, 200 points, Locomotive BASIC.png]]
[[File:Monte Carlo, 5000 points, Locomotive BASIC.png]]
 
=={{header|Logo}}==
<langsyntaxhighlight lang="logo">
to square :n
output :n * :n
Line 1,769 ⟶ 2,081:
show sim 100000 10000 ; 3.145
show sim 1000000 10000 ; 3.140828
</syntaxhighlight>
</lang>
 
=={{header|LSL}}==
To test it yourself; rez a box on the ground, and add the following as a New Script.
(Be prepared to wait... LSL can be slow, but the Servers are typically running thousands of scripts in parallel so what do you expect?)
<langsyntaxhighlight LSLlang="lsl">integer iMIN_SAMPLE_POWER = 0;
integer iMAX_SAMPLE_POWER = 6;
default {
Line 1,795 ⟶ 2,107:
llOwnerSay("Done.");
}
}</langsyntaxhighlight>
{{out}}
<pre>Estimating Pi (3.141593)
Line 1,808 ⟶ 2,120:
 
=={{header|Lua}}==
<langsyntaxhighlight lang="lua">function MonteCarlo ( n_throws )
math.randomseed( os.time() )
 
Line 1,824 ⟶ 2,136:
print( MonteCarlo( 100000 ) )
print( MonteCarlo( 1000000 ) )
print( MonteCarlo( 10000000 ) )</langsyntaxhighlight>
{{out}}
<pre>3.1436
Line 1,833 ⟶ 2,145:
=={{header|Mathematica}}/{{header|Wolfram Language}}==
We define a function with variable sample size:
<langsyntaxhighlight Mathematicalang="mathematica">MonteCarloPi[samplesize_Integer] := N[4Mean[If[# > 1, 0, 1] & /@ Norm /@ RandomReal[1, {samplesize, 2}]]]</langsyntaxhighlight>
Example (samplesize=10,100,1000,....10000000):
<langsyntaxhighlight Mathematicalang="mathematica">{#, MonteCarloPi[#]} & /@ (10^Range[1, 7]) // Grid</langsyntaxhighlight>
gives back:
<pre>10 3.2
Line 1,845 ⟶ 2,157:
10000000 3.14134</pre>
 
<langsyntaxhighlight Mathematicalang="mathematica">monteCarloPi = 4. Mean[UnitStep[1 - Total[RandomReal[1, {2, #}]^2]]] &;
monteCarloPi /@ (10^Range@6)</langsyntaxhighlight>
 
A less elegant way to solve the problem, is to imagine a (well-trained) monkey, throwing a number of darts at a dartboard.
Line 1,853 ⟶ 2,165:
 
We create a function ''MonkeyDartsPi'', which can take a variable number of throws as input:
<langsyntaxhighlight Wolframlang="wolfram Languagelanguage">MonkeyDartsPi[numberOfThrows_] := (
xyCoordinates = RandomReal[{0, 1}, {numberOfThrows, 2}];
InsideCircle = Length[Select[Total[xyCoordinates^2, {2}],#<=1&]] ;
4*N[InsideCircle / Length[xyCoordinates],1+Log10[numberOfThrows]])</langsyntaxhighlight>
 
We do several runs with a larger number of throws each time, increasing by powers of 10.
<langsyntaxhighlight Wolframlang="wolfram Languagelanguage">Grid[Table[{n, MonkeyDartsPi[n]}, {n, 10^Range[7]} ], Alignment -> Left]</langsyntaxhighlight>
 
We see that as the number of throws increases, we get closer to the value of Pi:
Line 1,876 ⟶ 2,188:
 
Minimally Vectorized:
<langsyntaxhighlight MATLABlang="matlab">function piEstimate = monteCarloPi(numDarts)
 
%The square has a sides of length 2, which means the circle has radius
Line 1,892 ⟶ 2,204:
 
end
</syntaxhighlight>
</lang>
 
Completely Vectorized:
<langsyntaxhighlight MATLABlang="matlab">function piEstimate = monteCarloPi(numDarts)
piEstimate = 4*sum( sum(rand(numDarts,2).^2,2) <= 1 )/numDarts;
 
end</langsyntaxhighlight>
 
{{out}}
<langsyntaxhighlight MATLABlang="matlab">>> monteCarloPi(7000000)
 
ans =
 
3.141512000000000</langsyntaxhighlight>
 
=={{header|Maxima}}==
<langsyntaxhighlight Maximalang="maxima">load("distrib");
approx_pi(n):= block(
[x: random_continuous_uniform(0, 1, n),
Line 1,918 ⟶ 2,230:
4*cin/n);
float(approx_pi(100));</langsyntaxhighlight>
 
=={{header|MAXScript}}==
Line 1,939 ⟶ 2,251:
 
=={{header|МК-61/52}}==
<syntaxhighlight lang="text">П0 П1 0 П4 СЧ x^2 ^ СЧ x^2 +
1 - x<0 15 КИП4 L0 04 ИП4 4 *
ИП1 / С/П</langsyntaxhighlight>
 
''Example:'' for n = ''1000'' the output is ''3.152''.
 
=={{header|Nim}}==
<langsyntaxhighlight lang="nim">import math, random
 
randomize()
Line 1,958 ⟶ 2,270:
for n in [10e4, 10e6, 10e7, 10e8]:
echo pi(n)</langsyntaxhighlight>
{{out}}
<pre>3.15336
Line 1,966 ⟶ 2,278:
 
=={{header|OCaml}}==
<langsyntaxhighlight lang="ocaml">let get_pi throws =
let rec helper i count =
if i = throws then count
Line 1,977 ⟶ 2,289:
else
helper (i+1) count
in float (4 * helper 0 0) /. float throws</langsyntaxhighlight>
Example:
# get_pi 10000;;
Line 1,992 ⟶ 2,304:
=={{header|Octave}}==
 
<langsyntaxhighlight lang="octave">function p = montepi(samples)
in_circle = 0;
for samp = 1:samples
Line 2,007 ⟶ 2,319:
disp(montepi(l));
l *= 10;
endwhile</langsyntaxhighlight>
 
Since it runs slow, I've stopped it at the second iteration, obtaining:
Line 2,015 ⟶ 2,327:
=== Much faster implementation ===
 
<langsyntaxhighlight lang="octave">
function result = montepi(n)
result = sum(rand(1,n).^2+rand(1,n).^2<1)/n*4;
endfunction
</syntaxhighlight>
</lang>
 
=={{header|PARI/GP}}==
<langsyntaxhighlight lang="parigp">MonteCarloPi(tests)=4.*sum(i=1,tests,norml2([random(1.),random(1.)])<1)/tests;</langsyntaxhighlight>
A hundred million tests (about a minute) yielded 3.14149000, slightly more precise (and round!) than would have been expected. A million gave 3.14162000 and a thousand 3.14800000.
 
=={{header|Pascal}}==
{{libheader|Math}}
<langsyntaxhighlight lang="pascal">Program MonteCarlo(output);
 
uses
Line 2,056 ⟶ 2,368:
writeln (10**i, ' samples give ', MC_Pi(i):7:5, ' as pi.');
end.
</syntaxhighlight>
</lang>
{{out}}
<pre>:> ./MonteCarlo
Line 2,067 ⟶ 2,379:
 
=={{header|Perl}}==
<langsyntaxhighlight lang="perl">sub pi {
my $nthrows = shift;
my $inside = 0;
Line 2,080 ⟶ 2,392:
}
 
printf "%9d: %07f\n", $_, pi($_) for 10**4, 10**6;</langsyntaxhighlight>
{{out}}
<pre>
Line 2,088 ⟶ 2,400:
 
=={{header|Phix}}==
<!--<langsyntaxhighlight Phixlang="phix">(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">N</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">100</span>
Line 2,101 ⟶ 2,413:
<span style="color: #000000;">N</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">10</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<!--</langsyntaxhighlight>-->
{{out}}
<pre>
Line 2,113 ⟶ 2,425:
 
=={{header|PHP}}==
<syntaxhighlight lang="php"><?
<lang PHP><?
$loop = 1000000; # loop to 1,000,000
$count = 0;
Line 2,122 ⟶ 2,434:
}
echo "loop=".number_format($loop).", count=".number_format($count).", pi=".($count/$loop*4);
?></langsyntaxhighlight>
{{out}}
<pre>loop=1,000,000, count=785,462, pi=3.141848</pre>
 
=={{header|Picat}}==
Some general Monte Carlo simulators. <code>N</code> is the number of runs, <code>F</code> is the simulation function.
===Using while loop===
<syntaxhighlight lang="text">
sim1(N, F) = C =>
C = 0,
I = 0,
while (I <= N)
C := C + apply(F),
I := I + 1
end.</syntaxhighlight>
 
===List comprehension===
This is simpler, but slightly slower than using <code>while</code> loop.
<syntaxhighlight lang="picat">sim2(N, F) = sum([apply(F) : _I in 1..N]).</syntaxhighlight>
 
===Recursion===
<syntaxhighlight lang="picat">sim_rec(N,F) = S =>
sim_rec(N,N,F,0,S).
sim_rec(0,_N,_F,S,S).
sim_rec(C,N,F,S0,S) :-
S1 = S0 + apply(F),
sim_rec(C-1,N,F,S1,S).</syntaxhighlight>
 
===Test===
Of the three different MC simulators, <code>sim_rec/2</code> (using recursion) is slightly faster than the other two (<code>sim1/2</code> and <code>sim2/2</code>) which have about the same speed.
<syntaxhighlight lang="picat">go =>
foreach(N in 0..7)
sim_pi(10**N)
end,
nl.
 
% The specific pi simulation
sim_pi(N) =>
Inside = sim(N,pi_f),
MyPi = 4.0*Inside/N,
Pi = math.pi,
println([n=N, myPi=MyPi, diff=Pi-MyPi]).
 
% The simulation function:
% returns 1 if success, 0 otherwise
pi_f() = cond(frand()**2 + frand()**2 <= 1, 1, 0).</syntaxhighlight>
 
{{out}}
<pre>[n = 1,myPi = 4.0,diff = -0.858407346410207]
[n = 10,myPi = 3.2,diff = -0.058407346410207]
[n = 100,myPi = 3.12,diff = 0.021592653589793]
[n = 1000,myPi = 3.152,diff = -0.010407346410207]
[n = 10000,myPi = 3.1672,diff = -0.025607346410207]
[n = 100000,myPi = 3.13888,diff = 0.002712653589793]
[n = 1000000,myPi = 3.14192,diff = -0.000327346410207]
[n = 10000000,myPi = 3.1408988,diff = 0.000693853589793]</pre>
 
=={{header|PicoLisp}}==
<langsyntaxhighlight PicoLisplang="picolisp">(de carloPi (Scl)
(let (Dim (** 10 Scl) Dim2 (* Dim Dim) Pi 0)
(do (* 4 Dim)
Line 2,136 ⟶ 2,501:
 
(for N 6
(prinl (carloPi N)) )</langsyntaxhighlight>
{{out}}
<pre>3.4
Line 2,147 ⟶ 2,512:
=={{header|PowerShell}}==
{{works with|PowerShell|2}}
<langsyntaxhighlight lang="powershell">function Get-Pi ($Iterations = 10000) {
$InCircle = 0
for ($i = 0; $i -lt $Iterations; $i++) {
Line 2,163 ⟶ 2,528:
| Add-Member -PassThru NoteProperty Pi $Pi `
| Add-Member -PassThru NoteProperty "% Difference" $Diff
}</langsyntaxhighlight>
This returns a custom object with appropriate properties which automatically enables a nice tabular display:
<pre>PS Home:\> 10,100,1e3,1e4,1e5,1e6 | ForEach-Object { Get-Pi $_ }
Line 2,175 ⟶ 2,540:
100000 3,14712 0,1759409006731298209938938800
1000000 3,141364 0,0072782698142600895432451100</pre>
 
=={{header|PureBasic}}==
<lang PureBasic>OpenConsole()
Procedure.d MonteCarloPi(throws.d)
inCircle.d = 0
For i = 1 To throws.d
randX.d = (Random(2147483647)/2147483647)*2-1
randY.d = (Random(2147483647)/2147483647)*2-1
dist.d = Sqr(randX.d*randX.d + randY.d*randY.d)
If dist.d < 1
inCircle = inCircle + 1
EndIf
Next i
pi.d = (4 * inCircle / throws.d)
ProcedureReturn pi.d
EndProcedure
 
PrintN ("'built-in' #Pi = " + StrD(#PI,20))
PrintN ("MonteCarloPi(10000) = " + StrD(MonteCarloPi(10000),20))
PrintN ("MonteCarloPi(100000) = " + StrD(MonteCarloPi(100000),20))
PrintN ("MonteCarloPi(1000000) = " + StrD(MonteCarloPi(1000000),20))
PrintN ("MonteCarloPi(10000000) = " + StrD(MonteCarloPi(10000000),20))
 
PrintN("Press any key"): Repeat: Until Inkey() <> ""
</lang>
{{out}}
<pre>'built-in' #PI = 3.14159265358979310000
MonteCarloPi(10000) = 3.17119999999999980000
MonteCarloPi(100000) = 3.14395999999999990000
MonteCarloPi(1000000) = 3.14349599999999980000
MonteCarloPi(10000000) = 3.14127720000000020000
Press any key</pre>
 
=={{header|Python}}==
Line 2,216 ⟶ 2,547:
 
One use of the "sum" function is to count how many times something is true (because True = 1, False = 0):
<langsyntaxhighlight lang="python">>>> import random, math
>>> throws = 1000
>>> 4.0 * sum(math.hypot(*[random.random()*2-1
Line 2,231 ⟶ 2,562:
for q in [0,1]]) < 1
for p in xrange(throws)) / float(throws)
3.1415666400000002</langsyntaxhighlight>
 
===As a program using a function===
<langsyntaxhighlight lang="python">
from random import random
from math import hypot
Line 2,252 ⟶ 2,583:
for n in [10**4, 10**6, 10**7, 10**8]:
print "%9d: %07f" % (n, pi(n))
</syntaxhighlight>
</lang>
 
===Faster implementation using Numpy===
<langsyntaxhighlight lang="python">
import numpy as np
 
n = input('Number of samples: ')
print np.sum(np.random.rand(n)**2+np.random.rand(n)**2<1)/float(n)*4
</syntaxhighlight>
</lang>
 
=={{header|Quackery}}==
Line 2,266 ⟶ 2,597:
{{trans|Forth}}
 
<langsyntaxhighlight Quackerylang="quackery"> [ $ "bigrat.qky" loadfile ] now!
 
[ [ 64 bit ] constant
Line 2,280 ⟶ 2,611:
swap 20 point$ echo$ cr ] is trials ( n --> )
 
' [ 10 100 1000 10000 100000 1000000 ] witheach trials</langsyntaxhighlight>
 
{{out}}
Line 2,292 ⟶ 2,623:
 
=={{header|R}}==
<langsyntaxhighlight Rlang="r"># nice but not suitable for big samples!
monteCarloPi <- function(samples) {
x <- runif(samples, -1, 1) # for big samples, you need a lot of memory!
Line 2,318 ⟶ 2,649:
print(monteCarloPi(1e4))
print(monteCarloPi(1e5))
print(monteCarlo2Pi(1e7))</langsyntaxhighlight>
 
=={{header|Racket}}==
<langsyntaxhighlight lang="racket">#lang racket
 
(define (in-unit-circle? x y) (<= (sqrt (+ (sqr x) (sqr y))) 1))
Line 2,352 ⟶ 2,683:
;; to see how it looks as a decimal we can exact->inexact it
(let ((mc (monte-carlo 10000000 1000000 random-point-in-2x2-square in-unit-circle? passed:samples->pi)))
(printf "exact = ~a~%inexact = ~a~%(pi - guess) = ~a~%" mc (exact->inexact mc) (- pi mc)))</langsyntaxhighlight>
{{out}}
<pre>1000000 samples of 10000000: 785763 passed -> 785763/250000
Line 2,369 ⟶ 2,700:
A little more Racket-like is the use of an iterator (in this case '''for/fold'''),
which is clearer than an inner function:
<langsyntaxhighlight lang="racket">#lang racket
(define (in-unit-circle? x y) (<= (sqrt (+ (sqr x) (sqr y))) 1))
;; Good idea made in another task that:
Line 2,396 ⟶ 2,727:
;; to see how it looks as a decimal we can exact->inexact it
(let ((mc (monte-carlo/2 10000000 1000000 random-point-in-unit-square in-unit-circle? passed:samples->pi)))
(printf "exact = ~a~%inexact = ~a~%(pi - guess) = ~a~%" mc (exact->inexact mc) (- pi mc)))</langsyntaxhighlight>
 
[Similar output]
Line 2,404 ⟶ 2,735:
{{works with|rakudo|2015-09-24}}
We'll consider the upper-right quarter of the unitary disk centered at the origin. Its area is <math>\pi \over 4</math>.
<syntaxhighlight lang="raku" perl6line>my @random_distances = ([+] rand**2 xx 2) xx *;
 
sub approximate_pi(Int $n) {
Line 2,413 ⟶ 2,744:
say "$_ iterations: ", approximate_pi $_
for 100, 1_000, 10_000;
</syntaxhighlight>
</lang>
{{out}}
<pre>Monte-Carlo π approximation:
Line 2,422 ⟶ 2,753:
We don't really need to write a function, though. A lazy list would do:
 
<syntaxhighlight lang="raku" perl6line>my @pi = ([\+] 4 * (1 > [+] rand**2 xx 2) xx *) Z/ 1 .. *;
say @pi[10, 1000, 10_000];</langsyntaxhighlight>
 
=={{header|REXX}}==
A specific─purpose commatizer function is included to format the number of iterations.
<langsyntaxhighlight lang="rexx">/*REXX program computes and displays the value of pi÷4 using the Monte Carlo algorithm*/
numeric digits 20 /*use 20 decimal digits to handle args.*/
parse arg times chunk digs r? . /*does user want a specific number? */
Line 2,460 ⟶ 2,791:
exit 0 /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
commas: procedure; arg _; do k=length(_)-3 to 1 by -3; _=insert(',',_,k); end; return _</langsyntaxhighlight>
{{out|output|text=&nbsp; when using the default input:}}
<pre>
Line 2,479 ⟶ 2,810:
 
=={{header|Ring}}==
<langsyntaxhighlight lang="ring">
decimals(8)
see "monteCarlo(1000) = " + monteCarlo(1000) + nl
Line 2,492 ⟶ 2,823:
t = (4 * n) / t
return t
</syntaxhighlight>
</lang>
Output:
<pre>
Line 2,498 ⟶ 2,829:
monteCarlo(10000) = 3.00320000
monteCarlo(100000) = 2.99536000
</pre>
 
=={{header|RPL}}==
≪ 0
1 3 PICK '''START'''
RAND SQ RAND SQ + 1 < +
'''NEXT'''
SWAP / 4 *
≫ '<span style="color:blue">MCARL</span>' STO
 
100 <span style="color:blue">MCARL</span>
1000 <span style="color:blue">MCARL</span>
10000 <span style="color:blue">MCARL</span>
100000 <span style="color:blue">MCARL</span>
{{out}}
<pre>
4: 3.2
3: 3.084
2: 3.1684
1: 3.14154
</pre>
 
=={{header|Ruby}}==
<langsyntaxhighlight lang="ruby">def approx_pi(throws)
times_inside = throws.times.count {Math.hypot(rand, rand) <= 1.0}
4.0 * times_inside / throws
Line 2,508 ⟶ 2,859:
[1000, 10_000, 100_000, 1_000_000, 10_000_000].each do |n|
puts "%8d samples: PI = %s" % [n, approx_pi(n)]
end</langsyntaxhighlight>
{{out}}
<pre> 1000 samples: PI = 3.2
Line 2,517 ⟶ 2,868:
 
=={{header|Rust}}==
<langsyntaxhighlight Rustlang="rust">extern crate rand;
 
use rand::Rng;
Line 2,549 ⟶ 2,900:
println!("{:9}: {:<11} dev: {:.5}%", samples, estimate, deviation);
}
}</langsyntaxhighlight>
{{out}}
<pre>Real pi: 3.141592653589793
Line 2,560 ⟶ 2,911:
 
=={{header|Scala}}==
<langsyntaxhighlight lang="scala">object MonteCarlo {
private val random = new scala.util.Random
 
Line 2,586 ⟶ 2,937:
}
}
}</langsyntaxhighlight>
{{out}}
<pre>10000 simulations; pi estimation: 3.1492
Line 2,595 ⟶ 2,946:
 
=={{header|Seed7}}==
<langsyntaxhighlight lang="seed7">$ include "seed7_05.s7i";
include "float.s7i";
 
Line 2,620 ⟶ 2,971:
writeln(" 10000000: " <& pi( 10000000) digits 5);
writeln("100000000: " <& pi(100000000) digits 5);
end func;</langsyntaxhighlight>
 
{{out}}
Line 2,633 ⟶ 2,984:
=={{header|SequenceL}}==
First solution is serial due to the use of random numbers. Will always give the same result for a given n and seed
<syntaxhighlight lang="sequencel">
<lang sequenceL>
import <Utilities/Random.sl>;
import <Utilities/Conversion.sl>;
Line 2,657 ⟶ 3,008:
result when n < 0 else
monteCarloHelper(n - 1, yRand.Generator, newResult);
</syntaxhighlight>
</lang>
 
The second solution will run in parallel. It will also always give the same result for a given n and seed. (Note, the function monteCarloHelper is the same in both versions).
 
<syntaxhighlight lang="sequencel">
<lang sequenceL>
import <Utilities/Random.sl>;
import <Utilities/Conversion.sl>;
Line 2,687 ⟶ 3,038:
result when n < 0 else
monteCarloHelper(n - 1, yRand.Generator, newResult);
</syntaxhighlight>
</lang>
 
=={{header|Sidef}}==
<langsyntaxhighlight lang="ruby">func monteCarloPi(nthrows) {
4 * (^nthrows -> count_by {
hypot(1.rand(2) - 1, 1.rand(2) - 1) < 1
Line 2,698 ⟶ 3,049:
for n in [1e2, 1e3, 1e4, 1e5, 1e6] {
printf("%9d: %07f\n", n, monteCarloPi(n))
}</langsyntaxhighlight>
{{out}}
<pre>
Line 2,707 ⟶ 3,058:
1000000: 3.142344
</pre>
 
=={{header|SparForte}}==
As a structured script.
<syntaxhighlight lang="ada">#!/usr/local/bin/spar
pragma annotate( summary, "monte" )
@( description, "A Monte Carlo Simulation is a way of approximating the" )
@( description, "value of a function where calculating the actual value is" )
@( description, "difficult or impossible. It uses random sampling to define" )
@( description, "constraints on the value and then makes a sort of 'best" )
@( description, "guess.'" )
@( description, "" )
@( description, "Write a function to run a simulation like this with a" )
@( description, "variable number of random points to select. Also, show the" )
@( description, "results of a few different sample sizes. For software" )
@( description, "where the number pi is not built-in, we give pi to a couple" )
@( description, "of digits: 3.141592653589793238462643383280 " )
@( see_also, "http://rosettacode.org/wiki/Monte_Carlo_methods" )
@( author, "Ken O. Burtch" );
pragma license( unrestricted );
 
pragma restriction( no_external_commands );
 
procedure monte is
function pi_estimate (throws : positive) return float is
inside : natural := 0;
begin
for throw in 1..throws loop
if numerics.random ** 2 + numerics.random ** 2 <= 1.0 then
inside := @ + 1;
end if;
end loop;
return 4.0 * float (inside) / float (throws);
end pi_estimate;
 
begin
? " 1_000:" & strings.image (pi_estimate ( 1_000))
@ " 10_000:" & strings.image (pi_estimate ( 10_000))
@ " 100_000:" & strings.image (pi_estimate ( 100_000))
@ " 1_000_000:" & strings.image (pi_estimate ( 1_000_000));
end monte;</syntaxhighlight>
 
=={{header|Stata}}==
<langsyntaxhighlight lang="stata">program define mcdisk
clear all
quietly set obs `1'
Line 2,725 ⟶ 3,117:
 
. mcdisk 100000000
3.1416253</langsyntaxhighlight>
 
=={{header|Swift}}==
{{trans|JavaScript}}
<langsyntaxhighlight Swiftlang="swift">import Foundation
 
func mcpi(sampleSize size:Int) -> Double {
Line 2,754 ⟶ 3,146:
println(mcpi(sampleSize: 1000000))
println(mcpi(sampleSize: 10000000))
println(mcpi(sampleSize: 100000000))</langsyntaxhighlight>
{{out}}
<pre>
Line 2,767 ⟶ 3,159:
 
=={{header|Tcl}}==
<langsyntaxhighlight lang="tcl">proc pi {samples} {
set i 0
set inside 0
Line 2,781 ⟶ 3,173:
foreach runs {1e2 1e4 1e6 1e8} {
puts "$runs => [pi $runs]"
}</langsyntaxhighlight>
result
<pre>PI is approx 3.141592653589793
Line 2,791 ⟶ 3,183:
 
=={{header|Ursala}}==
<langsyntaxhighlight Ursalalang="ursala">#import std
#import flo
 
mcp "n" = times/4. div\float"n" (rep"n" (fleq/.5+ sqrt+ plus+ ~~ sqr+ minus/.5+ rand)?/~& plus/1.) 0.</langsyntaxhighlight>
Here's a walk through.
* <code>mcp "n" = </code>... defines a function named <code>mcp</code> in terms of a dummy variable <code>"n"</code>, which will be the number of iterations used in the simulation
Line 2,812 ⟶ 3,204:
* The result of the division is quadrupled by <code>times/4.</code>.
test program:
<langsyntaxhighlight Ursalalang="ursala">#cast %eL
 
pis = mcp* <10,100,1000,10000,100000,1000000></langsyntaxhighlight>
{{out}}
<pre><
Line 2,827 ⟶ 3,219:
{{trans|Kotlin}}
{{libheader|Wren-fmt}}
<langsyntaxhighlight ecmascriptlang="wren">import "random" for Random
import "./fmt" for Fmt
 
var rand = Random.new()
Line 2,850 ⟶ 3,242:
Fmt.print("$9d -> $10.8f -> $6.4f", n, pi, err)
n = n * 10
}</langsyntaxhighlight>
 
{{out}}
Line 2,866 ⟶ 3,258:
 
=={{header|XPL0}}==
<langsyntaxhighlight XPL0lang="xpl0">code Ran=1, CrLf=9;
code real RlOut=48;
 
Line 2,886 ⟶ 3,278:
RlOut(0, MontePi( 1_000_000)); CrLf(0);
RlOut(0, MontePi(100_000_000)); CrLf(0);
]</langsyntaxhighlight>
 
{{out}}
Line 2,897 ⟶ 3,289:
 
=={{header|zkl}}==
<langsyntaxhighlight lang="zkl">fcn monty(n){
inCircle:=0;
do(n){
Line 2,904 ⟶ 3,296:
}
4.0*inCircle/n
}</langsyntaxhighlight>
Or, in a more functional style (using a reference for state info):
<langsyntaxhighlight lang="zkl">fcn monty(n){
4.0 * (1).pump(n,Void,fcn(r){
x:=(0.0).random(1); y:=(0.0).random(1);
Line 2,912 ⟶ 3,304:
r
}.fp(Ref(0)) ).value/n;
}</langsyntaxhighlight>
{{out}}
<pre>
9,476

edits