A single function for implementing Newton forward and backward polynomial interpolation

Question

On various websites, I've come across a lot of realizations of Newton polynomial interpolation that use separate functions for forward and backward interpolation, respectively. Such solutions seem ugly to me because the functions look pretty much the same.

For a single Newton interpolation function, I treat the backward interpolation case as forward interpolation with the data points written in reverse order. Usually the data points come already sorted in ascending order by x values, so the reverse order means simply going through the data points from right to left over an interval.

With all these in mind, my Newton interpolation function in C is:

/*
Error codes:
-2 - malloc() error
-1 - range check error
 0 - no error (success)
*/

// The necessary typedefs

typedef struct {
    int points;
    double *xList, *yList;
} DataTable;

typedef struct {
    double *coeff;
    int degree;
    double *xList;
} NewtonPoly;

// I've excluded the code for checking if the data table is filled properly.

int getNewtonPoly(DataTable *table, int nodes, NewtonPoly *poly)
{
    int n = abs(nodes);
    int points = table->points;
    if (n > points || n == 0)
        return -1;

    double *result = (double *)malloc(n * sizeof(double));
    if (result == NULL)
        return -2;
    double *tmp = (double *)malloc(n * sizeof(double));
    if (tmp == NULL)
        return -2;
    double *xList = (double *)malloc(n * sizeof(double));
    if (xList == NULL)
        return -2;

    if (nodes > 0)
        for (int i = 0; i < n; ++i) {
            xList[i] = table->xList[i];
            tmp[i] = table->yList[i];
        }
    else
        for (int i = points - 1, j = 0; i >= points - n; --i, ++j) {
            xList[j] = table->xList[i];
            tmp[j] = table->yList[i];
        }

    result[0] = tmp[0];
    for (int i = 1; i < n; ++i) {
        for (int j = 0; j <= n - i - 1; ++j)
            tmp[j] = (tmp[j + 1] - tmp[j]) / (xList[j + i] - xList[j]);
        result[i] = tmp[0];
    }

    poly->degree = n - 1;
    poly->coeff = result;
    poly->xList = xList;
    free(tmp);
    return 0;
}

Are there some better approaches to this problem?

Answer 1

Don't cast the result from malloc()

Instead, include <stdlib.h> so that the compiler knows it returns a void*. Also, it's a good practice to use the object name as argument to sizeof, so its type can be changed in one place (e.g. if you wanted to use float or long double instead). It's idiomatic to test a pointer isn't null by using the ! operator on it.

double *result = malloc(sizeof result * n);
if (!result)
    return -2;

Be careful about memory leaks

I don't see a corresponding free() of result or xList - we should be clear that this is the caller's responsibility. It's also easily missed that the function won't free an existing poly->coeff or poly->xList before overwriting - I think that's a big magnet for inviting leaks.