165 lines
9.1 KiB
Plaintext
165 lines
9.1 KiB
Plaintext
[section:minimax Minimax Approximations and the Remez Algorithm]
|
|
|
|
The directory `libs/math/minimax` contains an interactive command-line driven
|
|
program for the generation of minimax approximations using the Remez
|
|
algorithm. Both polynomial and rational approximations are supported,
|
|
although the latter are tricky to converge: it is not uncommon for
|
|
convergence of rational forms to fail. No such limitations are present
|
|
for polynomial approximations which should always converge smoothly.
|
|
|
|
It's worth stressing that developing rational approximations to functions
|
|
is often not an easy task, and one to which many books have been devoted.
|
|
To use this tool, you will need to have a reasonable grasp of what the Remez
|
|
algorithm is, and the general form of the approximation you want to achieve.
|
|
|
|
Unless you already familar with the Remez method, you should first read the
|
|
[link math_toolkit.remez brief background article explaining the principles behind the Remez algorithm].
|
|
|
|
The program consists of two parts:
|
|
|
|
[variablelist
|
|
[[main.cpp][Contains the command line parser, and all the calls to the Remez code.]]
|
|
[[f.cpp][Contains the function to approximate.]]
|
|
]
|
|
|
|
Therefore to use this tool, you must modify f.cpp to return the function to
|
|
approximate. The tools supports multiple function approximations within
|
|
the same compiled program: each as a separate variant:
|
|
|
|
NTL::RR f(const NTL::RR& x, int variant);
|
|
|
|
Returns the value of the function /variant/ at point /x/. So if you
|
|
wish you can just add the function to approximate as a new variant
|
|
after the existing examples.
|
|
|
|
In addition to those two files, the program needs to be linked to
|
|
a [link math_toolkit.high_precision.use_ntl patched NTL library to compile].
|
|
|
|
Note that the function /f/ must return the rational part of the
|
|
approximation: for example if you are approximating a function
|
|
/f(x)/ then it is quite common to use:
|
|
|
|
[expression f(x) = g(x)(Y + R(x))]
|
|
|
|
where /g(x)/ is the dominant part of /f(x)/, /Y/ is some constant, and
|
|
/R(x)/ is the rational approximation part, usually optimised for a low
|
|
absolute error compared to |Y|.
|
|
|
|
In this case you would define /f/ to return [role serif-italic f(x)/g(x)] and then set the
|
|
y-offset of the approximation to /Y/ (see command line options below).
|
|
|
|
Many other forms are possible, but in all cases the objective is to
|
|
split /f(x)/ into a dominant part that you can evaluate easily using
|
|
standard math functions, and a smooth and slowly changing rational approximation
|
|
part. Refer to your favourite textbook for more examples.
|
|
|
|
Command line options for the program are as follows:
|
|
|
|
[variablelist
|
|
[[variant N][Sets the current function variant to N. This allows multiple functions
|
|
that are to be approximated to be compiled into the same executable.
|
|
Defaults to 0.]]
|
|
[[range a b][Sets the domain for the approximation to the range \[a,b\], defaults
|
|
to \[0,1\].]]
|
|
[[relative][Sets the Remez code to optimise for relative error. This is the default
|
|
at program startup. Note that relative error can only be used
|
|
if f(x) has no roots over the range being optimised.]]
|
|
[[absolute][Sets the Remez code to optimise for absolute error.]]
|
|
[[pin \[true|false\]]["Pins" the code so that the rational approximation
|
|
passes through the origin. Obviously only set this to
|
|
/true/ if R(0) must be zero. This is typically used when
|
|
trying to preserve a root at \[0,0\] while also optimising
|
|
for relative error.]]
|
|
[[order N D][Sets the order of the approximation to /N/ in the numerator and /D/
|
|
in the denominator. If /D/ is zero then the result will be a polynomial
|
|
approximation. There will be N+D+2 coefficients in total, the first
|
|
coefficient of the numerator is zero if /pin/ was set to true, and the
|
|
first coefficient of the denominator is always one.]]
|
|
[[working-precision N][Sets the working precision of NTL::RR to /N/ binary digits. Defaults to 250.]]
|
|
[[target-precision N][Sets the precision of printed output to /N/ binary digits:
|
|
set to the same number of digits as the type that will be used to
|
|
evaluate the approximation. Defaults to 53 (for double precision).]]
|
|
[[skew val]["Skews" the initial interpolated control points towards one
|
|
end or the other of the range. Positive values skew the
|
|
initial control points towards the left hand side of the
|
|
range, and negative values towards the right hand side.
|
|
If an approximation won't converge (a common situation)
|
|
try adjusting the skew parameter until the first step yields
|
|
the smallest possible error. /val/ should be in the range
|
|
\[-100,+100\], the default is zero.]]
|
|
[[brake val][Sets a brake on each step so that the change in the
|
|
control points is braked by /val%/. Defaults to 50,
|
|
try a higher value if an approximation won't converge,
|
|
or a lower value to get speedier convergence.]]
|
|
[[x-offset val][Sets the x-offset to /val/: the approximation will
|
|
be generated for `f(S * (x + X)) + Y` where /X/ is the
|
|
x-offset, /S/ is the x-scale
|
|
and /Y/ is the y-offset. Defaults to zero. To avoid
|
|
rounding errors, take care to specify a value that can
|
|
be exactly represented as a floating point number.]]
|
|
[[x-scale val][Sets the x-scale to /val/: the approximation will
|
|
be generated for `f(S * (x + X)) + Y` where /S/ is the
|
|
x-scale, /X/ is the x-offset
|
|
and /Y/ is the y-offset. Defaults to one. To avoid
|
|
rounding errors, take care to specify a value that can
|
|
be exactly represented as a floating point number.]]
|
|
[[y-offset val][Sets the y-offset to /val/: the approximation will
|
|
be generated for `f(S * (x + X)) + Y` where /X/
|
|
is the x-offset, /S/ is the x-scale
|
|
and /Y/ is the y-offset. Defaults to zero. To avoid
|
|
rounding errors, take care to specify a value that can
|
|
be exactly represented as a floating point number.]]
|
|
[[y-offset auto][Sets the y-offset to the average value of f(x)
|
|
evaluated at the two endpoints of the range plus the midpoint
|
|
of the range. The calculated value is deliberately truncated
|
|
to /float/ precision (and should be stored as a /float/
|
|
in your code). The approximation will
|
|
be generated for `f(x + X) + Y` where /X/ is the x-offset
|
|
and /Y/ is the y-offset. Defaults to zero.]]
|
|
[[graph N][Prints N evaluations of f(x) at evenly spaced points over the
|
|
range being optimised. If unspecified then /N/ defaults
|
|
to 3. Use to check that f(x) is indeed smooth over the range
|
|
of interest.]]
|
|
[[step N][Performs /N/ steps, or one step if /N/ is unspecified.
|
|
After each step prints: the peek error at the extrema of
|
|
the error function of the approximation,
|
|
the theoretical error term solved for on the last step,
|
|
and the maximum relative change in the location of the
|
|
Chebyshev control points. The approximation is converged on the
|
|
minimax solution when the two error terms are (approximately)
|
|
equal, and the change in the control points has decreased to
|
|
a suitably small value.]]
|
|
[[test \[float|double|long\]][Tests the current approximation at float,
|
|
double, or long double precision. Useful to check for rounding
|
|
errors in evaluating the approximation at fixed precision.
|
|
Tests are conducted at the extrema of the error function of the
|
|
approximation, and at the zeros of the error function.]]
|
|
[[test \[float|double|long\] N] [Tests the current approximation at float,
|
|
double, or long double precision. Useful to check for rounding
|
|
errors in evaluating the approximation at fixed precision.
|
|
Tests are conducted at N evenly spaced points over the range
|
|
of the approximation. If none of \[float|double|long\] are specified
|
|
then tests using NTL::RR, this can be used to obtain the error
|
|
function of the approximation.]]
|
|
[[rescale a b][Takes the current Chebeshev control points, and rescales them
|
|
over a new interval \[a,b\]. Sometimes this can be used to obtain
|
|
starting control points for an approximation that can not otherwise be
|
|
converged.]]
|
|
[[rotate][Moves one term from the numerator to the denominator, but keeps the
|
|
Chebyshev control points the same. Sometimes this can be used to obtain
|
|
starting control points for an approximation that can not otherwise be
|
|
converged.]]
|
|
[[info][Prints out the current approximation: the location of the zeros of the
|
|
error function, the location of the Chebyshev control points, the
|
|
x and y offsets, and of course the coefficients of the polynomials.]]
|
|
]
|
|
|
|
[endsect] [/section:minimax Minimax Approximations and the Remez Algorithm]
|
|
|
|
[/
|
|
Copyright 2006 John Maddock and Paul A. Bristow.
|
|
Distributed under the Boost Software License, Version 1.0.
|
|
(See accompanying file LICENSE_1_0.txt or copy at
|
|
http://www.boost.org/LICENSE_1_0.txt).
|
|
]
|