umontreal.iro.lecuyer.functionfit.PolInterp Maven / Gradle / Ivy

Go to download

Show more of this group Show more artifacts with this name
Show all versions of ssj Show documentation

SSJ is a Java library for stochastic simulation, developed under the direction of Pierre L'Ecuyer, in the Département d'Informatique et de Recherche Opérationnelle (DIRO), at the Université de Montréal. It provides facilities for generating uniform and nonuniform random variates, computing different measures related to probability distributions, performing goodness-of-fit tests, applying quasi-Monte Carlo methods, collecting (elementary) statistics, and programming discrete-event simulations with both events and processes.

The newest version!



/*
 * Class:        PolInterp
 * Description:  polynomial that interpolates through a set of points
 * Environment:  Java
 * Software:     SSJ 
 * Copyright (C) 2001  Pierre L'Ecuyer and Université de Montréal
 * Organization: DIRO, Université de Montréal
 * @author       
 * @since

 * SSJ is free software: you can redistribute it and/or modify it under
 * the terms of the GNU General Public License (GPL) as published by the
 * Free Software Foundation, either version 3 of the License, or
 * any later version.

 * SSJ is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.

 * A copy of the GNU General Public License is available at
   GPL licence site.
 */

package umontreal.iro.lecuyer.functionfit;

import java.io.Serializable;
import umontreal.iro.lecuyer.functions.Polynomial;

import umontreal.iro.lecuyer.util.PrintfFormat;
import cern.colt.matrix.DoubleMatrix2D;
import cern.colt.matrix.impl.DenseDoubleMatrix2D;
import cern.colt.matrix.linalg.Algebra;


/**
 * Represents a polynomial that interpolates through a set of points. More
 * specifically, let 
 * (x₀, y₀),…,(x_n, y_n) be a set of points and
 * p(x) the constructed polynomial of degree n. Then, for 
 * i = 0,…, n,
 * 
 * p(x_i) = y_i.
 * 
 */
public class PolInterp extends Polynomial implements Serializable {
   private static final long serialVersionUID = -710451931485296501L;
   private static final Algebra alg = new Algebra ();
   private double[] x;
   private double[] y;



   /**
    * Constructs a new polynomial interpolating through the given points
    *  (x[0], y[0]), ..., (x[n], y[n]). This constructs a polynomial of
    *  degree n from n+1 points.
    * 
    * @param x the x coordinates of the points.
    * 
    *    @param y the y coordinates of the points.
    * 
    *    @exception NullPointerException if x or y are null.
    * 
    *    @exception IllegalArgumentException if the lengths of x and y are different,
    *                or if less than two points are specified.
    * 
    */
   public PolInterp (double[] x, double[] y) {
      super (getCoefficients (x, y));
      this.x = x.clone ();
      this.y = y.clone ();
   }


   /**
    * Computes and returns the coefficients the polynomial interpolating
    *  through the given points (x[0], y[0]), ..., (x[n], y[n]). 
    *  This polynomial has degree n and there are n+1 coefficients.
    * 
    * @param x the x coordinates of the points.
    * 
    *    @param y the y coordinates of the points.
    * 
    *    @return the coefficients the interpolating polynomial.
    * 
    */
   public static double[] getCoefficients (double[] x, double[] y) {
      if (x.length != y.length)
         throw new IllegalArgumentException (
               "x and y must have the same length");
      if (x.length <= 1)
         throw new IllegalArgumentException ("At least two points are needed");
      final DoubleMatrix2D u = new DenseDoubleMatrix2D (x.length, x.length);
      for (int i = 0; i < x.length; i++) {
         double v = 1;
         for (int j = 0; j < x.length; j++) {
            u.setQuick (i, j, v);
            v *= x[i];
         }
      }
      final DoubleMatrix2D yMat = new DenseDoubleMatrix2D (x.length, 1);
      yMat.viewColumn (0).assign (y);
      final DoubleMatrix2D bMat = alg.solve (u, yMat);
      return bMat.viewColumn (0).toArray ();
   }


   /**
    * Returns the x coordinates of the interpolated points.
    * 
    * @return the x coordinates of the interpolated points.
    * 
    */
   public double[] getX() {
      return x.clone ();
   }


   /**
    * Returns the y coordinates of the interpolated points.
    * 
    * @return the y coordinates of the interpolated points.
    * 
    */
   public double[] getY() {
      return y.clone ();
   }


   /**
    * Makes a string representation of a set of points.
    * 
    * @param x the x coordinates of the points.
    * 
    *    @param y the y coordinates of the points.
    * 
    *    @return the string representing the points.
    * 
    */
   public static String toString (double[] x, double[] y) {
      final StringBuilder sb = new StringBuilder ("Points: ");
      for (int i = 0; i < x.length; i++) {
         if (i > 0)
            sb.append (", ");
         final String xstr = PrintfFormat.format (8, 3, 3, x[i]);
         final String ystr = PrintfFormat.format (8, 3, 3, y[i]);
         sb.append ('(').append (xstr).append (", ").append (ystr).append (')');
      }
      return sb.toString ();
   }

   @Override



   /**
    * Calls {@link #toString( ) toString}(double[], double[]) with the
    *  associated  points.
    * 
    * @return a string containing the points.
    * 
    */
   public String toString() {
      return toString (x, y);
   }


   public PolInterp clone() {
      final PolInterp p = (PolInterp) super.clone ();
      p.x = x.clone ();
      p.y = y.clone ();
      return p;
   }
}