## Name:

remez computes the minimax of a function on an interval.

## Library names:

sollya_obj_t sollya_lib_remez(sollya_obj_t, sollya_obj_t, sollya_obj_t, ...) sollya_obj_t sollya_lib_v_remez(sollya_obj_t, sollya_obj_t, sollya_obj_t,                                 va_list)

## Usage:

remez(f, n, range, w, quality, bounds) : (function, integer, range, function, constant, range) -> function remez(f, L, range, w, quality, bounds) : (function, list, range, function, constant, range) -> function

## Parameters:

• f is the function to be approximated
• n is the degree of the polynomial that must approximate f
• L is a list of integers or a list of functions and indicates the basis for the approximation of f
• range is the interval where the function must be approximated
• w (optional) is a weight function. Default is 1.
• quality (optional) is a parameter that controls the quality of the returned polynomial p, with respect to the exact minimax p*. Default is 1e-5.
• bounds (optional) is a parameter that allows the user to make the algorithm stop earlier, whenever a given accuracy is reached or a given accuracy is proved unreachable. Default is [0, +Inf].

## Description:

• remez computes an approximation of the function f with respect to the weight function w on the interval range. More precisely, it searches p such that ||p*w-f|| is (almost) minimal among all p of a certain form. The norm is the infinity norm, e.g. ||g|| = max {|g(x)|, x in range}.
• If w=1 (the default case), it consists in searching the best polynomial approximation of f with respect to the absolute error. If f=1 and w is of the form 1/g, it consists in searching the best polynomial approximation of g with respect to the relative error.
• If n is given, p is searched among the polynomials with degree not greater than n. If L is given and is a list of integers, p is searched as a linear combination of monomials X^k where k belongs to L. In the case when L is a list of integers, it may contain ellipses but cannot be end-elliptic. If L is given and is a list of functions g_k, p is searched as a linear combination of the g_k. In that case L cannot contain ellipses. It is the user responsibility to check that the g_k are linearly independent over the interval range. Moreover, the functions w*g_k must be at least twice differentiable over range. If these conditions are not fulfilled, the algorithm might fail or even silently return a result as if it successfully found the minimax, though the returned p is not optimal.
• The polynomial is obtained by a convergent iteration called Remez' algorithm (and an extension of this algorithm, due to Stiefel). The algorithm computes a sequence p1, ..., pk, ... such that ek = ||pk*w-f|| converges towards the optimal value e. The algorithm is stopped when the relative error between ek and e is less than quality.
For this reason, the returned polynomial usually has dyadic coefficients. However, one should not rely on that property. Especially, a noticeable exception to that statement is when we manage to detect that f/w exactly simplifies to a polynomial of the required degree, in which case we return it as it is without further rounding.
• The optional argument bounds is an interval [satisfying_err, target_err] with the following behavior:
• if, during the algorithm, we manage to prove that target_err is unreachable, we stop the algorithm returning the last computed polynomial.
• if, during the algorithm, we obtain a polynomial with an error smaller than satisfying_err, we stop the algorithm returning that polynomial.
• otherwise we loop until we find an optimal polynomial with the required quality, as usual.
Examples of use:
[0, +Inf] (compute the optimal polynomial with the required quality)
[target_err] (stops as soon as a polynomial achieving target_err is obtained or as soon as such a polynomial is proved not to exist).
[0, target_err] (finds the optimal polynomial, but provided that its error is smaller than target_err).
[satisfying_err, +Inf] (stops as soon as a polynomial achieving satisfying_err is obtained. If such a polynomial does not exist, returns the optimal polynomial).

## Example 1:

> p = remez(exp(x),5,[0;1]);
> degree(p);
5
> dirtyinfnorm(p-exp(x),[0;1]);
1.1295698151096148707171193829266077607222634589363e-6

## Example 2:

> p = remez(1,[|0,2,4,6,8|],[0,Pi/4],1/cos(x));
> canonical=on!;
> p;
0.99999999994393732180959690352543887130348096061124 - 0.49999999571556857768772053063721544670949467222259 * x^2 + 4.1666613233473633009941059480570275870113220089059e-2 * x^4 - 1.3886529147145693651355523880319714051047635695061e-3 * x^6 + 2.4372679177224179934800328511009205218114284220126e-5 * x^8

## Example 3:

> p1 = remez(exp(x),5,[0;1],default,1e-5);
> p2 = remez(exp(x),5,[0;1],default,1e-10);
> p3 = remez(exp(x),5,[0;1],default,1e-15);
> dirtyinfnorm(p1-exp(x),[0;1]);
1.1295698151096148707171193829266077607222634589363e-6
> dirtyinfnorm(p2-exp(x),[0;1]);
1.12956980227478675612619255125474525171079325793124e-6
> dirtyinfnorm(p3-exp(x),[0;1]);
1.12956980227478675612619255125474525171079325793124e-6

## Example 4:

> L = [|exp(x), sin(x), cos(x)-1, sin(x^3)|];
> g = (2^x-1)/x;
> p1 = remez(g, L, [-1/16;1/16]);
> p2 = remez(g, 3, [-1/16;1/16]);
> dirtyinfnorm(p1 - g, [-1/16;1/16]);
9.8841323829271038137685646777951687620288462194746e-8
> dirtyinfnorm(p2 - g, [-1/16;1/16]);
2.54337800593461418356437401152248866818783932027105e-9

## Example 5:

> f = sin(x);
> I = [-3b-5;-1b-1074];
> time(popt = remez(1, [|1, 3, 4, 5, 7, 8, 9|], I, 1/f));
0.185799695000000000000000000000000000002291314966024
> time(p1 = remez(1, [|1, 3, 4, 5, 7, 8, 9|], I, 1/f, default, [0, 1b-73]));
0.134645933000000000000000000000000000002249871361033
> time(p2 = remez(1, [|1, 3, 4, 5, 7, 8, 9|], I, 1/f, default, [3b-72, +@Inf@]));
0.150942799000000000000000000000000000010601473971874
> dirtyinfnorm(popt/f-1, I);
2.06750931454112835098093903810531156576504665659064e-22
> dirtyinfnorm(p1/f-1, I);
2.49711266837493110470637913808914046704452778960875e-22
> dirtyinfnorm(p2/f-1, I);
5.4567247553615435246376977231253834265248756996947e-22
> 1b-73;
1.05879118406787542383540312584955245256423950195312e-22
> 3b-72;
6.3527471044072525430124187550973147153854370117187e-22
Go back to the list of commands