5
$\begingroup$

Points are randomly scattered inside the unit square, some fall within the unit circle with probability $P=\pi/4$.

so $P$ is approximated by the fraction $$P\approx \frac{\text{Number of red points}}{\text{Number of all points}}$$ this leads $$\pi \approx 4\frac{\text{Number of red points}}{\text{Number of all points}}$$ (see following image)

There is a code for this:

tinyColor[color_, point_] := {PointSize[Small], color, Point[point]} 
colorChoose[point_] := 
If[Norm[point] <= 1, tinyColor[Red, point], tinyColor[Blue, point]] 
darts = RandomReal[{0, 1}, {40000, 2}]; 
coloredDarts =ParallelMap[colorChoose, darts];
insides = Map[Boole[Norm[#] <= 1] &, darts]; 
piapprox = Accumulate[insides]/Range[Length[darts]]
inner = Select[darts, Norm[#] <= 1 &];
outer = Select[darts, Norm[#] > 1 &];

Show[Plot[Sqrt[1 - x^2], {x, 0, 1}, Filling -> Axis, AspectRatio -> 1,
PlotLabel -> n == Length[darts] TildeTilde[π, 4.0*piapprox[[-1]]]], 
ListPlot[{inner, outer}, 
PlotStyle -> {{PointSize[Tiny], Red}, {PointSize[Tiny], Blue}}, 
ImageSize -> {500, 500}]]

enter image description here

I tried to simplify this problem:

pts = RandomPoint[Rectangle[], 40000];

ListPlot[pts, AspectRatio -> 1, PlotStyle -> Blue]

enter image description here

The problem is following:

How can I split set of points pts into two parts, "inside the circle " and "outside the circle"?

$\endgroup$
  • 1
    $\begingroup$ You should also investigate RegionMember, e.g., rf = RegionMember[Disk[]]; rf[darts] $\endgroup$ – Carl Woll Feb 10 '17 at 22:26
7
$\begingroup$

You can also use Select:

ptsin = Select[pts, Norm[#] < 1 &];

N[Length[ptsin]/Length[pts]]*4
(* 3.1496 *)
$\endgroup$
7
$\begingroup$

How can I split set of points pts into two parts, "inside the circle " and "outside the circle"?

{in, out} = SortBy[GatherBy[pts, Norm[#] < 1 &], Norm[#[[1, 1]]] &];

ListPlot[{in, out}, AspectRatio -> 1, PlotStyle -> {Red, Blue}]

Mathematica graphics

$\endgroup$
  • $\begingroup$ What's the purpose of SortBy? $\endgroup$ – anderstood Feb 4 '17 at 18:27
  • $\begingroup$ The approximation of $\pi$ can be recovered with {in, out} = GatherBy[pts, Norm[#] < 1 &]; 4*Length[in]/(Length[out] + Length[in]) // N. $\endgroup$ – anderstood Feb 4 '17 at 18:30
  • 2
    $\begingroup$ @anderstood, depending on the Norm of the first element in lst, the the first list in the output produced by GatherBy may be the "inside sublist" or the "outside sublist". Sorting the output of GatherBy makes the first list the "inside" one. $\endgroup$ – kglr Feb 4 '17 at 18:32
  • $\begingroup$ Comment to my own comment: it should be {in, out} = SortBy[GatherBy[pts, Norm[#] < 1 &], Norm[#[[1, 1]]] &], cf kglr's comment above. $\endgroup$ – anderstood Feb 4 '17 at 18:55
3
$\begingroup$

Playing with Norm as shown in other answers:

pts = RandomReal[1, {40000, 2}];

4` True/(True + False) /. CountsBy[pts, Norm[#] < 1 &]

3.1474

That doesn't help with drawing the graphic however.

$\endgroup$
2
$\begingroup$

You can compute the norm of the point and verify if it is inside the circle.

At = 4;
d = 2;
totalpoints = 40000;
pts = RandomReal[{-1, 1}, {totalpoints, 2}];
pointsinsidecircle = Select[pts, Norm[#] < 1 &];
counter = Length[pointsinsidecircle];
approxpi = (4. At counter/totalpoints)/d^2;
Print["approx \[Pi] = ", approxpi]
ListPlot[{pts, pointsinsidecircle}, AspectRatio -> 1, 
 PlotStyle -> {Blue, Red}]

enter image description here

(*approx \[Pi] = 3.1328*)
$\endgroup$
  • $\begingroup$ You could you Norm directly, and avoid AppendTo which is very slow. Using Select is probably a much faster option. $\endgroup$ – anderstood Feb 4 '17 at 18:24
  • $\begingroup$ @anderstood i have changed the answer. Thank you. $\endgroup$ – Diogo Feb 4 '17 at 18:30
  • $\begingroup$ Tip: Don't save graphics as JPEGs, they become fuzzy and lose color fidelity. Use PNGs. $\endgroup$ – Rahul Feb 10 '17 at 21:19
2
$\begingroup$

A slightly different approach:

inside = Pick[pts, Map[# ∈ Disk[] &, pts], True];
outside = Complement[pts, inside];

Also as pointed in a comment above:

inside = Pick[#, RegionMember[Disk[]][#], True] &@pts
outside = Complement[pts, inside];
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.