Version 10 (modified by sam, 9 years ago) (diff)

--

Remez tutorial 1/5: exp(x) the quick way

This is a hands-on example of the Lol Remez toolkit.

In this section we are going to approximate the exp(x) function using a polynomial.

Getting started

If you do not have the full Lol Engine source code, download and unpack the latest LolRemez tarball.

The file you should edit is remez.cpp.

What does Remez do?

Given a function f and a range [a,b], the Remez algorithm looks for the polynomial P(x) that minimises the error value E:

\[\max_{x \in [a,b]}{\big\vert f(x) - P(x)\big\vert} = E\]

We will see an example right now.

Source code

#include "lol/math/real.h"
#include "lol/math/remez.h"

using lol::real;
using lol::RemezSolver;

real f(real const &x) { return exp(x); }

int main(int argc, char **argv)
{
    RemezSolver<4, real> solver;
    solver.Run(-1, 1, f, 30);
    return 0;
}

What does this mean?

  • we declare function f which returns the exponential of x: this is the function we want to approximate.
  • we create a RemezSolver object for 4th-degree polynomials and real numbers.
  • we run the solver on the [-1,1] range, approximating function f for 30 iterations.

Compilation

If you are using LolRemez, just put the above source code in remez.cpp and type:

make

Execution

To launch the test, type:

./remez

After all the iterations the output should be as follows:

Final error: 5.462771976237482581009771665937582411463e-4
Polynomial estimate:
x**0*1.000090756764725753887362987792025308996
+x**1*9.973086551667860566788019540269306006270e-1
+x**2*4.988332174505582284710918757571761729419e-1
+x**3*1.773462612793916519454714108029230813767e-1
+x**4*4.415666059995979611944324860870682575219e-2

Using the results

The above results can be used in a more CPU-friendly implementation such as the following one:

double fastexp(double x)
{
    const double a0 = 1.000090756764725753887362987792025308996;
    const double a1 = 9.973086551667860566788019540269306006270e-1;
    const double a2 = 4.988332174505582284710918757571761729419e-1;
    const double a3 = 1.773462612793916519454714108029230813767e-1;
    const double a4 = 4.415666059995979611944324860870682575219e-2;

    return a0 + x * (a1 + x * (a2 + x * (a3 + x * a4)));
}

Analysing the results

Plotting the real exponential function and our fastexp function gives the following curves:

The curves are undistinguishable. Actually they differ by no more than 5.462772e-4, which is the value the ./remez output gave.

It can be verified on the following error curve:

Conclusion

You should now be all set up for your own minimax polynomial computation!

Please report any trouble you may have had with this document to sam@hocevar.net. You may then carry on to the next section: switching to relative error.

Attachments (2)

Download all attachments as: .zip