Tutorial 3 of 5: changing variables for simpler polynomials

In the previous section, we introduced function $g(x)$ and saw that lolremez looked for a polynomial $P(x)$ and the smallest error value $E$ such that:

$$\max_{x \in [-1, 1]}{\dfrac{\big\lvert f(x) - P(x)\big\rvert}{\big\lvert g(x)\big\rvert}} = E$$

Where $g(x)$ is a function used as a weight to minimise relative error.

We will now see that $g(x)$ can be used for other purposes.

Odd and even functions

An odd function is such that in all points, $f(-x) = -f(x)$. Similarly, an even function is such that $f(-x) = f(x)$.

An example of an odd function is $\sin(x)$. An example of an even function is $\cos(x)$.

A reasonable expectation is that approximating an odd function over an interval centred at zero will lead to an odd polynomial. But due to the way the Remez exchange algorithm works, this is not necessarily the case. We can help the algorithm a little.

The example of $\sin(x)$

Suppose we want to approximate $\sin(x)$ over $[-\pi/2, \pi/2]$:

$$\max_{x \in [-\pi/2, \pi/2]}{\big\vert\sin(x) - P(x)\big\vert} = E$$

We know $\sin(x)$ is an odd function, so instead we look for an odd polynomial. Odd polynomials are polynomials that only use odd powers of $x$, such as $x$, $x^3$, $x^5$…

We notice that an odd polynomial $P(x)$ can always be rewritten as $xQ(x^2)$ where $Q$ has the same coefficients. Our minimax constraint becomes:

$$\max_{x \in [-\pi/2, \pi/2]}{\big\vert\sin(x) - xQ(x^2)\big\vert} = E$$

Another observation is that replacing $x$ with $-x$ has no effect on the expression in the left. We can therefore restrict the range to $[0, \pi/2]$:

$$\max_{x \in [0, \pi/2]}{\big\vert\sin(x) - xQ(x^2)\big\vert} = E$$

Now since $x$ is always positive, we can at last perform our change of variable, replacing $x^2$ with $y$ (and $x$ with $\sqrt{y}$, since $x$ is positive), and changing the range accordingly:

$$\max_{y \in [0, \pi^2/4]}{\big\lvert\sin(\sqrt{y}) - \sqrt{y}Q(y)\big\rvert} = E$$

This is almost something lolremez can deal with, if only there wasn’t this $\sqrt{y}$ in front of $Q(y)$. Fortunately we can simply divide everything by $\sqrt{y}$:

$$\max_{y \in [0, \pi^2/4]}{\dfrac{\bigg\lvert\dfrac{\sin(\sqrt{y})}{\sqrt{y}} - Q(y)\bigg\rvert}{\dfrac{1}{|\sqrt{y}|}}} = E$$

That’s it! We have our minimax equation, and we can feed it to lolremez.

Invocation

lolremez --degree 4 --range 1e-50:pi*pi/4 "sin(sqrt(x))/sqrt(x)" "1/sqrt(x)"

Note: since $f$ is not defined at zero, we simply reduce the range to $[1e-50, \pi^2/4]$ instead of $[0, \pi^2/4]$.

After all the iterations the output should be as follows:

/* Approximation of f(x) = sin(sqrt(x))/sqrt(x)
 * with weight function g(x) = 1/sqrt(x)
 * on interval [ 9.9999999999999999999e-51, 2.4674011002723396547 ]
 * with a polynomial of degree 4. */
float f(float x)
{
    float u = 2.4609388329975758276e-6f;
    u = u * x + -1.9698112762578435577e-4f;
    u = u * x + 8.3298294559966612167e-3f;
    u = u * x + -1.6666236485125293496e-1f;
    return u * x + 9.9999788400553332261e-1f;
}

We can therefore write the corresponding C++ function:

float fastsin(float x)
{
    return x * f(x * x);
}

Analysing the results

The obtained polynomial needs 5 constants for a maximum error of about 3.3381e-9 over $[-\pi/2, \pi/2]$. The real 9th degree minimax polynomial on $[-\pi/2, \pi/2]$ needs 10 constants for a maximum error of about… 3.3381e-9!

Again, the error curves match almost perfectly. It means that the odd coefficients in the minimax approximation of an even function play very little role.

Conclusion

You should now be able to give hints to the Remez solver through changes of variables, and find polynomials with fewer constants!

Please report any trouble you may have had with this document to [email protected]. You may then carry on to the next section: fixing lower-order parameters.