Re: Random points in triangle

*To*: mathgroup at smc.vnet.net*Subject*: [mg111644] Re: Random points in triangle*From*: Ray Koopman <koopman at sfu.ca>*Date*: Sun, 8 Aug 2010 07:23:43 -0400 (EDT)*References*: <i3jc3s$cuh$1@smc.vnet.net>

On Aug 7, 3:21 am, Bill Rowe <readn... at sbcglobal.net> wrote: > On 8/6/10 at 6:55 AM, stev... at ROADRUNNER.COM (S. B. Gray) wrote: > >> I was looking for a simple way to place random points inside a >> triangle with uniform distribution. Here's a good way: >> newtri := Module[{x}, >> ptri = RandomReal[{-5, +5}, {3, 2}]; tredg = Subsets[ptri, {2}]; >> ] >> newpts[nump_] := Module[{wts}, >> inpoints = {}; Do [ wts = RandomReal[GammaDistribution[1, 2], 3]; >> wts = wts/Total[wts]; newin = Total[ptri*wts]; inpoints = >> Append[inpoints, newin], {nump}]; >> ] >> shotri := Module[{x}, >> Graphics[{Blue, Line[tredg], Red, Point[inpoints]}, ImageSize -> 500] >> ] >> The same idea works for points in a tetrahedron; they will be >> uniformly distributed if you use args such as >> GammaDistribution[.6,.1]. > > It seems to me you are making something that should be > reasonably simple complex and doing so in a manner where you > aren't truly getting your stated desired result. Here, I take > your "random points inside a triangle with uniform distribution" > to mean points uniformly distributed over a specified triangle. > > GammaDistribution[1, 2] is the same as ExponentialDistribution[2] > > The sum of n independent items drawn from > ExponentialDistribution[b] will be distributed as > GammaDistribution[n, b] > > The ratio x/(x+y) will have a beta distribution when x and y are > independent values drawn from gamma distributions. > > A uniform distribution can be regarded as a special case of the > beta distribution with a = b = 1. > > So, your variable newwin which is the product of two independent > things each coming from a beta distribution will certainly not > be uniformly distributed. Try this a few times. It looks pretty uniform to me. p = RandomReal[NormalDistribution[0,1],{3,2}]; r = #.p / Total[#,{2}] & @ RandomReal[ExponentialDistribution[1],{5000,3}]; ListPlot[{p,r},PlotRange->All,AspectRatio->1,Frame->True,Axes->None, PlotStyle->{{Red,PointSize[.02]},{Black,PointSize[.005]}}] > However, the actual distribution for > newwin *might* be close enough to a uniform distribution for > your purpose. I haven't taken the time to work out the actual > distribution of newwin. > > But rather than going through all of the above and trying to > work out the exact distribution for newwin, there is a simpler > approach which is guaranteed to give points uniformly > distributed over a given shape. > > To select random points over some planar shape and ensure the > points are uniformly distributed over that shape, simply > generate a pair of uniform deviates using RandomReal[{-1,1}, 2]. > That gives you uniformly distributed points over a square. Now > all you need do to get uniform points over your desired shape is > to drop points that lie outside your desired shape. > > For simplicity, assume your desire is to have points distributed > uniformly over a unit circle. Then the values returned by: > > Cases[RandomReal[{-1, 1},{1000,2}],_?(Norm@#<=1&)] > > will be a set of points uniformly distributed over a circle with > radius 1 > > This can be readily see by doing > > ListPlot[Cases[RandomReal[{-1, 1}, {1000, 2}], _?(Norm@# <= 1 &)], > AspectRatio -> 1] > > The same idea can be extended for solids. For example, > > Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)] > > will be a set of points uniformly distributed through out the > volume of a sphere with radius 1. > > Here, I chose a sphere and circle as examples because it is very > simple to determine whether a randomly selected point lies > inside or not. But the idea of doing this selection is valid no > matter what shape is desired. > > Also, by using Cases, I am not controlling how many points I get > inside the desired area/volume. If you need to have a specific > number of points inside your desired shape, you could do > something like: > > In[12]:= Table[pt = RandomReal[1, 2]; > While[Norm[pt] > 1, pt = RandomReal[1, 2]]; pt, {5}] > > Out[12]= {{0.337431, 0.0873031}, {0.25394, 0.927761}, {0.566057, > 0.391587}, {0.497942, 0.441549}, {0.0179655, 0.854459}} > > which gives you a list of 5 points inside or on the circle with > radius 1. > > And of course it should go with out saying, larger or smaller > shapes can easily be obtained by multiplying your selected > points by the appropriate scale factor. > > There is one other aspect of this approach that should be > mentioned. There are more random numbers being generated than > are actually used. If the desired shape is a very small fraction > of unit square, there will be a large number of the generated > points that are simply thrown away, causing a loss of > efficiency. That is, doing > > Cases[RandomReal[{-10, 10},{100000,3}],_?(Norm@#<=1&)] > > will generate approximately the same number of points in a > circle with radius 1 as > > Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)] > > at the cost of doing about 100 times as much work. > > As long as the area/volume of the desired shape is reasonably > close to a square/cube then this shouldn't be much of a problem. > In fact, this method may out perform your approach since it is > takes less computation to generate a set of uniform deviates > that it does to generated deviates from other distributions.

**Re: Suggestions**

**Re: A new graphic user interface**

**Re: Random points in triangle**

**Re: Random points in triangle**