Relay-Version: version B 2.10 5/3/83; site utzoo.UUCP
Posting-Version: version B 2.10.2 9/18/84; site utcsri.UUCP
Path: utzoo!utcsri!hofbauer
From: hofbauer@utcsri.UUCP (John Hofbauer)
Newsgroups: net.micro
Subject: Re: Arctan approximation for generating pi
Message-ID: <1867@utcsri.UUCP>
Date: Sun, 29-Dec-85 00:32:04 EST
Article-I.D.: utcsri.1867
Posted: Sun Dec 29 00:32:04 1985
Date-Received: Sun, 29-Dec-85 02:13:36 EST
References: <1086@brl-tgr.ARPA>
Organization: CSRI, University of Toronto
Lines: 68

>    1.  For numeric approximations which people are actually going to have to
>        use in their work, avoid power series if possible.  Deadly slow, even
>        if, as some will, they converge sufficiently within 5 terms or fewer.
> 
>    The best way I know of to get around both of these is to use
>    continued fractions.  And most fortunately, there is a very nice
>    continued fraction form for arctan which requires perhaps 3 floating
>    point divisions, and 2 floating-point additions.  Furthermore, it
>    gives a very pleasant number of significant digits.  Incredibly few
>    operations for a very good result.

Of course there are many better types of approximations than power series:
continued fractions (as mentioned), rational functions, Chebyshev series,
minimax, to name a few. But which one do you choose? It depends on what
your objectives are. Most frequently you want to compute an approximation
to some desired accuracy. This is the case for builtin functions were you
want sin, for example, to single or double precision. Since you know what
the desired accuracy is at the time you construct the algorithm you usually
pick the type of approximation which will minimize the number of operations
and best preserves numerical stability. It's fascinating to look at the
construction of the builtin function library which comes with every compiler.
The variety of approximations used is astounding; double precision version
of exponential, say, may use an entirely different approximation than the
single precision version simply to save one add or multiply. The problem
is that if the accuracy requirement is changed a new approximation must
be derived from scratch. So, if you want to play around with
generating pi to various number of digits you've got problems,
since unlike with power series it is usually not possible
to generate the next term in the approximation dynamically. Therefore
the power series approximation is not all that bad provided you expand
about a point very close to zero, as some people suggested.

>    2.  Don't raise values to exponents unless it's unavoidable.  The
>        implementations used in most languages use EXP and LN, which are
>        themselves approximations, resulting in representation error,
>        round-off error, and amplification of these in the
>        multiplications which produce the final result.  Wirth knew what

It is true that a**x will be computed as exp(x*ln a), x real, [well known
formula] but a**n will usually be implemented as a*a*...*a (n-1 multiplies)
n integer, at least that's how IBM Fortran compilers do it. But then
in generating the next term in a power series the previous term is
multiplied by x and, usually, divided by some constant so the exponentiation
problem should never arise in practice. Incidently, I once heard a
great horror story how a**x = exp(x*ln a) almost sunk an engineering
Ph.D. thesis here at UofT. It happened nearly 20 years ago before the
builtin functions we all rely on got as good as they are now.

>        Wirth knew what he was doing when he refused to put exponentiation
>	into Pascal.

Wirth probably didn't include exponentiation in Pascal because of
the weird arithmetic on the CDC6600 computer, which was the machine
on which Pascal was first implemented. This is also why Pascal has
packed arrays. The 6600 was a number cruncher and so did not have
byte addressing. Wirth published a paper in Software: Practice and
Experience in 1971 which goes into the problems of implementing
Pascal on the 6600, among other things. It makes your hair stand on
end. Talk about lousy architecture! Under certain circumstances you
could multiply a number by one and change the last few bits!

>    PS.  Like Sheldon Meth, I am very intrigued by the apparently special
> 	relationship between atan(1/5) and atan(1/239).  But is it really true
> 	than pi, a transcendental number, is EQUAL to a linear combination of
> 	these two, or does it only lie within an acceptable tolerance?

What's so unusual about a transcendental number being the sum of two other
transcendental numbers?