Re: circular random number distribution

["Robert G. Brown" <rgb at phy dot duke dot edu>]

> > > ,----
> > > | do {
> > > |     x = gsl_ran_flat(rng, -1, 1);
> > > |     y = gsl_ran_flat(rng, -1, 1);
> > > | } while(sqrt(x*x + y*y) > 1.);
> > > `----
> >
> > In practice this might be faster than the technique described above.
>
> Well, probably. I'll test it on some machines. Does anybody have some
> experience with respect to performance here?
>
> > Try also
> >
> > phi = gsl_ran_flat(rng, 0, 2.0*PI);
> > r = gsl_ran_flat(rng, 0, 1);
> > x = r*cos(phi);
> > y = r*sin(phi);
>
> That comes down to the same problem as above. The correct solution,
> given at http://mathworld.wolfram.com/DiskPointPicking.html, is
> ,----
> | phi = gsl_ran_flat(rng, 0, 2.0*PI);
> | r = gsl_ran_flat(rng, 0, 1);
> | x = sqrt(r)*cos(phi);
> | y = sqrt(r)*sin(phi);

I'm probably late with this (sorry, busy weekend) but:  Your
accept-reject method is almost certainly much faster than using ANYTHING
with transcendental calls in it.  In fact, your sqrt call is redundant
and only makes the routine take longer.  There are likely ways you can
improve your accept-reject efficiency and MAYBE speed the program up a
bit, but even validating them would be costly in human time.

The basic problem is that generalized (hyper)spherical coordinates tend
not to be uniformly distributed on the manifold they cover -- they pack
together at the poles of a sphere, for example.  This actually
complicates doing whole classes of integration on these surfaces, as
covering them with a rectangular grid is just like trying to draw their
surface inside a rectagle -- the result is a projection that distorts
(and VASTLY overweights) the polar regions at the expense of the
equator.  So you encounter exactly the same thing when trying to pick
(theta,phi) randomly and uniformly on the surface of a sphere -- you
cannot just pick them in the ranges 0-PI or 0-2*PI...

  rgb

> `----
>
> Greetings,
> Jochen
> --
> Einigkeit und Recht und Freiheit                http://www.Jochen-Kuepper.de
>     LibertÃ, ÃgalitÃ, Fraternità                GnuPG key: CC1B0B4D
>         (Part 3 you find in my messages before fall 2003.)

Attachment: pgp00000.pgp
Description: PGP signature