JPolynome.hh

Go to the documentation of this file.

#ifndef __JMATH__JPOLYNOME__
#define __JMATH__JPOLYNOME__

#include <istream>
#include <ostream>
#include <vector>
#include <cmath>
#include <limits>
#include <stdarg.h>


/**
 * \author mdejong
 */

namespace JMATH {}
namespace JPP { using namespace JMATH; }

namespace JMATH {

  /**
   * Polynome function object.
   *
   * Evaluates function, derivative and integral values.
   */
  class JPolynome : 
    public std::vector<double>
  {
  public:
    /**
     * Default constructor.
     */
    JPolynome()
    {}


    /**
     * Constructor.
     *
     * \param  n          degree of polynome (followed by comma separated argument list)
     */
    JPolynome(int n, ...)
    {
      va_list  ap;

      va_start(ap, n);
    
      for (double x; n != 0; --n) {

        x = va_arg(ap, double);  // next argument

        push_back(x);
      }

      va_end(ap);
    }


    /**
     * Constructor.
     *
     * \param  __begin    begin of data
     * \param  __end      end   of data
     */
    template<class T>
    JPolynome(T __begin,
              T __end)
    {
      for (T i = __begin; i != __end; ++i)
        push_back(*i);
    }


    /**
     * Function value.
     *
     * \param  x          abscissa value
     * \return            function value
     */
    double getValue(const double x) const
    {
      double y = 0.0;
      double z = 1.0;

      for (const_iterator i = begin(); i != end(); ++i, z *= x) {
        y += (*i) * z;
      }

      return y;
    }


    /**
     * Derivative value.
     *
     * \param  x          abscissa value
     * \return            derivative value
     */
    double getDerivative(const double x) const
    {
      double y = 0.0;

      if (!empty()) {

        double z = 1.0;
        int    n = 1;

        for (const_iterator i = begin(); ++i != end(); ++n, z *= x) {
          y += (*i) * z * n;
        }
      }

      return y;
    }


    /**
     * Integral value.
     *
     * \param  x          abscissa value
     * \return            integral value
     */
    double getIntegral(const double x) const
    {
      double y = 0.0;
      double z = x;
      int    n = 1;

      for (const_iterator i = begin(); i != end(); ++i, ++n, z *= x) {
        y += (*i) * z / n;
      }

      return y;
    }


    /**
     * Function value.
     *
     * \param  x          abscissa value
     * \return            function value
     */
    double operator()(const double x) const
    {
      return getValue(x);
    }


    /**
     * Derivative function.
     *
     * \return            derivative function
     */
    JPolynome getDerivative() const
    {
      JPolynome buffer;

      int n = 1;

      for (const_iterator i = begin(); ++i != end(); ++n) {
        buffer.push_back(n * (*i));
      }

      return buffer;
    }


    /**
     * Integral function.
     *
     * \return            integral function
     */
    JPolynome getIntegral() const
    {
      JPolynome buffer(1, 0.0);

      int n = 1;

      for (const_iterator i = begin(); i != end(); ++i, ++n) {
        buffer.push_back((*i) / (double) n);
      }

      return buffer;
    }

    
    /**
     * Equality.
     *
     * \param  P        polynome
     * \param  eps      numerical precision
     * \return          true if polynomes identical; else false
     */
    bool equals(const JPolynome& P,
                const double     eps = std::numeric_limits<double>::min()) const
    {
      if (this->size() == P.size()) {

        for (const_iterator p = this->begin(), q = P.begin(); p != this->end(); ++p, ++q) {
          if (fabs(*p - *q) > eps) {
            return false;
          }
        }

        return true;
        
      } else {

        return false;
      }
    }

    
    /**
     * Read polynome from input.
     *
     * \param  in         input stream
     * \param  object     polynome
     * \return            input stream
     */
    friend inline std::istream& operator>>(std::istream& in, JPolynome& object)
    {
      for (double x; in >> x; ) {
        object.push_back(x);
      }

      return in;
    }


    /**
     * Write polynome to output.
     *
     * \param  out        output stream
     * \param  object     polynome
     * \return            output stream
     */
    friend inline std::ostream& operator<<(std::ostream& out, const JPolynome& object)
    {
      for (JPolynome::const_iterator i = object.begin(); i != object.end(); ++i) {
        out << ' ' << *i;
      }

      return out;
    }
  };
}

#endif