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
  • 1
    $\begingroup$ You should also investigate RegionMember, e.g., rf = RegionMember[Disk[]]; rf[darts] $\endgroup$
    – Carl Woll
    Feb 10, 2017 at 22:26

5 Answers 5

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$
4
  • $\begingroup$ What's the purpose of SortBy? $\endgroup$
    – anderstood
    Feb 4, 2017 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, 2017 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, 2017 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, 2017 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$
3
  • $\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, 2017 at 18:24
  • $\begingroup$ @anderstood i have changed the answer. Thank you. $\endgroup$
    – Stratus
    Feb 4, 2017 at 18:30
  • $\begingroup$ Tip: Don't save graphics as JPEGs, they become fuzzy and lose color fidelity. Use PNGs. $\endgroup$
    – user484
    Feb 10, 2017 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 and acknowledge that you have read and understand our privacy policy and code of conduct.

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