(2002-01-25) Basic Functions
What are some of the common and "special" numerical functions?

Well, people have been inventing special functions ad nauseam.
The list is quite literally endless, but we may attempt the beginning of a classification
for those functions which are common enough to have a universally accepted name.
Let's start with the truly elementary functions:

Polynomial functions:
The value y is obtained from the variable x using only a finite
number of additions and/or multiplications involving given constants.
The simplest such functions are the null function
(zero value, infinite degree)
and the other constant functions
(degree 0; y = a ¹ 0).
Next are linear functions (degree 1, very rarely called monic;
y = ax+b),
quadratic functions (degree 2;
y = ax^{2}+bx+c),
cubic functions,
quartic (or biquadratic), quintic, sextic (rarely hexic),
etc.
Specific qualifiers
are virtually unused for polynomials beyond degree 6:
degree 7 is heptic rather than septic [sic!]
degree 8 is octic, 9 is nonic, 10 is decic.
We're told that some people have called degree 100 hectic.

Rational functions:
The functions you obtain when division is allowed as well.
A rational function is the quotient of two polynomials. The simplest of these is the
reciprocal function y = 1/x.

Algebraic functions: The term applies to any function for which the value y and
the variable x are algebraically related, which is to say that there's a two-variable
polynomial P such that P(x,y) = 0.
Only a few of these functions have an explicit name and/or symbol.
The most notable is the
square root function y = Öx.
Like the square root function, most algebraic functions can only be defined continuously
over the complex plane as multivalued functions.
Alternately, such functions may
be construed as univalued (ordinary) functions of a variable whose domain is
a so-called Riemann surface for which several points may have the same projection
on the complex plane.

Elementary transcendental functions:
The simplest of these is y = e^{x},
the (natural)
exponential function y = exp(x), a function which is equal to its
own derivative (all other such functions are proportional to it).
The exponential function is defined (univalued) over the entire complex plane,
and so are the other transcendental functions which may be defined directly in terms of it,
including trigonometric functions (circular functions) or
hyperbolic functions like:
y = sin(x) = ½ (exp(ix) - exp(-ix))
y = sh(x) = ½ (exp(x) - exp(-x))
Modern usage is to consider only the 3 preferred trigonometric functions
(sine, cosine and tangent)
whereas their 3 reciprocals (cosecant, secant and cotangent) are being deprecated.
A similar remark applies to the 3 preferred hyperbolic functions (sh, ch, th), whose
reciprocals are rarely used, if ever.
Also classified as elementary transcendental functions are the inverses of
the above, starting with the (natural) logarithm function, y = ln(x),
which is the inverse of the exponential (x = exp(y)).
If continuity is required in the realm of complex numbers,
the logarithm function may only be defined as a multivalued function.
The same thing is true of the inverse trigonometric functions (arcsin, arccos, arctg)
or the inverse hyperbolic functions, which complete the modern list of
elementary functions.

Important named combinations
of elementary functions._{ }
One example is the Gudermannian or hyperbolic amplitude,
named after the German mathematician
Christoph
Gudermann (1798-1852):

The fact that gd is odd
is clear from the latter expression (but obfuscated by the former one).
The derivative of gd(x) is 1/ch(x),
and the inverse of the Gudermannian is a primitive of 1/cos(x)...
In a Mercator conformal map,
the distance to the equator of a point at latitude gd(u) is proportional
to u._{ }

sinc(x) = sin(x)/x
is the so-called^{ }sine cardinal or sampling function
(which arises to express Fourier transtorms of
rectangular functions).

Gamma function (see below).
Arguably, the most common special function,
or the least "special" of them.
The other transcendental functions listed below are called "special" because
you could conceivably avoid some of them by staying away from many
specialized mathematical topics.
On the other hand, the Gamma function y = G(x)
is most difficult to avoid.

Elliptic functions and elliptic integrals.

Exponential & logarithm integral (Ei, li),
sine & cosine integral (si, ci).

Riemann's Zeta function.
A simple function with a
rich structure,
best known for its nontrivial zeroes (its trivial zeroes are the even
negative integers).
Infinitely many of these have been shown (by G.H. Hardy)
to have a real part of ½,
and billions of them have actually been found on that critical line,
but it's still not known whether all of them are there,
as Bernhard Riemann
(1826-1866) first conjectured in 1859.
The far-reaching implications of this statement,
known as the Riemann Hypothesis,
make it the
most important unproved mathematical proposition
of our times.

... ...

As advertised, the list is endless...

(2006-11-20) The Trigonometric Circle Diagram
Memorize 3 basic trigonometric functions (out of 6)

Modern usage retains mostly the three trigonometric functions depicted at right:
sine (sin) cosine (cos) and tangent (tan or tg).
The cotangent function (cot or cotg) may be used for
the reciprocal of the tangent.

The names of the reciprocals of the sine and cosine functions are deprecated.
If you must know, the secant (sec) is the reciprocal of the cosine and the
cosecant (csc or cosec) is the reciprocal of the sine function. For the record:

tan q

=

sin q / cos q

sec q

=

1 / cos q

cot q

=

cos q / sin q

csc q

=

1 / sin q

Each of those 6 trigonometric functions is the ratio of two sides in a right triangle where
one of the acute angles is specified. However, it's much better to consider the
so-called
trigometric circle of unit radius depicted above,
which generalizes those functions smoothly
to obtuse and/or negative values of the angle q
(e.g., the cosine of an obtuse angle is negative).

Just memorize the above pictorial definitions of
the 3 main trigonometric functions and ignore the rest
(unless you expect a quizz at school, that is).

(2002-01-29) Trigonometry in a Nutshell
What are the basic laws used with trigonometric functions
to obtain all the elements of a triangle when only some of them are known?

A planar triangle is determined by 3 independent quantities.
These could be the 3 sides (a,b,c,
each of which is less than or equal to the sum of the other two),
1 angle and 2 sides, or 1 side and 2 angles.
Instead of a side, you may be given some other length related to the triangle; for example,
the radius r of the circumscribed circle
(see figure at right).
The angular data is usually (but not always) expressed directly in terms of the
inside angles between sides
(the angles a, b, g, which add up to
p).

The Law of Sines :

2 r =

a

=

b

=

c

sin a

sin b

sin g

The Law of Cosines :

c^{2} = a^{2} + b^{2}
- 2 ab cos g

The Law of Tangents :

a - b

=

tg (a/2 - b/2)

a + b

tg (a/2 + b/2)

The Sine Law may well be the most useful of the three.
It delivers directly any missing quantity in all basic cases not covered
by the next paragraph...

If the 3 sides are given, the Law of Cosines gives you any angle
you may want (via its cosine).
The Law of Cosines will also give you the missing side (c)
when a, b and g are given
(the "SAS case", in high-school parlance).
Having the three sides, you could then obtain either of the missing angles by using
the Law of Cosines again.
However, it's more elegant and more direct to compute
the missing angles with the Law of Tangents
(especially, if you do not care about the value of the missing side c):
Since you know a/2+b/2
(it's equal to p/2-g/2),
the Law of Tangents gives you the (tangent of the) angle
a/2-b/2.
and the missing angles are simply the sum and the
difference of a/2+b/2 and a/2-b/2.

Proofs :

The Law of Sines can be proved by remarking that,
if O is the center of the circumscribed circle,
one may consider an isoceles triangle like OBC which has two sides of length
r forming an angle 2a.
The length of the base (a) is twice the side opposite to an angle
a
in a right triangle of hypotenuse r.

The sine law is traditionally attributed to
Abu Nasr Mansur (970-1036)
who either discovered it himself or got it from
Abu'l-Wafa (940-998).
Abu Nasr Mansur was the mathematical mentor and partner of
Al-Biruni (973-1048).
Historians now agree that
Nasir al-Tusi
(1201-1274) was mistaken when he attributed that law to the astronomer
Al-Khujandi (940-1000)
who never made such a claim himself.

The Law of Cosines is best proved (or memorized) in the modern context of
vector algebra:
Consider two vectors U (from A to C) and V (from A to B),
the vector U-V
goes from B to C and its square (U-V)^{2} is
|U|^{2} - 2 U.V + |V|^{2}.
The scalar product ("dot product") U.V is equal to the cosine of the angle
formed by the two vectors multiplied by the product of their magnitudes.

More traditionally, the law of cosines can be construed
to be equivalent to Heron's formula for the surface area of a triangle,
which appears in Book I of the Metrica by
Heron of Alexandria
(c. AD 10-75).
However, it would be quite a stretch to make it go back this far, since Heron
didn't know what a cosine was...

The Law of Tangents was first stated around 1580 by
François
Viète (Viette, or Vieta).
It has always been the least popular of the three,
and it's not always found in textbooks...
It may be proved by expressing the ratio
(a - b)/(a + b) in terms of
sin a and sin b,
using the Law of Sines.
The result is then immediately obtained from the following identity for the denominator and
from its counterpart for the numerator
(replace b by -b).

sin a + sin b
= 2
sin(a/2 + b/2) cos(a/2 - b/2)

(2003-06-20) Spherical Trigonometry
What are some basic relations applicable to spherical triangles?

A spherical triangle
is a figure on the surface of a sphere of radius R,
featuring three sides which are arcs of great circles
(a "great circle" is the intersection of the sphere with a plane containing the
sphere's center).
Each such figure divides the surface of the sphere into two parts,
whose areas add up to
4pR^{2 }.
Unless otherwise specified, the smaller part is usually considered the "inside"
of the triangle, but this need not always be so...

The study of spherical triangles is often called spherical trigonometry
and is about as ancient as the simpler planar trigonometry summarized
above.

The internal angles of a spherical triangle always add up to more than a
flat angle (of p radians).
Expressed in radians, the difference
(denoted e, with
0< e <4p)
is usually called the spherical excess, a term coined around 1626
by the French-born Dutch mathematician
Albert
Girard (1595-1632),
who showed that the surface area of a spherical triangle is
simply equal to:

e R^{ 2} =
( (a + b + g) - p ) R^{ 2}

There are some striking similarities between the two kinds of trigonometries,
including the spherical Law of Sines of Abu'l-Wafa
(940-998):

sin a

=

sin b

=

sin c

sin a

sin b

sin g

In this, a, b and c are the angular "lengths" of the
sides (as seen from the sphere's center); they are the
curvilinear distances along the
great arcs, using R as a unit.
The spherical excess may also be expressed in terms of these quantities
and the semiperimeters = ½ (a+b+c),
using the spherical equivalent of Hero's formula:

This beautiful formula is named after the Swiss mathematician
Simon
L'Huilier (1750-1840) who was once a teacher of Charles-François Sturm in Geneva.
(His last name is sometimes also spelled "L'Huillier" or "Lhuilier".)

(V. R. of India. 2000-10-16)
Let cos A + cos B = 2p and sin A + sin B = 2q .
Prove that tan A/2 + tan B/2 = 2q / (p^{2} + q^{2} + p)

Let u be tan(A/2) and v be tan(B/2).
We have cos(A) = (1-u^{2})/(1+u^{2}),
sin(A) = 2u/(1+u^{2}),
and similar relations for B and v. Therefore:

Expressions such as these, which are
symmetrical with respect to u and v, may be expressed
in terms of the sum X=u+v and the product Y=uv.
For example u^{2}+v^{2} is X^{2}-2Y and
(1+u^{2})(1+v^{2}) is 1+X^{2}-2Y+Y^{2} or
X^{2}+(1-Y)^{2}.
The above two relations thus become:

Adding or subtracting these two after multiplying each by either (1-Y) or X
(and removing the nonzero factor X^{2}+(1-Y)^{2} which turns up)
greatly simplifies this system of equations, which boils down to a linear system:

pX = q(1-Y) and p (1-Y)+ q X = (1+Y)

This may be rewritten

p X + q Y = q and q X - (p+1) Y = (1-p)

Solving for X gives X [p(p+1)+q^{2}] = q(p+1)+q(1-p) = 2q,
which is the desired relation:

Those two results are equivalent to the statement that the two roots of the
following quadratic equation in t are tan(A/2)
and tan(B/2) :

( p^{ 2} + q^{ 2} + p ) t^{ 2}
- 2 q t +
( p^{ 2} + q^{ 2} - p )
= 0

brentw (Brent Watts
of Hickory, NC. 2001-03-11)
Show that | sin (x + iy) |^{ 2} =
sin^{ 2}(x) + sinh^{ 2}(y)

The well-known relation exp(x+iy)=exp(x)(cos(y)+ i sin(y))
may be turned into a definition of the cosine and sine function
(since x and y need not be real in this).
In particular, exp(iz)=cos(z)+ i sin(z),
so sin(z)=(exp(iz)-exp(-iz))/2i. Therefore:

sin(x+iy) = [e^{ix-y} - e^{-ix+y }] / 2i
= [e^{-y}(cos(x)+ i sin(x)) - exp^{y}(cos(x)- i sin(x))]/2i
= [cos(x)(e^{-y} - e^{y }) + i sin(x)(e^{-y} + e^{y })] / 2i

So far, we did not assume that x and y were real, now we do:
|z|^{ 2} is the sum of the squares of
the real and imaginary parts of z.
When z is the last of the above expressions, this translates into

This expression can be idenfified with the given one by noticing that:

sin^{ 2 }(x) = 1/2 - cos(2x)/2

sinh^{ 2 }(y) =
(e^{ 2y} + e^{-2y })/4 - 1/2

Therefore, the entire expression is indeed
sin^{2}(x) + sh^{2}(y),
as advertised.

(J. S. of Canada. 2000-10-15)
How do you solve these equations to exact values for x?

1) If ln(x-2)-3 = ln(x+1) then ln((x-2)/(x+1)) = ln(e^{3})
so we must have
(x-2) = (x+1)e^{3}
and x can only be equal to
(2+e^{3})/(1-e^{3}).
Now, however, this value of x happens to be negative
(it's about -1.157) which makes it unacceptable,
since both (x-2) and (x+1) should be positive
(or else you can't take their logarithm).
Therefore, the original equation does not have any solutions at all!

2) Rewrite
sin(2x)sin(x)+cos(x) = 0 as 2 cos(x)sin(x)sin(x) + cos(x) = 0,
or cos(x)[2sin^{2}(x)+1] = 0.
As the second factor cannot be zero,
this equation boils down to cos(x) = 0,
which has infinitely many solutions of the form
x = (k+½)p,
where k is any integer (positive or not).

FlyingHellfish
(2003-07-28)
In a broken calculator, only the 6 functions shown
at right are available.
Can any positive rational number be obtained from an initial 0?

Yes. Actually, Ö[p/q]
can be obtained for any positive integers p and q, since:

If p < q, then
Ö[p/q]

= sin arctan Ö[p/(q-p)]

If p > q, then
Ö[p/q]

= tan arccos Ö[q/(p+q)]
= tan arccos sin arctan Ö[q/p]
= tan arccos sin arctan sin arctan Ö[q/(p-q)]

Using whichever relation is relevant, we may thus reduce any case to a simpler one,
until we're faced with p = q,
which we solve by pushing cos once.

There's (almost) no need to say that the above shows that all positive rationals
can be so obtained, since each of them is the square root of its square.

For example, if we wish to obtain 5/8, we observe that it's the square root of
25/64, which is the sin arctan of the
square root of 25/39, itself the sin arctan of the
square root of 25/14, itself the tan arccos sin arctan sin arctan of the
square root of 14/11, itself the tan arccos sin arctan sin arctan of the
square root of 11/3, itself the tan arccos sin arctan sin arctan of the
square root of 3/8, itself the sin arctan of the
square root of 3/5, itself the sin arctan of the
square root of 3/2, itself the tan arccos sin arctan sin arctan of the
square root of 2, itself the tan arccos sin arctan sin arctan of (the
square root of) 1,
which is, of course, the cos of 0...
Only 39 keys to press on that broken calculator.

It's irrelevant whether the calculator works in degrees or in radians,
since we only use trigonometric functions on angles obtained from
inverse trigonometric functions, except for the initial zero angle
(either 0 degrees or 0 radians).

( John of Garland, TX. 2000-11-19)
How are the values of trigonometric functions calculated?
For example, how do we determine that sin(32°) = 0.52991?

Basically, the following relation is used:
sin(x) = x
- x^{3}/6
+ x^{5}/120
- x^{7}/5040
+ x^{9}/362880 - ...
+ (-1)^{k}
x^{2k+1}/(2k+1)! + ...

To use this for actual computations,
you've got to remember that x should be expressed in radians
(1° = p/180 rad).
In your example, x = 32 ° = 0.558505360638... rad.
The series "converges" very rapidly:

After 1 term, S = 0.55850536063818
After 2 terms, S = 0.52946976180816
After 3 terms, S = 0.52992261296708
After 4 terms, S = 0.52991924970365
After 5 terms, S = 0.52991926427444
After 6 terms, S = 0.52991926423312
After 7 terms, S = 0.52991926423332
(no change at this precision after this)

Your computer and/or calculator uses this along with a technique called
economization
(the most popular of which is the Chebyshev economization)
which allows a polynomial of high degree (or any reasonable function)
to be very well approximated by a polynomial of lower degree.

In the case of the sine function, the convergence of the above series is so good
that economization only saves you a couple of multiplications for a given precision.
In some other cases (like the atan function), it is quite indispensable.

Footnote: about atan: The atan function has a nice
Chebyshev expansion which allows one to bypass the intermediates
step of a so-called Taylor expansion like the above.
This is rather fortunate because the convergence of atan's
Taylor expansion is quite lousy when x is close to 1.
Modern atan routines use an economized polynomial for x between 0 and 1,
and reduce the computation of atan(x) to that of atan(1/x)
when x is above 1.
See the following article for more details...

(2000-11-19) What is Chebyshev economization?

Over a finite interval, it is always possible to approximate a continuous function
with arbitrary precision by a polynomial of sufficiently high degree.
In some cases
[one example is the sine function in the previous article]
truncation of the function's Taylor series works well enough.
In other cases, the Taylor series may either converge too slowly or not at all
(the function may not be analytic or, if it is analytic, the radius of
convergence of its Taylor series may be too small to cover comfortably
the desired interval).

If a good polynomial approximation of the continuous real function
f (x)
is desired over a finite interval, the following approach may be used and
is in fact the most popular one.
We may consider without loss of generality
that the desired range of x
is [-1,1]
(if it's not, a linear change of variable will make it so).
Thus, a new variable q
(whose range is [0,p] )
can be introduced via the relation
cos q = x.
Either variable is a decreasing function of the other.

The fundamental remark is that cos(nq) is a polynomial
function of cos(q).
In fact, either of the following relations defines a polynomial of degree n
known as the Chebyshev polynomial [of the first kind] of degree n.
The symbol "T" is conventionally used for these because
of alternate transliterations from Russian, like Tchebycheff
or Tchebychev which are a better match for the Russian pronounciation
(the spellings "Chebychev" and "Tchebyshev" also appear).

The trigonometric formula
cos(n+2)x = 2 cos x cos(n+1)x - cos nx
translates into a simple recurrence relation which makes Chebyshev polynomials
very easy to tabulate: T_{n+2}(x) =
2x T_{n+1}(x) - T_{n}(x)

T_{0}(x)

=

1

T_{1}(x)

=

x

T_{2}(x)

=

-1

+2x^{2}

T_{3}(x)

=

-3x

+4x^{3}

T_{4}(x)

=

1

-8x^{2}

+8x^{4}

T_{5}(x)

=

5x

-20x^{3}

+16x^{5}

T_{6}(x)

=

-1

+18x^{2}

-48x^{4}

+32x^{6}

T_{7}(x)

=

-7x

+56x^{3}

-112x^{5}

+64x^{7}

T_{8}(x)

=

1

-32x^{2}

+160x^{4}

-256x^{6}

+128x^{8}

We must remark prominently that, if y^{2} = x^{2}-1
(y need not be real ), then:

T_{n}(x) = [ (x+y)^{n} + (x-y)^{n} ] / 2

This is a consequence of de Moivre's relation
(with x = cos q
and y = i sin q ):

[ cos q + i sin q ]^{ n}
=
exp(i q)^{ n} =
exp(i nq) =
cos nq + i sin nq

Now, f (cos q)
is clearly an even
function of q which is
continuous when f is.
As such, it has a tame
Fourier expansion which contains only cosines and
translates into the so-called Chebyshev-Fourier expansion of f(x):

^{ }f (cos q) =
c_{o }/2 +

¥

å

n=1

^{ }c_{n }cos(nq)
therefore: f (x) = c_{o }/2 +

¥

å

n=1

^{ }c_{n}T_{n}(x)

The last expression is a series which is always convergent.
For "infinitely smooth" functions, it converges exponentially fast
(as a function of n, the coefficient has to be smaller than the reciprocal of
a polynomial of degree k+1, for any k,
or else the Fourier series of the k-th derivative of
f (cos q) would not converge).
This is much more than what can be said about a
Taylor power series...
A truncated Fourier-Chebyshev series is thus expected to give
a much better approximation than a Taylor series truncated to the same order.

What is known as Chebyshev economization is often limited to the
following dubious technique: Take a good polynomial approximant
with many terms (possibly coming from a Taylor expansion) and express it
as a linear combination of Chebyshev polynomials
(whose coefficients may be obtained from the inversion formula below).
This expression may be truncated at some low order to
obtain a good approximation as a polynomial of lower degree.

A better approach, whenever possible, is to compute the exact
Chebyshev expansion of the target function and to truncate that
in order to obtain a good approximation by a polynomial of low degree...
The following inversion formula can be used for to obtain the
Chebychev expansion
[watch out for the explicit halving of c_{0 }]
of an analytic function given by its Taylor expansion:

f (x) =

¥

å

n=0

a_{n} x^{ n }
= ½ c_{o} +

¥

å

n=1

^{ }c_{n}T_{n}(x)

^{ }c_{n} = 2

¥

å

p=0

ì î

2p+n p

ü þ

a_{ 2p+n}

2^{2p+n}

x

=

T_{1}(x)

2 x^{2}

=

1

+ T_{2}(x)

4 x^{3}

=

3 T_{1}(x)

+ T_{3}(x)

8 x^{4}

=

3

+ 4 T_{2}(x)

+ T_{4}(x)

16 x^{5}

=

10 T_{1}(x)

+ 5 T_{3}(x)

+ T_{5}(x)

32 x^{6}

=

10

+ 15 T_{2}(x)

+ 6 T_{4}(x)

+ T_{6}(x)

2^{ n-1} x^{ n} =

_{(n-1)/2}

å

k=0

ì î

n _{ }k_{ }

ü þ

T_{n-2k }(x) +

½ C(n,^{ n}/_{2 }) if n is even.
0 if n is odd.

The above complete inversion formula (infinite sum)
is occasionally handy, but one may also
always obtain the coefficients c_{n} via the
Euler formulas, which give:

c_{n} =

2

ó õ

^{ 1} _{-1}

f (x) T_{n}(x)

_{ } dx =

2

ó õ

^{ p } _{0}

_{ }f (cos q)
cos(nq) dq

p

p

Ö

1-x^{2}

In at least one (important) case, we may even obtain
the Chebyshev expansion directly by algebraic methods...
Consider the arctangent function, which gives the angle in radians between
-p/2 and p/2
whose tangent equals its given [real] argument.
That function is variously abbreviated
Arctg (Int'l/European), arctan (US), atg or atan (computerese).
The following relation is true for small enough arguments.
[It's true modulo p for unrestricted arguments,
because of the formula giving tg(a+b) as (u+v)/(1-uv) if
u and v are the respective tangents of a and b.]
This may thus be considered an algebraic relation
between formal power series:

Arctg( (u+v)/(1-uv) ) = Arctg(u) + Arctg(v)

With this in mind, we may as well use this formal identity for the complex numbers
u = k [x+iÖ(1-x^{2
})]
and
v = k [x-iÖ(1-x^{2 })], so that
2k^{n} T_{n}(x) =
(u^{ n} + v^{ n}).
This turns the RHS of the above identity directly into a Chebyshev expansion where the
coefficient c_{n} is simply the coefficient of the arctangent power series
multiplied by 2k^{n}.
On the other hand, the LHS becomes
Arctg(2kx/(1-k^{2 })).
If we let k be Ö2-1,
this boils down to Arctg(x) and we have:

That's [almost] all there is to it:
We got the Chebyshev expansion at very little cost!
How good is the convergence of this series?
Well, we may first remark that it converges even if the magnitude of
x exceeds unity.
More precisely, when x is larger than 1, T_{n}(x) is asymptotically
equal to half the n-th power of
x+Ö(x^{2}-1) ,
a quantity which equals the reciprocal of
Ö2-1
when x is Ö2.
Therefore, the series converges if and only if the magnitude of x is
less than (or equal to) Ö2.

More importantly, when the magnitude of x is not more than 1,
a partial sum approximates the whole thing with an error smaller than the coefficient
of the first discarded term.
Suppose we want to use this to find a polynomial approximant of the arctangent function
at a precision of about 13 significant digits
(we need it only over the interval [-1,1], as we may obtain the arctangent of x for x>1
as p/2 minus the arctangent of 1/x).
We find that for 2n+1=31, the relevant coefficient is about
0.88 10^{-13}
so that the corresponding term is just about small enough to be dropped.
The method will thus give the desired precision with an odd polynomial of degree 29,
whose value can be computed using 16 multiplications and 14 additions.
A similar accuracy would require about 10 000 000 000 000
operations with the "straight" Taylor series...
Some economization, indeed!

The above "formal" computation gives the same results as the (unambiguous) relevant
Euler formula for the coefficients of the Chebyshev expansion
of the arctangent function.
This may puzzle a critical reader, since the whole thing seems to work as long as
the quantity 2k/(1-k^{2 }) is equal to unity,
and this quadratic condition is true not only when k is
Ö2-1, but also for the alternate root
-(Ö2+1) as well.
This latter value, however, leads to a formal Chebyshev series which diverges for
any value of x...

(Mark Barnes, UK.
2000-10-24)
What can you
tell
me about the
Gamma function?
I can work out values for G(x)
if x is integral or x is an integer plus one half.
How can I calculate values for
G(x),
if x is some other value,
like 2.8 or 67/9?
What actually is the function?

The following intimidating definitions of the transcendental Gamma function hide
its simple nature:
G(z+1) is merely the generalization of the factorial function
(z!) to all real or complex values of the number z [besides negative integers].

Euler integral of the 2^{nd} kind
(valid only if Re(z)>0):
G(z) =
ò_{0}^{¥}e^{-t} t^{z-1} dt
G(z) =
ò_{1}^{¥}e^{-t} t^{z-1} dt
+ å (-1)^{n}/(n!(n+z))

Weierstrass's definition:^{ } (g
being the Euler-Mascheroni constant, namely
0.5772156649015328606065120900824024310421593359399235988...
):
G(z) = e^{-g z} / z
Õ e^{z/n}/(1+z/n)

G(z) has an elementary
expression only when z is either a positive integer n,
or a positive or negative half-integer
(½+n or ½-n):

G(n)

(n-1)!
G(1/2 + n)

Öp

(2n-1)!!
G(1/2 - n)

(-2)^{n }Öp

2^{n}

(2n-1)!!

In this, k! ("k factorial") is the product of all positive integers less than or equal to k,
whereas k!! ("k double-factorial")
is the product of all such integers which have the same parity as k, namely
k(k-2)(k-4)...
Note that k!, is undefined (¥) when k is a negative integer
(the G function is undefined at z = 0,-1,-2,-3,...
as it has a simple pole at z = -n with a residue of
(-1)^{n}/n! , for any natural integer n).
However, the double factorial k!! may also
be defined for negative odd values of k:
The expression (-2n-1)!! = -(-1)^{n} / (2n-1)!! ) may be obtained
through the recurrence relation (k-2)!! = k!! / k , starting with k=1.
In particular (-1)!! = 1, so that either of the above
formulas does give G(1/2) = Öp , with n=0.
(You may also notice that either relation holds for positive or negative values of n.)

When the real 2x is not an integer, we do not know any expression of
G(x) in terms of elementary functions:

The real [little known] gem which I have to offer about numerical values of the Gamma
function is the so-called "Lanczos approximation formula" [pronounced
"LAHN-tsosh" and named after the Hungarian mathematician
Cornelius
Lanczos (1893-1974), who published it in 1964].
Its form is quite specific to the Gamma function whose values
it gives with superb precision, even for complex numbers.
The formula is valid as long as Re(z)
[the real part of z] is positive.
The nominal accuracy, as I recall, is stated for
Re(z) > ½,
but it's a simple application of the "reflection formula" (given below)
to obtain the value for the rest of the complex plane with a similar accuracy.
The Lanczos formula makes the Gamma function almost as straightforward to compute
as a sine or a cosine. Here it is:

e(z) is a small error term whose value is bounded over
the half-plane described above.
The values of the coefficients C_{i}
depend on the choice of the integers p and n.
For p=5 and n=6, the formula gives a relative error less than
2.2´10^{-10}
with the following choice of coefficients:
C_{1}=76.18009173,
C_{2}= -86.50532033,
C_{3}=24.01409822,
C_{4}= -1.231739516,
C_{5}=0.00120858003,
and C_{6}= -0.00000536382.

I used this particular set of coefficients extensively for years
(other sources may be used for confirmation)
and stated so in my original article here.
This prompted Paul Godfrey
of Intersil Corp. to share a more precise set and his own method to compute any
such sets (without the fear of uncontrolled rounding errors).
Paul has kindly agreed to
let us post his (copyrighted) notes on the subject here.

Some of the fundamental properties of the Gamma function are:

Reflection formula:
G(z)G(1-z) =
p/sin(pz)

Recursion formula: G(1+z) =
zG(z)

Exact values (when n is an integer;
see above when n is negative):
G(n) = (n-1)! and
G(n+1/2) =
Öp (2n)! / (n!4^{n})

Multiplication formula (also called "duplication formula" when n=2):
G(nz) =
(2p)^{(1/2-n/2)}
n^{(nz-1/2)} [
G(z) + ... +
G(z+k/n) + ... +
G(z+(n-1)/n) ]

Other interesting remarks about the Gamma function include:

| G(ix) |^{ 2} =
p / (x sinh px ) for x real

Louis Vlemincq
(Belgium._{ } 2004-02-19; e-mail)
How is the equation t + ln(t) = T ln( I / i ) solved for t and i ?

Taking the exponential of both sides makes it easy to solve for i:

t e^{t} = [I / i]^{ T}

i = I / ( t e^{ t} )^{ 1/T}

To solve for t, you have to use Lambert's W function,
one of the more common "special" functions presented above:
Apply W to both sides of the first of the above equations.
Since, by definition, W(t exp(t)) is equal to t, we obtain:

t = W( [I / i]^{ T })

This solution is valid for positive values of t (the original equation does not make
sense for negative ones). By itself, the equation
x = t exp(t) has 2 real solutions for t when x is
between -1/e and 0 and no real solution when x is less than -1/e.