Arithmetic Approximations of the Perimeter of an Ellipse
David W. Cantrell
approximations are computationally efficient, using only arithmetic
operations, and have |rel. error| < 10^-5.
There are now two articles on the web which present extensive surveys of
approximations of the perimeter of an ellipse. One article is Gerard
Michon's at <http://www.numericana.com/answer/ellipse.htm>; the other,
which appeared quite recently, is Stanislav Sykora's at
<http://www.ebyte.it/library/docs/math05a/EllipsePerimeterApprox05.html>.
Near the end of the latter, Sykora says "Those who wish to pursue this
topic further should look for an efficient formula using only the four
algebraic operations (if possible, avoiding even square-root) with a
maximum error below 10 ppm. If would be also nice if such a formula were
exact for both the circle (y=1) and the degenerate flat ellipse (y=0)."
Here I provide two approximations which meet those criteria.
With a and b denoting the lengths of the semiaxes of the ellipse, both
approximations have the form [for viewing in a fixed-width font]
2 ( a b ) 2
p (a + b) + q a b + r ( ----- )
a b ( a + b )
4(a + b) - ----- * ----------------------------------
a + b 2 ( a b ) 2
(a + b) + s a b + t ( ----- )
( a + b )
or [written linearly]
4(a+b) - ab/(a+b) * (p(a+b)^2 + qab + r(ab/(a+b))^2)/((a+b)^2 + sab + t(ab/(a+b))^2)
where p, q, r, s and t are parameters to be determined.
Of course, there are many other ways this form could be expressed, some of
which might be preferable for aesthetic or computational reasons. [For
example, the form can easily be put in terms of the three Pythagorean
means: Replace each occurrence of ab/(a+b) by H/2, where H is the harmonic
mean of a and b. Then similarly, replace each occurrence of (a+b) and ab by
2A and G^2, where A and G are the arithmetic and geometric means, resp.]
If eccentricity is 1 (so that the ellipse is degenerate, a=0 or b=0), the
form gives the perimeter precisely. At the other extreme, when eccentricity
is 0 (so that the ellipse is a circle, a=b), if we require that the form
also give the perimeter precisely, it is simple to solve for one parameter
in terms of the others. For example, solving for r gives
r = 4 [ (4 - pi)(4s + t + 16) - (4p + q) ]
which we shall use henceforth. Our approximations will therefore be
extreme-perfect (i.e., give the perimeter precisely at both extremes of
eccentricity). Values of the other parameters (p, q, s and t) then remain
to be determined according to other criteria.
Historically, in approximating the ellipse perimeter, attention has most
often been given to making the approximation as accurate as possible for
near-circular ellipses. Here, however, we wish |relative error| to be small
over _all_ eccentricities. As such, an obvious way to determine values of
the four other parameters is to minimize worst |rel. error|. Doing this
numerically, we find
p = 3.982901, q = 66.71674, s = 18.31287 and t = 23.39728
(and accordingly r = 56.2007)
and that our form of approximation then has |rel. error| < 5.2 * 10^-6.
This approximation meets all of Sykora's criteria. [However, despite the
presence of a square root, I still prefer my second approximation
in "Two New Approximations, in a Certain Form, for the Perimeter of an
Ellipse" (sci.math, 2004 May 23)
at <http://groups.google.com/group/sci.math/msg/618a2c0f319dd8e4>
since it is neater and slightly more accurate.]
Looking at the value of p given above, one might ask "Should that not be
exactly 4 perhaps?" Answer: No, not if we really wish to minimize
worst |rel. error| for an extreme-perfect approximation in the form being
considered. However, if we do choose p = 4 exactly, the approximation then
becomes more accurate for near-degenerate ellipses and, with the values of
the three remaining parameters adjusted appropriately, worst |rel. error|
is still less than 10^-5. Specifically, with
p = 4, q = 129.79130, s = 34.53115 and t = 63.29008
(and accordingly r = 163.3562),
our form of approximation then has |rel. error| < 9.4 * 10^-6. Readers may
be surprised that such a small change in the value of p induced such large
changes in the other parameter values, but I think this merely shows that
this form of approximation is quite sensitive to the values of the
parameters. This matter is related to the following
Disclaimer: My numerically approximated parameter values may differ
somewhat from the optimal ones! (In other words, their last few digits
may well be spurious.)
Nonetheless, for the two sets of parameter values given, the error bounds
stated above are correct. And "refinement" of the parameter values should
not yield any significant reduction in the error bounds.
David W. Cantrell
Paul Abbott
Indeed (see below). One possible exact symmetric formula for the
perimeter of an ellipse is the rapidly converging Gauss-Kummer series:
Pi (a + b) Hypergeometric2F1[-1/2, -1/2, 1, h]
where h is the "standard" (symmetric) variable
h = (a - b)^2 / (a + b)^2 = 1 - (4 a b)/(a + b)^2
Dividing through your approximation by (a + b), one can write it as a
[3,2] rational expression in h.
> If eccentricity is 1 (so that the ellipse is degenerate, a=0 or b=0), the
> form gives the perimeter precisely. At the other extreme, when eccentricity
> is 0 (so that the ellipse is a circle, a=b), if we require that the form
> also give the perimeter precisely, it is simple to solve for one parameter
> in terms of the others. For example, solving for r gives
>
> r = 4 [ (4 - pi)(4s + t + 16) - (4p + q) ]
>
> which we shall use henceforth. Our approximations will therefore be
> extreme-perfect (i.e., give the perimeter precisely at both extremes of
> eccentricity). Values of the other parameters (p, q, s and t) then remain
> to be determined according to other criteria.
There is a neater way of expressing such extreme-perfect approximations.
The linear extreme-perfect approximation, given by interpolation, is
(a + b) (4 h + (1 - h) Pi)
since eccentricity = 1 => h = 1 and eccentricity = 0 => h = 0. So, if
you subtract off and then add back on the linear extreme-perfect
approximation, you can express your approximation as
h (1 - h) (h - c/d)
(a + b) (4 h + (1 - h) Pi + d ------------------- )
t (h - e) (h - f)
where
{c -> -128 + 8 p + 32 Pi + q - 16 s + 4 Pi s,
d -> -64 + 4 p + 16 Pi + q - 16 s + 4 Pi s,
e -> 1 + 2 s/t + 2 Sqrt[s^2 - 4 t]/t,
f -> 1 + 2 s/t - 2 Sqrt[s^2 - 4 t]/t}
This form has several advantages:
[1] Clearly, it is extreme perfect (h (1 - h) = 0 at h = 0 and h = 1);
[2] Generalization is trivial (more numerator or denominator terms);
[3] The numerically determined parameters are ~1 and do not vary greatly;
[4] Location of zeros and poles in rational term is transparent.
> Historically, in approximating the ellipse perimeter, attention has most
> often been given to making the approximation as accurate as possible for
> near-circular ellipses. Here, however, we wish |relative error| to be small
> over _all_ eccentricities. As such, an obvious way to determine values of
> the four other parameters is to minimize worst |rel. error|. Doing this
> numerically, we find
>
> p = 3.982901, q = 66.71674, s = 18.31287 and t = 23.39728
>
> (and accordingly r = 56.2007)
>
> and that our form of approximation then has |rel. error| < 5.2 * 10^-6.
Your approximation reads
0.257903 (h - 1.3641)
(a + b) (4 h + (1 - h) Pi + h (1 - h) -------------------------- )
(h - 3.89452) (h - 1.23625)
Note that the numerically determined parameters are ~1 and that there
are no zeros or poles of the rational term in (0, 1).
> Disclaimer: My numerically approximated parameter values may differ
> somewhat from the optimal ones! (In other words, their last few digits
> may well be spurious.)
Actually, more than the last few digits.
> Nonetheless, for the two sets of parameter values given, the error bounds
> stated above are correct. And "refinement" of the parameter values should
> not yield any significant reduction in the error bounds.
This is not correct. For the parameters
{p -> 3.98454, q -> 77.5672, s -> 21.0894, t -> 31.0392}
the approximation reads
0.237042 (h - 1.29951)
(a + b) (4 h + (1 - h) Pi + h (1 - h) -------------------------- )
(h - 3.51262) (h - 1.20515)
and has |rel. error| < 4.3 * 10^-6.
Note that, dropping the extreme-perfect requirement, computing the best
[N,M] minimax approximation is straightforward. E.g., using Mathematica,
here is the best [3,2] approximant (after dividing out (a + b)):
<< NumericalMath`
approx = MiniMaxApproximation[Pi Hypergeometric2F1[-1/2, -1/2, 1, h],
{h, {0, 1}, 3, 2}][[2,1]]
(0.14756854 h^3 - 0.08368956 h^2 - 2.7172838 h + 3.1415843)/
(0.2370797 h^2 - 1.1150345 h + 1)
This formula has |rel. error| < 2.7 * 10^-6. One can reduce the |rel.
error| by considering more general [N,M] rational approximants.
Cheers,
Paul
_______________________________________________________________________
Paul Abbott Phone: 61 8 6488 2734
School of Physics, M013 Fax: +61 8 6488 1014
The University of Western Australia (CRICOS Provider No 00126G)
AUSTRALIA http://physics.uwa.edu.au/~paul
David W. Cantrell
BTW, this simple approximation (not written in terms of h, however) was
discussed by me in a previous article.
You're right, I'm embarassed to say!
> > Nonetheless, for the two sets of parameter values given, the error
> > bounds stated above are correct. And "refinement" of the parameter
> > values should not yield any significant reduction in the error bounds.
>
> This is not correct. For the parameters
>
> {p -> 3.98454, q -> 77.5672, s -> 21.0894, t -> 31.0392}
>
> the approximation reads
>
> 0.237042 (h - 1.29951)
> (a + b) (4 h + (1 - h) Pi + h (1 - h) -------------------------- )
> (h - 3.51262) (h - 1.20515)
>
> and has |rel. error| < 4.3 * 10^-6.
However, your parameter values are not quite optimal either. In fact we can
get |rel. error| < 3.5 * 10^-6 using
0.235493 (h - 1.28588)
(a + b) (4 h + (1 - h) Pi + h (1 - h) -------------------------- )
(h - 3.47665) (h - 1.19687)
Not without some trepidation, I'll now say again that "refinement" of these
parameter values should not yield any significant reduction in the error
bound.
> Note that, dropping the extreme-perfect requirement, computing the best
> [N,M] minimax approximation is straightforward.
Of course this can be done. But at least for me, the approximation above
is "nicer" than that below. Indeed I think there are no well known
approximations of the ellipse perimeter which fail to give the perimeter
exactly when eccentricity is 0.
Cheers,
David
David W. Cantrell
> Paul Abbott <pa...@physics.uwa.edu.au> wrote:
> > In article <20060112180301.278$a...@newsreader.com>,
> > David W. Cantrell <DWCan...@sigmaxi.org> wrote:
> > > There are now two articles on the web which present extensive surveys
> > > of approximations of the perimeter of an ellipse. One article is
> > > Gerard Michon's at <http://www.numericana.com/answer/ellipse.htm>;
> > > the other, which appeared quite recently, is Stanislav Sykora's at
> > > <http://www.ebyte.it/library/docs/math05a/EllipsePerimeterApprox05.html>
> > > With a and b denoting the lengths of the semiaxes of the ellipse
and with h = ((a - b)/(a + b))^2:
> > There is a neater way of expressing such extreme-perfect
> > approximations. The linear extreme-perfect approximation, given by
> > interpolation, is
> >
> > (a + b) (4 h + (1 - h) Pi)
>
> BTW, this simple approximation (not written in terms of h, however) was
> discussed by me in a previous article.
The next-most-complicated extreme-perfect approximation would seem to be of
the form
(a + b) (4 h + (1 - h) Pi + c h (1 - h))
where the constant factor c is to be determined. Choosing c so that
worst |rel. error| is minimized, we get c = -0.0920186 and the
approximation then has |rel. error| < 0.00069 . I mention this
approximation because (1) its accuracy is quite good considering its
simplicity and (2) I cannot recall having seen it mentioned before.
Of course, as Paul had indicated, if we replace c by rational functions in
h (possibly just polynomials) of higher and higher degree, we can reduce
worst |rel. error| more and more...
Cheers,
David Cantrell