Error and Complementary Error Functions

Reading Problems

Outline
Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2

Definitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

Gaussian function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
Error function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
Complementary Error function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

Relations and Selected Values of Error Functions . . . . . . . . . . . . . . . . . . . . . . . . 12

Numerical Computation of Error Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19

Rationale Approximations of Error Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

Assigned Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27

1
Background
The error function and the complementary error function are important special functions
which appear in the solutions of diffusion problems in heat, mass and momentum transfer,
probability theory, the theory of errors and various branches of mathematical physics. It
is interesting to note that there is a direct connection between the error function and the
Gaussian function and the normalized Gaussian function that we know as the bell curve.
The Gaussian function is given as

2 /(2 2 )
G(x) = Aex

where is the standard deviation and A is a constant.

The Gaussian function can be normalized so that the accumulated area under the curve is
unity, i.e. the integral from to + equals 1. If we note that the definite integral

r

Z
ax2
e dx =
a

then the normalized Gaussian function takes the form

1 2 /(2 2 )
G(x) = ex
2

If we let

2
x2 1
t = and dt = dx
2 2 2

then the normalized Gaussian integrated between x and +x can be written as

x x
1
Z Z
2
G(x) dx = et dt
x x

or recognizing that the normalized Gaussian is symmetric about the yaxis, we can write

2
 x x  
2 x
Z Z
t2
G(x) dx = e dt = erf x = erf
x 0 2

and the complementary error function can be written as

2
Z
2
erfc x = 1 erf x = et dt
x

Historical Perspective

The normal distribution was first introduced by de Moivre in an article in 1733 (reprinted in
the second edition of his Doctrine of Chances, 1738 ) in the context of approximating certain
binomial distributions for large n. His result was extended by Laplace in his book Analytical
Theory of Probabilities (1812 ), and is now called the Theorem of de Moivre-Laplace.

Laplace used the normal distribution in the analysis of errors of experiments. The important
method of least squares was introduced by Legendre in 1805. Gauss, who claimed to have
used the method since 1794, justified it in 1809 by assuming a normal distribution of the
errors.

The name bell curve goes back to Jouffret who used the term bell surface in 1872 for a
bivariate normal with independent components. The name normal distribution was coined
independently by Charles S. Peirce, Francis Galton and Wilhelm Lexis around 1875 [Stigler].
This terminology is unfortunate, since it reflects and encourages the fallacy that everything
is Gaussian.

3
Definitions
1. Gaussian Function
The normalized Gaussian curve represents the probability distribution with standard
distribution and mean relative to the average of a random distribution.

1 2 /(2 2 )
G(x) = e(x)
2

This is the curve we typically refer to as the bell curve where the mean is zero and
the standard distribution is unity.

2. Error Function
The error function
equals twice the integral of a normalized Gaussian function between
0 and x/ 2.

x
2
Z
2
y = erf x = et dt for x 0, y [0, 1]
0

where
x
t=
2

3. Complementary Error Function

The complementary error function equals one minus the error function

2
Z
2
1 y = erfc x = 1 erf x = et dt for x 0, y [0, 1]
x

4. Inverse Error Function

x = inerf y

inerf y exists for y in the range 1 < y < 1 and is an odd function of y with a
Maclaurin expansion of the form

4

X
inverf y = cn y 2n1
n=1

5. Inverse Complementary Error Function

x = inerfc (1 y)

5
Theory
Gaussian Function
The Gaussian function or the Gaussian probability distribution is one of the most fundamen-
tal functions. The Gaussian probability distribution with mean and standard deviation
is a normalized Gaussian function of the form

1 2 /(2 2 )
G(x) = e(x) (1.1)
2

where G(x), as shown in the plot below, gives the probability that a variate with a Gaussian
distribution takes on a value in the range [x, x + dx]. Statisticians commonly call this
distribution the normal distribution and, because of its shape, social scientists refer to it as
the bell curve. G(x) has been normalized so that the accumulated area under the curve
between x + totals to unity. A cumulative distribution function, which totals
the area under the normalized distribution curve is available and can be plotted as shown
below.
GHxL DHxL

x x
-4 -2 2 4 -4 -2 2 4

Figure 2.1: Plot of Gaussian Function and Cumulative Distribution Function

When the mean is set to zero ( = 0) and the standard deviation or variance is set to unity
( = 1), we get the familiar normal distribution

1 2
G(x) = ex /2 dx (1.2)
2

which is shown in the curve below. The normal distribution function N (x) gives the prob-
ability that a variate assumes a value in the interval [0, x]

x
1
Z
2 /2
N (x) = et dt (1.3)
2 0

6
 NHxL
0.4

0.3

0.2

0.1

x
-4 -2 2 4

Figure 2.2: Plot of the Normalized Gaussian Function

Gaussian distributions have many convenient properties, so random variates with unknown
distributions are often assumed to be Gaussian, especially in physics, astronomy and various
aspects of engineering. Many common attributes such as test scores, height, etc., follow
roughly Gaussian distributions, with few members at the high and low ends and many in
the middle.

Computer Algebra Systems

Function Maple Mathematica

Probability Density Function statevalf[pdf,dist](x) PDF[dist, x ]

- frequency of occurrence at x

Cumulative Distribution Function statevalf[cdf,dist](x) CDF[dist, x ]

- integral of probability
density function up to x dist = normald[, ] dist = NormalDistribution[, ]
= 0 (mean) = 0 (mean)
= 1 (std. dev.) = 1 (std. dev.)

Potential Applications

1. Statistical Averaging:

7
Error Function
The error function is obtained by integrating the normalized Gaussian distribution.

x
2
Z
2
erf x = et dt (1.4)
0

where the coefficient in front of the integral normalizes erf () = 1. A plot of erf x over
the range 3 x 3 is shown as follows.

0.5
erfHxL

-0.5

-1
-3 -2 -1 0 1 2 3
x

Figure 2.3: Plot of the Error Function

The error function is defined for all values of x and is considered an odd function in x since
erf x = erf (x).

The error function can be conveniently expressed in terms of other functions and series as
follows:
 
1 1 2
erf x = ,x (1.5)
2
   
2x 1 3 2
2x x2 3 2
= M , , x = e M 1, , x (1.6)
2 2 2

2 X (1)n x2n+1
= (1.7)
n=0 n!(2n + 1)

where () is the incomplete gamma function, M() is the confluent hypergeometric function
of the first kind and the series solution is a Maclaurin series.

8
Computer Algebra Systems

Function Maple Mathematica

Error Function erf(x ) Erf[x ]

Complementary Error Function erfc(x ) Erfc[x ]

Inverse Error Function fslove(erf(x )=s) InverseErf[s]

Inverse Complementary fslove(erfc(x )=s) InverseErfc[s]

Error Function
where s is a numerical value and we solve for x

Potential Applications

1. Diffusion: Transient conduction in a semi-infinite solid is governed by the diffusion

equation, given as

2T 1 T
=
x2 t

where is thermal diffusivity. The solution to the diffusion equation is a function of

either the erf x or erfc x depending on the boundary condition used. For instance,
for constant surface temperature, where T (0, t) = Ts

T (x, t) Ts
 
x
= erfc
Ti Ts 2 t

9
complementary Error Function
The complementary error function is defined as

erfc x = 1 erf x
Z
2 2
= et dt (1.8)
x

1.5
erfHxL

0.5

-3 -2 -1 0 1 2 3
x

Figure 2.4: Plot of the complementary Error Function

and similar to the error function, the complementary error function can be written in terms
of the incomplete gamma functions as follows:

 
1 1 2
erfc x = ,x (1.9)
2

As shown in Figure 2.5, the superposition of the error function and the complementary error
function when the argument is greater than zero produces a constant value of unity.

Potential Applications

1. Diffusion: In a similar manner to the transient conduction problem described for the
error function, the complementary error function is used in the solution of the diffusion
equation when the boundary conditions are constant surface heat flux, where qs = q0

2q0 (t/)1/2 x2
   
q0 x x
T (x, t) Ti = exp erfc
k 4t k 2 t

10
 1 Erf x + Erfc x
0.8 Erf x

erfHxL
0.6
0.4
0.2 Erfc x

0 0.5 1 1.5 2 2.5 3

x

Figure 2.5: Superposition of the Error and complementary Error Functions


T 
and surface convection, where k = h[T T (0, t)]
x  x=0

 " !#
T (x, t) Ti h2 t
   
x hx x h t
= erfc exp + erfc +
T Ti 2 t k k2 2 t k

11
Relations and Selected Values of Error Functions

erf (x) = erf x erfc (x) = 2 erfc x

erf 0 = 0 erfc 0 = 1

erf = 1 erfc = 0

erf () = 1


Z
erfc x dx = 1/
0

Z
erfc2 x dx = (2 2)/
0

Ten decimal place values for selected values of the argument appear in Table 2.1.

Table 2.1 Ten decimal place values of erf x

x erf x x erf x

0.0 0.00000 00000 2.5 0.99959 30480

0.5 0.52049 98778 3.0 0.99997 79095
1.0 0.84270 07929 3.5 0.99999 92569
1.5 0.96610 51465 4.0 0.99999 99846
2.0 0.99532 22650 4.5 0.99999 99998

12
Approximations

Power Series for Small x (x < 2)

Since

2 x
2 x
(1)n t2n
Z Z
t2
X
erf x = e dt = dt (1.10)
0 0 n=0
n!

and the series is uniformly convergent, it may be integrated term by term. Therefore

2 X (1)n x2n+1
erf x = (1.11)
n=0 (2n + 1)n!

x3 x5 x7 x9
 
2 x
= + + (1.12)
1 0! 3 1! 5 2! 7 3! 9 4!

Asymptotic Expansion for Large x (x > 2)

Since

2 2 1
Z Z
t2 2
erfc x = e dt = et t dt
x x t

we can integrate by parts by letting

1 2
u = dv = et d dt
t
1 2
du = t2 dt v = et
2

therefore

 Z t2

1 t2
  
1 1e
Z Z
t2
e t dt = uv v du = e 2
dt
x t x x 2t x x 2 t

13
Thus

( 2
)
2 1 1
et
Z
x2
erfc x = e dt (1.13)
2x 2 x t2

Repeating the process n times yields

13 1 3 (2n 3)
 
1 x2
1 1 n1
erfc x = e + + (1) +
2 2 x 2x3 22 x5 2n1 x2n1
2
1 3 (2n 1)
et
Z
n
+(1) dt (1.14)
2n x t2n

Finally we can write

2
X 1 3 5 (2n 1)
xex erfc x = 1 + (1)n (1.15)
n=1
(2x2 )n

This series does not converge, since the ratio of the nth term to the (n1)th does not remain
less than unity as n increases. However, if we take n terms of the series, the remainder,

2
1 3 (2n 1)
et
Z
dt
2n x t2n

is less than the nth term because

2

et
dt
Z Z
x2
dt < e <
x t2n 0 t2n

We can therefore stop at any term taking the sum of the terms up to this term as an
approximation of the function. The error will be less in absolute value than the last term
retained in the sum. Thus for large x, erfc x may be computed numerically from the
asymptotic expansion.

x2
X 1 3 5 (2n 1)
xe erfc x = 1 + (1)n
n=1
(2x2 )n

1 13 135
= 1 + + (1.16)
2x2 (2x2 )2 (2x2 )3

14
Some other representations of the error functions are given below:

2 x2 X x2n+1
erf x = e (1.17)
n=0
(3/2)n
 
2x 1 3 2
= M , , x (1.18)
2 2
 
2x x2 3 2
= e M 1, , x (1.19)
2
 
1 1 2
= ,x (1.20)
2
 
1 1 2
erfc x = ,x (1.21)
2

The symbols and represent the incomplete gamma functions, and M denotes the con-
fluent hypergeometric function or Kummers function.

Derivatives of the Error Function

 x 
d 2 d 2
Z
2 t2
erf x = ex = e dt (1.22)
dx dx 0

Use of Leibnitz rule of differentiation of integrals gives:

d2 d 2 2 2 2

2
erfc x = ex = (2x)ex (1.23)
dx dx

d3
 
d 2 x2
2 2

3
erfc x = (2x)e = (4x2 2)ex (1.24)
dx dx

In general we can write

dn+1 n
2 x2
erf x = (1) H n (x) e (n = 0, 1, 2 . . .) (1.25)
dxn+1

where Hn (x) are the Hermite polynomials.

15
Repeated Integrals of the Complementary Error Function
Z
n
i erfc x = in1 erfc t dt n = 0, 1, 2, . . . (1.26)
x

where

2 2
i1 erfc x = ex (1.27)

i0 erfc x = erfc x (1.28)

Z
1
i erfc x = ierfc x = erfc t dt
x

1
= exp(x2 ) x erfc x (1.29)

Z
2
i erfc x = i erfc t dt
x
 
1 2
2 2
= (1 + 2x ) erfc x x exp(x )
4
1
= [erfc x 2x ierfc x] (1.30)
4

The general recurrence formula is

2nin erfc x = in2 erfc x 2xin1 erfc x (n = 1, 2, 3, . . .) (1.31)

Therefore the value at x = 0 is

1
in erfc 0   (n = 1, 0, 1, 2, 3, . . .) (1.32)
n
n
2 1+
2

It can be shown that y = in erfc x is the solution of the differential equation

d2 y dy
+ 2x 2ny = 0 (1.33)
dx2 dx

16
The general solution of

y 00 + 2xy 0 2ny = 0 x (1.34)

is of the form

y = Ain erfc x + Bin erfc (x) (1.35)

Derivatives of Repeated Integrals of the Complementary Error

Function

d
[in erfc x] = (1)n1 erfc x (n = 0, 1, 2, 3 . . .) (1.36)
dx

dn h x2
i 2
e erfc x = (1)n 2n n!ex in erfc x (n = 0, 1, 2, 3 . . .) (1.37)
dxn

17
Some Integrals Associated with the Error Function

x2
et
Z
dt = erf x (1.38)
0 t
Z x
t y

e dt = erf x (1.39)
0 2y
2 x2
1
et
Z
2
ex 1 {erf x}2
 
dt = (1.40)
0 1 + t2 2


et x
Z
dt = exy erfc ( xy) x>0 (1.41)
0 y+t x
2x

et
Z
2
dt = exy erfc ( xy) x > 0, y > 0 (1.42)
0 t2 + y2 2y

etx
Z
dt = exy erfc (xy) x > 0, y 6= 0 (1.43)
0 (t + y) t y

y
Z p
et x
erf ( yt) dt = (x + y)1/2 (x + y) > 0 (1.44)
0 x

1 2xy
Z p
et x erf ( y/t dt = e x > 0, y > 0 (1.45)
0 x
Z
erfc (t) dt = ierfc (a) + 2a = ierfc (a) (1.46)
a
Z a
erf (t) dt = 0 (1.47)
a
Z a
erfc(t) dt = 2a (1.48)
a

1
Z
ierfc (t) dt = i2 erfc (a) = + a i2 erfc (a) (1.49)
a 2
   
t+c a+c
Z
n n+1
i erfc dt = bi erfc (1.50)
a b b

18
Numerical Computation of Error Functions
The power series form of the error function is not recommended for numerical computations
when the argument approaches and exceeds the value x = 2 because the large alternat-
ing terms may cause cancellation, and because the general term is awkward to compute
recursively. The function can, however, be expressed as a confluent hypergeometric series.

 
2 x2
3 2
erf x = x e M 1, , x (1.51)
2

in which all terms are positive, and no cancellation can occur. If we write

X
erf x = b an 0x2 (1.52)
n=0

with

2x 2 x2
b = ex a0 = 1 an = an1 n1
(2n + 1)/2

then erf x can be computed very accurately (e.g. with an absolute error less that 109 ).
Numerical experiments show that this series can be used to compute erf x up to x = 5
to the required accuracy; however, the time required for the computation of erf x is much
greater due to the large number of terms which have to be summed. For x 2 an alternate
method that is considerably faster is recommended which is based upon the asymptotic
expansion of the complementary error function.

2
Z
2
erfc x = et dt
x

2
ex
 
1 1
= 2 Fo , 1, x (1.53)
x 2 x2

which cannot be used to obtain arbitrarily accurate values for any x. An expression that
converges for all x > 0 is obtained by converting the asymptotic expansion into a continued
fraction

19
 2 1
ex erfc x = x>0 (1.54)
1/2
x+
1
x+
3/2
x+
2
x+
5/2
x+
x + ...

which for convenience will be written as

2
ex
 
1 1/2 1 3/2 2
erfc x = x>0 (1.55)
x+ x+ x+ x+ x+

It can be demonstrated experimentally that for x 2 the 16th approximant gives erfc x
with an absolute error less that 109 . Thus we can write

2
ex
 
1 1/2 1 3/2 8
erfc x = x2 (1.56)
x+ x+ x+ x+ x

Using a fixed number of approximants has the advantage that the continued fraction can be
evaluated rapidly beginning with the last term and working upward to the first term.

20
Rational Approximations of the Error Functions (0 x < )
Numerous rational approximations of the error functions have been developed for digital
computers. The approximations are often based upon the use of economized Chebyshev
polynomials and they give values of erf x from 4 decimal place accuracy up to 24 decimal
place accuracy.

Two approximations by Hastings et al.11 are given below.

2
erf x = 1 [t(a1 + t(a2 + a3 t))] ex + (x) 0x (1.57)

where

1
t=
1 + px

and the coefficients are

p = 0.47047

a1 = 0.3480242

a2 = 0.0958798

a3 = 0.7478556

This approximation has a maximum absolute error of |(x)| < 2.5 105 .

Another more accurate rational approximation has been developed for example

2
erf x = 1 [t(a1 + t(a2 + t(a3 + t(a4 + a5 t))))] ex + (x) (1.58)

where

1
t=
1 + px

21
and the coefficients are

p = 0.3275911

a1 = 0.254829592

a2 = 0.284496736

a3 = 1.421413741

a4 = 1.453152027

a5 = 1.061405429

This approximation has a maximum absolute error of |(x)| < 1.5 107 .

22
Assigned Problems
Problem Set for Error and Due Date: February 12, 2004
Complementary Error Function

1. Evaluate the following integrals to four decimal places using either power series, asymp-
totic series or polynomial approximations:
Z 2 Z 0.002
x2 2
a) e dx b) ex dx
0 0.001

10
2 2
Z Z
x2 2
c) e dx d) ex dx
1.5 5

1.5   r  
1 2 1
Z Z
x2 x2
e) e dx f) e dx
1 2 1 2

2. The value of erf 2 is 0.995 to three decimal places. Compare the number of terms
required in calculating this value using:

a) the convergent power series, and

b) the divergent asymptotic series.

Compare the approximate errors in each case after two terms; after ten terms.

3. For the function ierfc(x) compute to four decimal places when x = 0, 0.2, 0.4, 0.8,
and 1.6.

4. Prove that

 
1 2
i) erf (x) = ,x
2

 
1 2
ii) erfc(x) = ,x
2
   
1 2
1 2
where ,x and ,x are the incomplete Gamma functions defined as:
2 2

23
 Z y
(a, y) = eu ua1 du
0

and

Z
(a, y) = eu ua1 du
y

5. Show that (x, t) = 0 erfc(x/2 t) is the solution of the following diffusion
problem:

2 1
= x 0, t > 0
x2 t

and

(0, t) = 0 , constant

(x, t) 0 as x

6. Given (x, t) = 0 erf x/2 t:

i) Obtain expressions for and at any x and all t > 0
t x

x
ii) For the function
2 0 x

showthat it has a maximum value when x/2 t = 1/ 2 and the maximum value
is 1/ 2e .

7. Given the transient point source solution valid within an isotropic half space

q
T = erfc(r/2 t), dA = r dr d
2kr

24
 derive the expression for the transient temperature rise at the centroid of a circular
area (a2 ) which is subjected to a uniform and constant heat flux q. Superposition of
point source solutions allows one to write

Z aZ 2
T0 = T dA
0 0

8. For a dimensionless time F o < 0.2 the temperature distribution within an infinite
plate L x L is given approximately by

T (, F o) Ts 1
 
1+
=1 erfc + erfc
T0 Ts 2 Fo 2 Fo

for 0 1 where = x/L and F o = t/L2 .

Obtain the expression for the mean temperature (T (F o) Ts )/(T0 Ts ) where

Z 1
T = T (, F o) d
0

The initial and surface plate temperature are denoted by T0 and Ts , respectively.

9. Compare the approximate short time (F o < 0.2) solution:

3
(2n 1) (2n 1) +
X  
n+1
(, F o) = 1 (1) erfc + erfc
n=1 2 Fo 2 Fo

and the approximate long time (F o > 0.2) solution

3
X 2(1)n+1 2
(, F o) = en F o cos(n )
n=1
n

with n = (2n 1)/2.

For the centerline ( = 0) compute to four decimal places (0, F o)ST and (0, F o)LT
for F o = 0.02, 0.06, 0.1, 0.4, 1.0 and 2.0 and compare your values with the exact
values given in Table 1.

25
Table 1: Exact values of (0, F o) for the Infinite Plate

Fo (0, F o)

0.02 1.0000
0.06 0.9922
0.10 0.9493
0.40 0.4745
1.0 0.1080
2.0 0.0092

26
References
1. Abramowitz, M. and Stegun, I.A., Handbook of Mathematical Functions, Dover,
New York, 1965.

2. Fletcher, A., Miller, J.C.P., Rosehead, L. and Comrie, L.J., An Index of

Mathematical Tables, Vols. 1 and 2, 2 edition, Addison-Wesley, Reading, Mass., 1962.

3. Hochsadt, H., Special Functions of Mathematical Physics, Holt, Rinehart and Win-
ston, New York, 1961.

4. Jahnke, E., Emdw, F. and Losch, F., Tables of Higher Functions, 6th Edition,
McGraw-Hill, New York, 1960.

5. Lebedev, A.V. and Fedorova, R.M., A Guide to Mathematical Tables, Pergamon

Press, Oxford, 1960.

6. Lebedev, N.N., Special Functions and Their Applications, Prentice-Hall, Englewood

Cliffs, NJ, 1965.

7. Magnus, W., Oberhettinger, F. and Soni, R.P., Formulas and Theorems for the
Functions of Mathematical Physics, 3rd Edition, Springer-Verlag, New York, 1966.

8. Rainville, E.D., Special Functions, MacMillan, New York, 1960.

9. Sneddon, I.N., Special Functions of Mathematical Physics and Chemistry, 2nd Edi-
tion, Oliver and Boyd, Edinburgh, 1961.

10. National Bureau of Standards, Applied Mathematics Series 41, Tables of Error
Function and Its Derivatives, 2nd Edition, Washington, DC, US Government Printing
Office, 1954.

11. Hastings, C., Approximations for Digital Computers, Princeton University Press,
Princeton, NJ, 1955.

27
Chebyshev Polynomials

Reading Problems

Differential Equation and Its Solution

The Chebyshev differential equation is written as

d2 y dy
(1 x2 ) x + n2 y = 0 n = 0, 1, 2, 3, . . .
dx2 dx

If we let x = cos t we obtain

d2 y
+ n2 y = 0
dt2

whose general solution is

y = A cos nt + B sin nt

or as

y = A cos(n cos1 x) + B sin(n cos1 x) |x| < 1

or equivalently

y = ATn (x) + BUn (x) |x| < 1

where Tn (x) and Un (x) are defined as Chebyshev polynomials of the first and second kind
of degree n, respectively.

1
If we let x = cosh t we obtain

d2 y
n2 y = 0
dt2

whose general solution is

y = A cosh nt + B sinh nt

or as

y = A cosh(n cosh1 x) + B sinh(n cosh1 x) |x| > 1

or equivalently
y = ATn (x) + BUn (x) |x| > 1

The function Tn (x) is a polynomial. For |x| < 1 we have

 p n
Tn (x) + iUn (x) = (cos t + i sin t)n = x + i 1 x2
 p n
Tn (x) iUn (x) = (cos t i sin t)n = x i 1 x2

from which we obtain

1 h p n  p n i
Tn (x) = x+i 1 x2 + x i 1 x2
2

For |x| > 1 we have

 p n
Tn (x) + Un (x) = ent = x x2 1
 p n
nt 2
Tn (x) Un (x) = e = x x 1

The sum of the last two relationships give the same result for Tn (x).

2
Chebyshev Polynomials of the First Kind of Degree n
The Chebyshev polynomials Tn (x) can be obtained by means of Rodrigues formula

(2)n n! p dn
Tn (x) = 1 x2 (1 x2 )n1/2 n = 0, 1, 2, 3, . . .
(2n)! dxn

The first twelve Chebyshev polynomials are listed in Table 1 and then as powers of x in
terms of Tn (x) in Table 2.

3
 Table 1: Chebyshev Polynomials of the First Kind

T0 (x) = 1

T1 (x) = x

T2 (x) = 2x2 1

T3 (x) = 4x3 3x

T4 (x) = 8x4 8x2 + 1

T5 (x) = 16x5 20x3 + 5x

T6 (x) = 32x6 48x4 + 18x2 1

T7 (x) = 64x7 112x5 + 56x3 7x

T8 (x) = 128x8 256x6 + 160x4 32x2 + 1

T9 (x) = 256x9 576x7 + 432x5 120x3 + 9x

T10 (x) = 512x10 1280x8 + 1120x6 400x4 + 50x2 1

T11 (x) = 1024x11 2816x9 + 2816x7 1232x5 + 220x3 11x

4
 Table 2: Powers of x as functions of Tn (x)

1 = T0

x = T1
1
x2 = (T0 + T2 )
2
1
x3 = (3T1 + T3 )
4
1
x4 = (3T0 + 4T2 + T4 )
8
1
x5 = (10T1 + 5T3 + T5 )
16
1
x6 = (10T0 + 15T2 + 6T4 + T6 )
32
1
x7 = (35T1 + 21T3 + 7T5 + T7 )
64
1
x8 = (35T0 + 56T2 + 28T4 + 8T6 + T8 )
128
1
x9 = (126T1 + 84T3 + 36T5 + 9T7 + T9 )
256
1
x10 = (126T0 + 210T2 + 120T4 + 45T6 + 10T8 + T10 )
512
1
x11 = (462T1 + 330T3 + 165T5 + 55T7 + 11T9 + T11 )
1024

5
Generating Function for Tn (x)
The Chebyshev polynomials of the first kind can be developed by means of the generating
function

1 tx X
= Tn (x)tn
1 2tx + t2 n=0

Recurrence Formulas for Tn (x)

When the first two Chebyshev polynomials T0 (x) and T1 (x) are known, all other polyno-
mials Tn (x), n 2 can be obtained by means of the recurrence formula

Tn+1 (x) = 2xTn (x) Tn1 (x)

The derivative of Tn (x) with respect to x can be obtained from

(1 x2 )Tn0 (x) = nxTn (x) + nTn1 (x)

Special Values of Tn (x)

The following special values and properties of Tn (x) are often useful:

Tn (x) = (1)n Tn (x) T2n (0) = (1)n

Tn (1) = 1 T2n+1 (0) = 0

Tn (1) = (1)n

6
Orthogonality Property of Tn (x)
We can determine the orthogonality properties for the Chebyshev polynomials of the first
kind from our knowledge of the orthogonality of the cosine functions, namely,

0 (m 6= n)
Z


cos(m) cos(n ) d = /2 (m = n 6= 0)
0



(m = n = 0)

Then substituting

Tn (x) = cos(n)

cos = x

to obtain the orthogonality properties of the Chebyshev polynomials:

0 (m 6= n)
Z 1

Tm (x) Tn (x) dx
= /2 (m = n 6= 0)
1 1 x2



(m = n = 0)

We observe that the Chebyshev polynomials form an orthogonal set on the interval 1
x 1 with the weighting function (1 x2 )1/2

Orthogonal Series of Chebyshev Polynomials

An arbitrary function f (x) which is continuous and single-valued, defined over the interval
1 x 1, can be expanded as a series of Chebyshev polynomials:

f (x) = A0 T0 (x) + A1 T1 (x) + A2 T2 (x) + . . .

X
= An Tn (x)
n=0

7
where the coefficients An are given by

1
1 f (x) dx
Z
A0 = n=0
1 1 x2

and

1
2 f (x) Tn (x) dx
Z
An = n = 1, 2, 3, . . .
1 1 x2

The following definite integrals are often useful in the series expansion of f (x):

1
dx x3 dx
1
Z Z
= = 0
1 1 x2 1 1 x2
1
x dx x4 dx
1
3
Z Z
= 0 =
1 1 x2 1 1 x2 8
1
x2 dx x5 dx
1
Z Z
= = 0
1 1 x2 2 1 1 x2

Chebyshev Polynomials Over a Discrete Set of Points

A continuous function over a continuous interval is often replaced by a set of discrete values
of the function at discrete points. It can be shown that the Chebyshev polynomials Tn (x)
are orthogonal over the following discrete set of N + 1 points xi , equally spaced on ,

2
i = 0, , , . . . (N 1) ,
N N N

where

xi = arccos i

We have

8


0 (m 6= n)
N 1


1 X 1
Tm (1)Tn (1)+ Tm (xi )Tn (xi )+ Tm (1)Tn (1) = N/2 (m = n 6= 0)
2 i=2
2



N (m = n = 0)

The Tm (x) are also orthogonal over the following N points ti equally spaced,

3 5 (2N 1)
i = , , , ...,
2N 2N 2N 2N

and

ti = arccos i

0 (m 6= n)
N


X
Tm (ti )Tn (ti ) = N/2 (m = n 6= 0)

i=1


N (m = n = 0)

The set of points ti are clearly the midpoints in of the first case. The unequal spacing of
the points in xi (N ti ) compensates for the weight factor

W (x) = (1 x2 )1/2

in the continuous case.

9
Additional Identities of Chebyshev Polynomials
The Chebyshev polynomials are both orthogonal polynomials and the trigonometric cos nx
functions in disguise, therefore they satisfy a large number of useful relationships.

The differentiation and integration properties are very important in analytical and numerical
work. We begin with

Tn+1 (x) = cos[(n + 1) cos1 x]

and

Tn1 (x) = cos[(n 1) cos1 x]

Differentiating both expressions gives

1 d[Tn+1 (x)] sin[(n + 1) cos1 x

=
(n + 1) dx 1 x2

and

1 d[Tn1 (x)] sin[(n 1) cos1 x

=
(n 1) dx 1 x2

Subtracting the last two expressions yields

1 d[Tn+1 (x)] 1 d[Tn1 (x)] sin(n + 1) sin(n 1)

=
(n + 1) dx (n 1) dx sin

or

0 0
Tn+1 (x) Tn1 (x) 2 cos n sin
= = 2Tn (x) n2
(n + 1) (n 1) sin

10
Therefore

T20 (x) = 4T1

T10 (x) = T0

T00 (x) = 0

We have the formulas for the differentiation of Chebyshev polynomials, therefore these for-
mulas can be used to develop integration for the Chebyshev polynomials:

 
1 Tn+1 (x) Tn1 (x)
Z
Tn (x)dx = +C n2
2 (n + 1) (n 1)
1
Z
T1 (x)dx = T2 (x) + C
4
Z
T0 (x)dx = T1 (x) + C

The Shifted Chebyshev Polynomials

For analytical and numerical work it is often convenient to use the half interval 0 x 1
instead of the full interval 1 x 1. For this purpose the shifted Chebyshev polynomials
are defined:

Tn (x) = Tn (2x 1)

Thus we have for the first few polynomials

T0 = 1

T1 = 2x 1

T2 = 8x2 8x + 1

T3 = 32x3 48x2 + 18x 1

T4 = 128x4 256x3 + 160x2 32x + 1

11
and the following powers of x as functions of Tn (x);

1 = T0
1
x = (T0 + T1 )
2
1
x2 = (3T0 + 4T1 + T2 )
8
1
x3 = (10T0 + 15T1 + 6T2 + T3 )
32
1
x4 = (35T0 + 56T1 + 28T2 + 8T3 + T4 )
128

The recurrence relationship for the shifted polynomials is:

Tn+1 (x) = (4x 2)Tn (x) Tn1

(x) T0 (x) = 1

or

1 1 1
xTn (x) =
Tn+1 (x) + Tn (x) +
Tn1 (x)
4 2 4

where

Tn (x) = cos n cos1 (2x 1) = Tn (2x 1)

 

Expansion of xn in a Series of Tn (x)

A method of expanding xn in a series of Chebyshev polynomials employes the recurrence
relation written as

1
xTn (x) = [Tn+1 (x) + Tn1 (x)] n = 1, 2, 3 . . .
2
xT0 (x) = T1 (x)

12
To illustrate the method, consider x4

4 2
x2 x
x = x (xT1 ) = [T2 + T0 ] = [T1 + T3 + 2T1 ]
2 4
1 1
= [3xT1 + xT3 ] = [3T0 + 3T2 + T4 + T2 ]
4 8
1 1 3
= T4 + T2 + T0
8 2 8

This result is consistent with the expansion of x4 given in Table 2.

Approximation of Functions by Chebyshev Polynomials

Sometimes when a function f (x) is to be approximated by a polynomial of the form

X
f (x) = an xn + EN (x) |x| 1
n=0

where |En (x)| does not exceed an allowed limit, it is possible to reduce the degree of the
polynomial by a process called economization of power series. The procedure is to convert
the polynomial to a linear combination of Chebyshev polynomials:

N
X N
X
n
an x = bn Tn (x) n = 0, 1, 2, . . .
n=0 n=0

It may be possible to drop some of the last terms without permitting the error to exceed
the prescribed limit. Since |Tn (x)| 1, the number of terms which can be omitted is
determined by the magnitude of the coefficient b.

The Chebyshev polynomials are useful in numerical work for the interval 1 x 1
because

1. |Tn (x)] 1 within 1 x 1

2. The maxima and minima are of comparable value.

13
 3. The maxima and minima are spread reasonably uniformly over the interval
1 x 1

4. All Chebyshev polynomials satisfy a three term recurrence relation.

5. They are easy to compute and to convert to and from a power series form.

These properties together produce an approximating polynomial which minimizes error in

its application. This is different from the least squares approximation where the sum of
the squares of the errors is minimized; the maximum error itself can be quite large. In
the Chebyshev approximation, the average error can be large but the maximum error is
minimized. Chebyshev approximations of a function are sometimes said to be mini-max
approximations of the function.

The following table gives the Chebyshev polynomial approximation of several power series.

14
 Table 3: Power Series and its Chebyshev Approximation

1. f (x) = a0

f (x) = a0 T0

2. f (x) = a0 + a1 x

f (x) = a0 T0 + a1 T1

3. f (x) = a0 + a1 x + a2 x2
   
a2 a2
f (x) = a0 + T0 + a1 T1 + T2
2 2

4. f (x) = a0 + a1 x + a2 x2 + a3 x3
       
a2 3a3 a2 a3
f (x) = a0 + T0 + a1 + T1 + T2 + T3
2 4 2 4

5. f (x) = a0 + a1 x + a2 x2 + a3 x3 + a4 x4
       
a2 a3 3a3 a2 a4 a3
f (x) = a0 + + T0 + a1 + T1 + + T2 + T3
2 8 4 2 2 8
 
a4
+ T4
8

6. f (x) = a0 + a1 x + a2 x2 + a3 x3 + a4 x4 + a5 x5
     
a2 3a4 3a3 5a5 a2 a4
f (x) = a0 + + T0 + a1 + + T1 + + T2
2 8 4 8 2 2
     
a3 5a5 a4 a5
+ + T3 + T4 + T5
4 16 8 16

15
 Table 4: Formulas for Economization of Power Series

x = T1
1
x2 = (1 + T2 )
2
1
x3 = (3x + T3 )
4
1
x4 = (8x2 1 + T4 )
8
1
x5 = (20x3 5x + T5 )
16
1
x6 = (48x4 18x2 + 1 + T6 )
32
1
x7 = (112x5 56x3 + 7x + T7 )
64
1
x8 = (256x6 160x4 + 32x2 1 + T8 )
128
1
x9 = (576x7 432x5 + 120x3 9x + T9 )
256
1
x10 = (1280x8 1120x6 + 400x4 50x2 + 1 + T10 )
512
1
x11 = (2816x9 2816x7 + 1232x5 220x3 + 11x + T11 )
1024

For easy reference the formulas for economization of power series in terms of Chebyshev are
given in Table 4.

16
Assigned Problems

Problem Set for Chebyshev Polynomials

1. Obtain the first three Chebyshev polynomials T0 (x), T1 (x) and T2 (x) by means of
the Rodrigues formula.

2. Show that the Chebyshev polynomial T3 (x) is a solution of Chebyshevs equation of

order 3.

3. By means of the recurrence formula obtain Chebyshev polynomials T2 (x) and T3 (x)
given T0 (x) and T1 (x).

4. Show that Tn (1) = 1 and Tn (1) = (1)n

5. Show that Tn (0) = 0 if n is odd and (1)n/2 if n is even.

6. Setting x = cos show that

1 h p n  p n i
Tn (x) = x + i 1 x2 + x i 1 x2
2

where i = 1.

7. Find the general solution of Chebyshevs equation for n = 0.

8. Obtain a series expansion for f (x) = x2 in terms of Chebyshev polynomials Tn (x),

3
X
2
x = An Tn (x)
n=0

9. Express x4 as a sum of Chebyshev polynomials of the first kind.

17