Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It only takes a minute to sign up.
Sign up to join this community
While I was having breakfast this morning, I thought: "I know how to put random points on a sphere or a torus, but how to place random points on an egg (with uniform probability with respect to any area of the egg)?"
Assume this equation (proposed by Nobuo Yamamoto) describes egg curve: $$ (x^2+y^2)^2 = ax^3 + \frac{3a}{10}xy^2, \tag{1} $$ where $0\leq x\leq a$, $a$ is the length of the axis of symmetry for the egg.
The egg is created by rotation of this curve around the axis of simmetry.
The problem as stated could mean to generate a random evenly distributed over the curve $$ \left( (x^2+y^2)^2 \right) = ax^3 + \alpha a xy^2 $$ in the $xy$ plane (with $\alpha = 3/10$), but you probably meant to ask to generate a random point evenly distributed over the surface of a 3-D egg, described in cylindrical coordinates as $$ \left( (x^2+\rho^2)^2 \right) = ax^3 + \alpha a x\rho^2 $$ where $\rho$ is the radial coordinate and $x$ the linear coordinate. No matter; the basic strategy will be very similar.
The strategy will base the point on an angle $\theta$ uniformly chosen on $(-\pi/2, \pi/2)$ and interpreted as an angle with base along the X axis and vertex at the origin. The line along this angle will intersect the egg curve at some point, and the important quantity is $\frac{ds}{d\theta}$ where $s$ is the arc length along the egg curve from the origin.
Having determined $\frac{ds}{d\theta}$ for the chosen $\theta$, we then generate a uniform variate $q$ on $(0,M)$ where $M$ is some fixed number greater than or equal to the maximum possible value of $\frac{ds}{d\theta}$. Now if $q > \frac{ds}{d\theta}$ we reject these to variates and try again. (Thus there is a greater likelihood of accepting a point in a region where a lot of arc length is packed into a small $\delta \theta$; this is how you get the points to be uniformly distributed on the egg curve.)
If the problem is the one restricted to the XY plane, then this algorithm is the full answer. In the 3-D case, you must generate two further uniform variates, $b$ on $(0,N)$ (where $N$ is some fixed number greater than or equal to the largest possible value of $\rho$) and $\phi$ on $(-\pi,\pi)$. $\phi$ will be the cylindrical azimuth coordinate, but you must first compare the radial coordinate $\rho$ with $b$; if $b > \rho$ you must reject these variates and start all over from the beginning. This rejection step reflects the fact that the area of a ring at radius $\rho$ is proportional to $\rho$; for example, very little of the total surface is very close to the egg end.
So the problem boils down to finding $\frac{ds}{d\theta}$ as a function of $\theta$, and finding suitable values of $M$ greater than the maximum $\frac{ds}{d\theta}$ and (if we are considering the surface of a 3D egg) $N$ greater than the largest vertical radius.
To find $\frac{ds}{d\theta}$ it is convenient to first write $$ u = \tan \theta = y/x $$ Then the curve is $$ x = a\frac{1+\alpha u^2}{(1+u^2)^2} $$ Now the usual arc-length formula will read $$ (ds)^2 = (dx)^2 + (dy)^2 = (dx)^2 + u^2(dx)^2 + x^2(du)^2 + 2ux(du)(dx) $$ and then $$ \frac{ds}{d\theta} = \sqrt{ \left( \frac{dx}{d\theta} \right)^2 \sec^2 \theta + x^2 \left( \frac{du}{d\theta} \right) ^2 + 2ux \frac{du}{d\theta} \frac{dx}{d\theta} } $$ But $$ \frac{dx}{d\theta} = \frac{du}{d\theta} \frac{dx}{du} = (1+u^2)\frac{dx}{du} $$ and $$ \frac{dx}{d\theta} = (1+u^2)\frac{dx}{du} = a \frac{2\alpha u (1+u^2)^2 - 2(1+u^2)(1+\alpha u^2)}{(1+u^2)^3} $$ So to get $\frac{ds}{d\theta}$ in terms of $\theta$, as needed, we plug that expression into $$ \frac{ds}{d\theta} = \sqrt{ \left( \frac{dx}{d\theta} \right)^2 \sec^2 \theta + x^2 \left( \frac{du}{d\theta} \right) ^2 + 2ux \frac{du}{d\theta} \frac{dx}{d\theta} } $$ and remember that $u = \tan \theta$ and $(1+u^2) = \frac{du}{d\theta} = \sec^2 \theta$. (The resulting expression is ugly, but trivial to cleanly produce in a program.)
As to the maximum value of $\frac{ds}{d\theta} $, this is too ugly to be worth finding, but it is fairly easy to show that it cannot exceed $$ M = 6 a \left( \alpha + \frac{1}{\alpha} \right) $$ which for our $\alpha = 0.3$ works out to $21.8a$. Using a larger $M$ than we have to merely impacts the efficiency of the algorithm.
Similarly, we need $N$ for the 3D case, and here it is obvious that $N = a$ is always greater than $\rho$; again, we could tighten this up for better efficiency.
So this solves both forms of the problem.
