2
$\begingroup$

I random chioce three points from unit disk.I want to calculate the probability exactly that they circumcircle is contained within the unit disk still.Of curse,we can use the Monte Carlo to estimate the result.But I try to calculate it exactly,So I use following code

Probability[
 If[And @@ 
   Thread[Norm /@ {x, y, z} <= 1], # + EuclideanDistance[#, x] &@
    RegionCentroid[Circumsphere[{x, y, z}]] < 
   1], {x \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}], 
  y \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}], 
  z \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}]}]

As you see,I get nothing.Can any suggestion give?


I suppose the Circumsphere donnot work in symbol case,so I write a basis function to solve the center of circumcircle.

center := 
  Function[{x1, y1, x2, y2, x3, 
    y3}, {((y3 - y1) (y2 y2 - y1 y1 + x2 x2 - x1 x1) + (y2 - 
        y1) (y1 y1 - y3 y3 + x1 x1 - x3 x3))/(
    2 (x2 - x1) (y3 - y1) - 
     2 (x3 - x1) (y2 - y1)), ((x3 - x1) (x2 x2 - x1 x1 + y2 y2 - 
        y1 y1) + (x2 - x1) (x1 x1 - x3 x3 + y1 y1 - y3 y3))/(
    2 (y2 - y1) (x3 - x1) - 2 (y3 - y1) (x2 - x1))}]

But it seem don't work still.

Probability[
 If[And @@ 
   Thread[Norm /@ {x, y, z} <= 
     1], # + EuclideanDistance[#, x] &@(center @@ Flatten[x, y, z]) < 
   1], {x \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}], 
  y \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}], 
  z \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}]}]
$\endgroup$
4
  • $\begingroup$ Your question is not clear. What does 'they circumcircle is contained within' mean? Moreover, have you checked what your code actually means or produces? e.g. RegionCentroid[Circumsphere[{{0, 0}, {1, 0}, {0, 1}}]] returns {1/2,1/2}, which is presumably not what you want. What is Norm[x] <= 1 && Norm[y] <= 1 && Norm[z] <= 1 going to do for you? $\endgroup$
    – wolfies
    Jul 13, 2016 at 16:54
  • $\begingroup$ @wolfies I have already fix the typo.the norm less 1 just want to restrict the random point within the unit disk $\endgroup$
    – yode
    Jul 13, 2016 at 17:11
  • $\begingroup$ @wolfies: my interpretation is that if he picks three random points uniformly distributed about the unit disk, what is the probability that the circumcircle of those three points is entirely within the unit disk? Recall that a sufficiently thin triangle can have a large circumcircle. $\endgroup$ Jul 13, 2016 at 17:17
  • $\begingroup$ @yode, please could you post final code. $\endgroup$
    – Dragutin
    Jul 14, 2016 at 11:59

1 Answer 1

3
$\begingroup$

I think I've got your idea, you simply want to know when all three points are lying uniformly in a unit disk, how to calculate the possibility of the Circumsphere lying fully in the unit sphere, right?

So, there's quite a lot mistakes in your code:

  1. What is If doing up here?

    Probability can only handle codes that'll generate True and False only, so using If is not a wise choice as it will return nothing if the the first part yield False. I suggest you to use conditioned possibility.

  2. Code will not do your desired job.

    This code will work better I suppose:

    f[x : {a_?NumberQ, b_}, y_, z_] := (And @@ Thread[Norm /@ {x, y, z} <= 1]) 
    && (Norm[#1] + #2 & @@ Circumsphere[{x, y, z}]) < 1;
    
    f[{.5, .5}, {2, 2}, {.3, .3}]
    
    Probability[f[x, y, z],
    {x \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}],
    y \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}], 
    z \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}]}]
    
  3. Possibility is simply not designed for these questions.

    The usage of Possibility is quite limited, check this:

    f[{x1_, x2_}] := x1 > .5;
    
    Probability[f[x], x \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}]]
    NProbability[f[x], x \[Distributed] UniformDistribution[{{-1, 1}, {-1, 1}}]]
    

    It will return nothing wihile it should return something like 1/4.

    I think we can open a new question and ask how to wrestle with this sort of complex definitions.


Edit 1

It seems that Probability will not work when we split those multidimension variables' components. so, simply change the method we use to define those functions can give you a desired result:

f[a : {ax_, ay_}, b : {bx_, by_}, 
   c : {cx_, cy_}] := (And @@ 
     Thread[Norm /@ {a, b, c} <= 1]) && (Norm[#1] + #2 & @@ 
      Circumsphere[{a, b, c}]) < 1;

NProbability[
 f[{ax, ay}, {bx, by}, {cx, cy}], {ax \[Distributed] 
   UniformDistribution[{-1, 1}], 
  ay \[Distributed] UniformDistribution[{-1, 1}], 
  bx \[Distributed] UniformDistribution[{-1, 1}], 
  by \[Distributed] UniformDistribution[{-1, 1}], 
  cx \[Distributed] UniformDistribution[{-1, 1}], 
  cy \[Distributed] UniformDistribution[{-1, 1}]}]

(*0.193461*)

But as I've said before, this is just the probability when point are distributed in a unit square, so the result should be multiplied by 4/Pi which yield the final result of

0.246322

I don't know whether it's correct as some error message was reported(NIntegrate converge slowly), but I think this shall be a correct way.


Edit 2

I used a simple Monte-Carlo method to proof that my code is correct:

n = 1000000;
N[Count[f @@@ RandomReal[{-1, 1}, {n, 3, 2}], True]/n]

(*0.19485*)

It seems that NProbability generate a proper result this time~


Edit 3

Post a new simpler code:

f[{ax_, ay_}, {bx_, by_}, {cx_, 
    cy_}] = (And @@ 
     Thread[Norm /@ {{ax, ay}, {bx, by}, {cx, cy}} <= 
       1]) && (Norm[#1] + #2 & @@ 
      Circumsphere[{{ax, ay}, {bx, by}, {cx, cy}}]) < 1;

Probability[
 f[{ax, ay}, {bx, by}, {cx, cy}], {ax, ay, bx, by, cx, 
   cy} \[Distributed] UniformDistribution[Array[{-1, 1} &, 6]]]

Using {ax, ay, bx, by, cx, cy} \[Distributed] UniformDistribution[Array[{-1, 1} &, 6]] can save you a lot of characters~ :P

Thus, the final result shall be:

0.24632

$\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.