In mathematics, approximation theory is concerned with how functions can best be approximated with simpler functions, and with quantitatively characterizing the errors introduced thereby. What is meant by best and simpler will depend on the application.
A closely related topic is the approximation of functions by generalized Fourier series, that is, approximations based upon summation of a series of terms based upon orthogonal polynomials.
One problem of particular interest is that of approximating a function in a computer mathematical library, using operations that can be performed on the computer or calculator (e.g. addition and multiplication), such that the result is as close to the actual function as possible. This is typically done with polynomial or rational (ratio of polynomials) approximations.
The objective is to make the approximation as close as possible to the actual function, typically with an accuracy close to that of the underlying computer's floating point arithmetic. This is accomplished by using a polynomial of high degree, and/or narrowing the domain over which the polynomial has to approximate the function.Narrowing the domain can often be done through the use of various addition or scaling formulas for the function being approximated. Modern mathematical libraries often reduce the domain into many tiny segments and use a low-degree polynomial for each segment.
Once the domain (typically an interval) and degree of the polynomial are chosen, the polynomial itself is chosen in such a way as to minimize the worst-case error. That is, the goal is to minimize the maximum value of
\midP(x)-f(x)\mid
+\varepsilon
-\varepsilon
\varepsilon
For example, the graphs shown to the right show the error in approximating log(x) and exp(x) for N = 4. The red curves, for the optimal polynomial, are level, that is, they oscillate between
+\varepsilon
-\varepsilon
To prove this is true in general, suppose P is a polynomial of degree N having the property described, that is, it gives rise to an error function that has N + 2 extrema, of alternating signs and equal magnitudes. The red graph to the right shows what this error function might look like for N = 4. Suppose Q(x) (whose error function is shown in blue to the right) is another N-degree polynomial that is a better approximation to f than P. In particular, Q is closer to f than P for each value xi where an extreme of P−f occurs, so
|Q(xi)-f(xi)|<|P(xi)-f(xi)|.
Q(xi)-f(xi)\le|Q(xi)-f(xi)|<|P(xi)-f(xi)|=P(xi)-f(xi),
f(xi)-Q(xi)\le|Q(xi)-f(xi)|<|P(xi)-f(xi)|=f(xi)-P(xi).
One can obtain polynomials very close to the optimal one by expanding the given function in terms of Chebyshev polynomials and then cutting off the expansion at the desired degree.This is similar to the Fourier analysis of the function, using the Chebyshev polynomials instead of the usual trigonometric functions.
If one calculates the coefficients in the Chebyshev expansion for a function:
f(x)\sim
infty | |
\sum | |
i=0 |
ciTi(x)
and then cuts off the series after the
TN
The reason this polynomial is nearly optimal is that, for functions with rapidly converging power series, if the series is cut off after some term, the total error arising from the cutoff is close to the first term after the cutoff. That is, the first term after the cutoff dominates all later terms. The same is true if the expansion is in terms of bucking polynomials. If a Chebyshev expansion is cut off after
TN
TN+1
TN+1
TN
In the graphs above, the blue error function is sometimes better than (inside of) the red function, but sometimes worse, meaning that it is not quite the optimal polynomial. The discrepancy is less serious for the exp function, which has an extremely rapidly converging power series, than for the log function.
Chebyshev approximation is the basis for Clenshaw–Curtis quadrature, a numerical integration technique.
The Remez algorithm (sometimes spelled Remes) is used to produce an optimal polynomial P(x) approximating a given function f(x) over a given interval. It is an iterative algorithm that converges to a polynomial that has an error function with N+2 level extrema. By the theorem above, that polynomial is optimal.
Remez's algorithm uses the fact that one can construct an N
th-degree polynomial that leads to level and alternating error values, given N+2 test points.Given N+2 test points
x1
x2
xN+2
x1
xN+2
\begin{align} P(x1)-f(x1)&=+\varepsilon\\ P(x2)-f(x2)&=-\varepsilon\\ P(x3)-f(x3)&=+\varepsilon\\ & \vdots\\ P(xN+2)-f(xN+2)&=\pm\varepsilon. \end{align}
The right-hand sides alternate in sign.
That is,
\begin{align} P0+P1x1+P2
2 | |
x | |
1 |
+P3
3 | |
x | |
1 |
+...+PN
N | |
x | |
1 |
-f(x1)&=+\varepsilon\\ P0+P1x2+P2
2 | |
x | |
2 |
+P3
3 | |
x | |
2 |
+...+PN
N | |
x | |
2 |
-f(x2)&=-\varepsilon\\ & \vdots \end{align}
Since
x1
xN+2
f(x1)
f(xN+2)
P0
P1
PN
\varepsilon
x1
xN+2
\varepsilon
The graph below shows an example of this, producing a fourth-degree polynomial approximating
ex
\varepsilon
The error graph does indeed take on the values
\pm\varepsilon
The second step of Remez's algorithm consists of moving the test points to the approximate locations where the error function had its actual local maxima or minima. For example, one can tell from looking at the graph that the point at −0.1 should have been at about −0.28. The way to do this in the algorithm is to use a single round of Newton's method. Since one knows the first and second derivatives of, one can calculate approximately how far a test point has to be moved so that the derivative will be zero.
Calculating the derivatives of a polynomial is straightforward. One must also be able to calculate the first and second derivatives of f(x). Remez's algorithm requires an ability to calculate
f(x)
f'(x)
f''(x)
After moving the test points, the linear equation part is repeated, getting a new polynomial, and Newton's method is used again to move the test points again. This sequence is continued until the result converges to the desired accuracy. The algorithm converges very rapidly. Convergence is quadratic for well-behaved functions—if the test points are within
10-15
10-30
Remez's algorithm is typically started by choosing the extrema of the Chebyshev polynomial
TN+1