I need to approximate $x \ln x$ on $[0,1]$ as a piecewise-linear function. If $P(x)$ is a piecewise-linear approximation, I want to minimize $$ \max_{0 \le x \le 1} |P(x) - x \ln x| \rightarrow \min_P. $$

I tried interpolation over evenly spread points $x_i = \frac{i}{n}$, for $i=0,1,\dotsc,n$. I also tried Chebyshev nodes together with additional nodes $x=0$ and $x=1$. The latter works much better. For my numerical problem, it was enough but I am interested in optimal solution (mathematical mind, you know :)

Is there a way to optimally spread points $x_0, x_1, \dotsc, x_n$ over $[0,1]$ for piecewise-linear interpolation/approximation?

Update To illustrate the difference between uniform and Chebyshev nodes, here is the error of approximation with $n=20$ points.

Uniform vs Chebyshev nodes

My guess is that for optimal nodes, the "bumps" should be of the same height.

Update 2

I have also the following idea of a numerical method.

First, let us see what is the maximum error of approximation between points $x_n$ and $x_{n+1}$.

error of approximation

We can find equation of an approximating (red) line $y = kx + b$, where $$ k = \frac{x_{n+1} \ln x_{n+1} - x_n \ln x_n}{x_{n+1} - x_n}, \quad b = -k x_n+x_n \ln x_n \,. $$

Next, for $x_n \le x \le x_{n+1}$, we can define the error of approximation at point $x$: $$ \delta(x) = k x + b - x \ln x \,. $$ To find the maximum error, we just solve $\delta'(x) = 0$: $$ 0 = (kx + b - x \ln x)' = k - \ln x - 1. $$ Therefore, the maximum error on $[x_n,x_{n+1}]$ is achieved at $x_c=e^{k-1}$ and it is equal to $$ d(x_n, x_{n+1}) = k x_c + b - x_c \ln x_c \,, $$ where $k$, $b$, and $x_c$ depend on $x_n$ and $x_{n+1}$ as described above. Note also that if we fix $x_n$, the $d(x_n, x_{n+1})$ is increasing function in $x_{n+1}$.

Now, the method itself. Assume we want to ensure that approximation error is not larger than $\epsilon$. Then the following iterative procedure can be applied to construct the list of points.

  1. Set $x_0 = 0$.
  2. For $n=0,1,2,\dotsc$ find (e.g. by binary search) $x_{n+1}$ that ensures $d(x_n, x_{n+1}) = \epsilon$.
  3. Stop when some $x_{n+1} \ge 1$. Change $x_{n+1}$ to $x_{n+1} = 1$ (not necessary but to make things nicer).

But is there an analytical solution?

  • $\begingroup$ For example, if we set three points $x_0=0, x_1, x_2=1$, then it seems the optimal position for $x_1$ is $0.2869781560930248$... $\endgroup$ – Yauhen Yakimenka Jan 29 '20 at 10:25

There is an analytical solution to the problem in the following sense:

Given a number $N$, the optimal interpolation points $x_0=0, x_1, ..., x_{N-1}=1$ are the roots of an $(N-2) \times (N-2)$ system of non-linear equations described below.

The main insight to construct this system is the following proposition (using the same notation as in the question):

The interpolation points that minimize the maximal error are achieved when each interval $[x_i, x_{i+1}]$ attains exactly the same maximal error $d(x_i, x_{i+1}) = k_i x_{c_i} + b_i - x_{c_i} \ln x_{c_i}$.

This can be proved using the fact that, since $x \ln x$ is convex, the maximal error on any interval can be reduced by reducing the length of the interval. Thus, if we assume that the maximal interpolation error is attained at some interval $[x_i, x_{i+1}]$, then if the error on the neighboring interval $[x_{i-1}, x_i]$ was smaller, we could have moved $x_i$ to the right and reduce the maximal interpolation error. Therefore the error on $[x_{i-1}, x_i]$ cannot be smaller and the same argument can be made further to prove that all intervals attain the maximal error.

So, the minimal maximal error is attained when all the intervals have the same maximal error.

Using this proposition, and the $d(.)$ function developed in the question, we can build a set of non-linear equations of the form:

$$ d(0, x_1) = d(x_1, x_2) \\ ...\\ d(0, x_1) = d(x_{N-2}, 1) $$

and use any non-linear solver to solve it (in the results below I used python's scipy.optimize.fsolve(), with an initial guess of equally spaced $x$-values).

For $N=3$, we get a single equation: $$ d(0, x_1) = d(x_1, 1) $$

and the root is (as @yauhen-yakimenka notes in his comment) $x_1 = 0.28697816$ with an optimal error of $0.10557336369191296$.

In the error plot below, we can see this optimal value attained in two places.

enter image description here

For $N=20$, we get the roots at:

$$ 0, 0.00369174, 0.01286418, 0.02747082, 0.04751028, 0.0729823, 0.10388683, 0.14022384, 0.18199332, 0.22919527, 0.28182968, 0.33989654, 0.40339587, 0.47232766, 0.5466919, 0.62648861y, 0.71171777, 0.80237939, 0.89847347, 1 $$

With the optimal error of $0.001358114333332766$, see error plot below.

enter image description here


Each piece of deviation from a linear interpolation resembles a parabola. So if we interpolate a function $f$ over a segment of width $w$, then the maximum deviation of the curve from the segment is roughly $(w/2)^2 f''/2$. To get a maximum deviation of $\epsilon$, we should choose $w=\sqrt{8\epsilon/f''}$. To be more precise, we should use a value for $f''$ in the middle of the parabola, whose location we can estimate from the previous two points of interpolation.

This leads to the following procedure: Let \begin{align} x_0 &=1 + 1/N\\ x_1 &= 1\\ x_{n+1} &= x_n - \sqrt{8\epsilon\big/ f''\!\left(x_n+\frac{x_n - x_{n-1}}2\right)}\\ \end{align} Then we interpolate with the points $\{0,x_{N-1},\ldots,x_1\}$.

Here $N=20$, $f(x)=x \log x$, and by trial and error we choose $\epsilon=.00141$. This gives the points $$\{0.0000, 0.0039, 0.0112, 0.0242, 0.0427, 0.0670, 0.0968, 0.1323, 0.1735, 0.2203, \\ \ \ 0.2727, 0.3308, 0.3945, 0.4638, 0.5388, 0.6195, 0.7057, 0.7976, 0.8951, 1.0000\}.$$ Linear interpolation on those points approximates $f$ with a maximum deviation of only $.00145$.

  • $\begingroup$ Hmm... actually, I hoped for the analytical solution :) I think I have an iterative method that is exact... I will test it and post it here if it works $\endgroup$ – Yauhen Yakimenka Jan 29 '20 at 20:25
  • $\begingroup$ See Update 2 in the question. What do you think? $\endgroup$ – Yauhen Yakimenka Jan 29 '20 at 21:34
  • 1
    $\begingroup$ If I set $\epsilon=0.00145$ in my numerical method, I also get 20 points: {0.000000,0.003942,0.013735,0.029329,0.050725,0.077920,0.110915,0.149711,0.194307,0.244702,0.300898,0.362893,0.430689,0.504284,0.583680,0.668875,0.759871,0.856666,0.959261,1}. So your method is rather good (at least, for this $\epsilon$). $\endgroup$ – Yauhen Yakimenka Jan 29 '20 at 21:47
  • $\begingroup$ For most $f$, there is no closed form for the maximum error on an interval, so I like having a method that doesn't depend on numerical maximization at each of those 20 steps. $\endgroup$ – Matt F. Jan 29 '20 at 22:11
  • $\begingroup$ True, it requires solving $-k = f'(x)$ $\endgroup$ – Yauhen Yakimenka Jan 29 '20 at 22:18

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.