Orthogonal Functions∗ May 9, 2003
1
Introduction
A particularly important class of applications of linear algebra within mathematics arises from treating functions as vectors, with the rules • (f + g)(x) = f (x) + g(x) • (αf )(x) = αf (x), where α is a scalar. Rb • f · g = a f (x)g(x) dx. Here, the interval [a, b] over which we integrate is regarded as fixed within any particular discussion. Two important examples are 1. [a, b] = [0, 1] and the functions in question are taken to be polynomials. 2. [a, b] = [0, 2π] and the functions are periodic, e.g. trigonometric polynomials. A set of functions on [a,b] is linearly independent if no linear combination of those functions is identically 0 on [a,b]. In particular, the powers of x are linearly independent on any interval, as are the functions sin(nx), cos(nx) and 1 on [0, 2π], for n and m positive integers. However, the latter set has an additional important property: they are orthogonal. That is Z 2π Z 2π Z 2π sin(nx) cos(mx) dx = sin(nx) sin(mx) dx = cos(nx) cos(mx) dx = 0 0
0 ∗
0
Preliminary version: subject to later expansion
1
unless m = n, in which case Z 2π sin(nx) cos(nx) dx = 0;
Z
2π 2
2π
Z
cos2 (nx) dx = π.
sin (nx) dx =
0
0
0
The powers of x are not orthogonal on any interval. However if we confine our attention to any particular interval, such as [0, 1], we can use the GramSchmidt orthogonalization algorithm to produce orthogonal polynomials. PnSuppose we have orthogonal functions {fi }0≤i≤n , and a function g = i=0 ai fi , which is a linear combination of the functions fi . Then Z
b
g(x)fj (x) dx = a
n X
b
Z ai
Z fi (x)fj (x) dx = aj
a
i=0
so that
b
fj (x)2 dx,
a
Rb aj =
g(x)fj (x) dx , Rb 2 dx f (x) j a
a
and the coefficients aj can be recovered by integration.
2
Complex Exponentials and the Discrete Fourier Transform
In dealing with periodic functions on [0, 2π], it is particularly useful to consider the complex exponential functions einx as n varies over the integers. From the power series x
e =
∞ X xn 0
n
=1+x+
cos(x) = 1 −
x2 x3 + + ··· , 2 3
x2 x4 x6 + − + ··· , 2 4 6
and
x3 x5 x7 sin(x) = x − + − + ··· , 3 5 7 we recover the identity eix = cos(x) + i sin(x), 2
which is also consistent with the trigonometric addition formulae. Because we are considering functions that take complex values (although only as functions of the real variable x), we are, for the first time, considering a complex vector space. This means that scalars are allowed to be complex numbers, and also requires an adjustment to the dot product. We recall that if z = x + iy, then z¯ = x − iy, and set Z 2π f ·g = f (x)g(x) dx. 0
With this definition, einx · eimx =P2π if n = m, and 0 otherwise, since = e−inx . In particular, if f (x) = an einx , then Z 2π imx f ·e = f (x)e−imx = 2πam ,
einx
0
and we can recover the coefficients am by integration, as in the strictly real case. Let us now consider how to best approximate the coefficients am for 0 ≤ |m| ≤ N − 1 when we know the values f ( 2πk ) for a real valued function f . N From the fact that f is real values, it follows that a−m = am for all m, so that we need only determine the coefficients am for m ≥ 0. The next observation is that we can approximate the integral that gives us 2πam by the Riemann sum N −1 2π X 2πk −2πimk f( )e N , N k=0 N so that
N −1 N −1 1 X 2πk 1 X 2πk −2πimk am ≈ f( )e N = f( )WNmk , N k=0 N N k=0 N −2πi
where we define WN to be e N . If we set v to be the vector whose k t h component is f ( 2π(k−1) ) and MN to be the complex N × N symmetric matrix N (m−1)(k−1) whose (m, k)th entry is WN , the components of MN v are approximations to 2πN am 0 ≤ m ≤ N − 1. This approximation is called the discrete Fourier transform. It tends to be rather crude for values of m close to N .
3
3
The Fast Fourier Transform
In order to appreciate what follows, one must think of N as very large (perhaps as large as 106 ), so that the N 2 multiplications that are necessary to compute the discrete Fourier transform will take significant time. We will see now that there is a trick that allows us to do the computation with many fewer multiplications if N is a power of 2. The assumption is that multiplication is much more time consuming than addition, so that it is much more efficient to compute a(b + c + d) then ab + ac + ad. It will be more appropriate for what follows to shift the components and indices so that they go from 0 to N − 1, vk = f ( 2π(k ) and the (m, k)th N (m)(k) component of MN is WN . We make the observation that, for N even, N
WN2 = W N , and WN2 = −1. We write veven for the vector formed from the 2 even components of v taken in order, and vodd for the vector formed from the odd components of v, also taken in order. From this it follows that the first N components of MN v are the components of the vector of length N2 whose 2 j t h component, (ending with j = N2 − 1) is (M N veven )j + WNj (M N vodd )j , 2
while the (j +
N th ) 2
2
component of MN v (starting with j = 0) is (M N veven )j − WNj (M N vodd )j . 2
2
If we write µ(N ) for the minimum number of multiplications required to compute MN v, it follows from the foregoing that, for N even, µ(N ) ≤ 2µ(
N N )+ . 2 2
We are not counting the multiplications required to compute the powers of WN , because these can be stored and used repeatedly, and we are really interested in the speed for which the computation can be carried out for fixed N and constantly changing v. Since M2 =
1 1 , 1 −1
it follows that µ(2) ≤ 1 and, from the foregoing by induction, that µ(2p ) ≤ p 2p−1 . 4
4
Problems 1. Let f (x) be the function on [0, 1] defined by f (x) = x,
1 0≤x≤ ; 2
f (x) = 1 − x,
1 ≤ x ≤ 1. 2
Find the orthogonal projection of f on the space of quartic polynomials. Begin by finding an orthogonal basis for the space of quartic polynomials by applying Gram-Schmidt to the functions {1, x, x2 , x3 , x4 }. You are encouraged to use MATLAB for this purpose. 2. Let f (x) = π − |π − x| for 0 ≤ x ≤ 2π. (a) Find the Fourier coefficients am = direct integration for |m| ≤ 4.
1 2π
R 2π 0
f (x)e−imx dx of f by
(b) Use the values of f at multiples of π2 and the discrete Fourier transform to approximate the same coefficients.
5