// *=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=* // ** Copyright UCAR (c) 1990 - 2016 // ** University Corporation for Atmospheric Research (UCAR) // ** National Center for Atmospheric Research (NCAR) // ** Boulder, Colorado, USA // ** BSD licence applies - redistribution and use in source and binary // ** forms, with or without modification, are permitted provided that // ** the following conditions are met: // ** 1) If the software is modified to produce derivative works, // ** such modified software should be clearly marked, so as not // ** to confuse it with the version available from UCAR. // ** 2) Redistributions of source code must retain the above copyright // ** notice, this list of conditions and the following disclaimer. // ** 3) Redistributions in binary form must reproduce the above copyright // ** notice, this list of conditions and the following disclaimer in the // ** documentation and/or other materials provided with the distribution. // ** 4) Neither the name of UCAR nor the names of its contributors, // ** if any, may be used to endorse or promote products derived from // ** this software without specific prior written permission. // ** DISCLAIMER: THIS SOFTWARE IS PROVIDED "AS IS" AND WITHOUT ANY EXPRESS // ** OR IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED // ** WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. // *=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=*=* /*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*/ // RCS info // $Author: dixon $ // $Locker: $ // $Date: 2016/03/03 18:46:08 $ // $Id: NewtonPtFunc.cc,v 1.8 2016/03/03 18:46:08 dixon Exp $ // $Revision: 1.8 $ // $State: Exp $ /**-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-**/ /********************************************************************* * NewtonPtFunc.cc: class implementing a function defined by a series * of points. The function value for points is * calculated by using Newton polynomials. * * RAP, NCAR, Boulder CO * * January 2000 * * Nancy Rehak * *********************************************************************/ #include #include #include #include #include #include using namespace std; /********************************************************************** * Constructors */ NewtonPtFunc::NewtonPtFunc(map< double, double, less > function, const double tolerance) : PtFunction(function), _tolerance(tolerance), _degree(0) { // Do nothing } /********************************************************************** * Destructor */ NewtonPtFunc::~NewtonPtFunc(void) { } /********************************************************************** * computeFunctionValue() - Compute the function value using Newton * polynomial interpolation. */ double NewtonPtFunc::computeFunctionValue(const double x) { double *table; double delta; double m = 1.0; double p; // Allocate space for the table table = new double[_function.size()]; // Set the initial interpolating polynomial and first element // of the divided-difference table to f(x0). p = table[0] = (*_function.begin()).second; _degree = 0; map< double, double, less >::iterator k_iter; map< double, double, less >::iterator k_1_iter; map< double, double, less >::iterator i_iter; unsigned int k; int i; for (k = 0, k_iter = _function.begin(); k < _function.size() - 1; ++k, ++k_iter) { k_1_iter = k_iter; ++k_1_iter; ++_degree; table[k + 1] = (*k_1_iter).second; // Compute the next diagonal placed in the // divided-difference table to obtain the next coefficient. for (i = k, i_iter = k_iter; i >= 0; --i, --i_iter) table[i] = (table[i + 1] - table[i])/((*k_1_iter).first - (*i_iter).first); // Compute the next pk(x), one degree larger than the // last, and determine whether the specified tolerance will // be met. m = m * (x - (*k_iter).first); delta = table[0] * m; p = p + delta; if (fabs(delta) < _tolerance) { delete [] table; return p; } } /* endfor - k */ // The specified tolerance was never met for differences // between successive interpolating polynomials. delete [] table; return DBL_MAX; } /********************************************************************** * print() - Print the function information to the given stream. */ void NewtonPtFunc::print(FILE *stream) { fprintf(stream, "Newton Polynomial Interpolation Function:\n"); fprintf(stream, "\n"); fprintf(stream, " tolerance = %f\n", _tolerance); fprintf(stream, "\n"); _printFunction(stream); } void NewtonPtFunc::print(ostream &stream) { stream << "Newton Polynomial Interpolation Function:" << endl; stream << endl; stream << " tolerance = " << _tolerance << endl; stream << endl; _printFunction(stream); } /********************************************************************** * Private Member Functions * **********************************************************************/