Subversion Repositories programming

/* Copyright: Ira W. Snyder
 * Start Date: 2005-11-13
 * End Date: 2005-11-13
 * License: Public Domain
 *
 * Changelog Follows:
 *
 * 2005-11-13
 * - Implemented Functions 1-3 from the Project #3 Handout.
 * - Implemented Algorithm 7 (LaGrange) from Notes Set #3.
 * - Implemented Algorithm 8-9 (Newton Divided Diff) from Notes Set #3.
 *
 */

#include <cstdio>
#include <cmath>
using namespace std;

#define X_MIN = -5.0
#define X_MAX = 5.0

double func1 (double x)
{
    return 1.0 / (1.0 + pow(x, 2));
}

double func2 (double x)
{
    // NOTE: M_E is the value of e defined by the math.h header
    return (pow(M_E, x) + pow(M_E, -x)) / 2.0;
}

double func3 (double x)
{
    return pow(x, 14) + pow(x, 10) + pow(x, 4) + x + 1;
}

double LaGrangeMethod (double *x, double *y, int n, double point)
{
    double value = 0;
    double term  = 0;
    int i, j;

    for (i=0; i<n; i++)
    {
        term = y[i];

        for (j=0; j<n; j++)
        {
            if (i != j)
                term *= (point - x[j])/(x[i] - x[j]);
        }

        value += term;
    }

    return value;
}

double NewtonMethod (double *x, double *y, int n, double point)
{
    int i, j;
    double a[n];

    for (i=0; i<n; i++)
        a[i] = y[i];

    for (i=1; i<n; i++)
        for (j=n-1; j>=i; j--)
            a[j] = (a[j] - a[j-1]) / (x[j] - x[j-i]);

    // At this point, all of the a[i]'s have been calculated,
    // using a Divided Difference Table.

    double xterm = 1.0;
    double value = 0.0;

    for (i=0; i<n; i++)
    {
        value += a[i] * xterm;
        xterm *= (point - x[i]);
    }

    return value;
}

int main (void)
{
    double x[] = {-2, -1, 0, 2, 3};
    double y[] = {15, -3, -5, 15, 85};

    int size = 5;
    double point = 2.0;

    printf ("LaGrange = %e\n", LaGrangeMethod(x, y, size, point));
    printf ("Newton = %e\n", NewtonMethod(x, y, size, point));

    return 0;
}