How to obtain finite difference, which is continuous

Question

I want to calculate a finite difference (something like this SO Post). My data is as follows: I have x-values that are powers of two (4, 8, 16, 32 and 64). Corresponding to them are y-values, such that

$$Y=f(X)$$

where $f$ is a monotonic function. My interest is in calculating finite difference such as

$$\text{Diff}(x) = (y_2-y_1) /(x_2-x_1)$$ at different values of X, where X is in the domain (i.e. $4\le x \le 64$). Now my problem is that when I take piecewise linear approximation, Diff(x) changes a lot between 4-to-8 and 8-to-16. This is undesirable. If I try to fit quadratic equation, by taking three points at a time, and then calculate finite difference, it becomes negative at some points, which is not physically acceptable, since $f$ is monotonic function. In other words, I fitted $ax^2+bx+c$ and thus $\text{Diff}(X)=2ax+b$. So, the problem is that $2ax+b$ is not positive at all points in the domain.

Can someone point-out a method to solve this problem. I need to use C/C++ for solve this, I cannot use matlab. I just have five X values, so solution of any complexity is acceptable. Thanks for the help.

Answer 1

Instead of fitting a general polynomial, you way want to consider using piece-wise cubic splines to interpolate the function in the interval between each gridpoint. one good property of cubic splines is that they minimize the bending energy of the fitted curve, which may lead to fewer radical oscillations (which may be the cause of the negative values you experience with the quadratic approach you previously tried).

You can obtain a piecewise cubic spline by proposing a separate cubic polynomial of the form

$f_i(x)=a_ix^3+b_ix^2+c_ix+d_i$ in each interval $i$

whose endpoints are the gridpoints you selected. Then, you can impose the following conditions: 1. All functions must agree with the data at each of the nodes 2. The functions must be continuous at each of the node points 3. The first and second derivatives of each function must continuous at each of the node points

You can get two extra equations by imposing the so-called "natural spline" condition, that is the second derivative must be zero at the first and last node. This should give you enough information to solve for each of the coefficients of each cubic polynomial in each interval.

Answer 2

A simple and accurate approach in this case is to use standard differences in $\log$ space.

$$ \begin{align*} s(x) &= \log x \\ \partial_x f(x) &= \partial_s f(e^{s(x)}) \partial_x s(x) \\ &= \frac 1 x \partial_s f(e^{s(x)}) \\ \end{align*} $$

Since $s(x)$ is equally spaced, it is easy to approximate the derivative on the bottom line using any standard differencing.

Answer 3

Fitting by a parametrized family of functions (like degree $n$-polynomials, where there are $n+1$ free parameters corresponds to minimizing an error term on a parameter space. If the error term is a least squares approach, then this leads to a quadratic program.

Quadratic Programming

If you want your function to be non-negative, or to be a monotonic function (it is not clear what you mean), then you arrive a contrained minimization, which is another topic on its own.

Constraint programming

In your case, a low order polynomial should suffice. The first wikipedia article should lead to plenty of ways to solve your problem.