// This class does the Newton Interpolation
import java.lang.Math;

class Interpolator {
    // Constructor for no points
    public Interpolator() {
        pts     = 0;
        alloced = 1;
        x       = new double[1];
        ddt     = new double[1][];
    }

    // Copy constructor of sorts, with a twist. This twist is used for
    // the copy-and-add constructor.
    public Interpolator(Interpolator frm, int extras) {
        int i,j;

        pts = frm.pts;

        if (extras < 0)
            alloced = frm.pts;
        else
            alloced = frm.pts + extras;

        x   = new double[alloced];
        ddt = new double[alloced][];

        for (i=0;i<pts;i++) x[i] = frm.x[i];
        for (i=0;i<pts;i++) ddt[i] = new double[i+1];
        for (i=0;i<pts;i++)
            for (j=0;j<=i;j++)
                ddt[i][j] = frm.ddt[i][j];
    }

    // Copy constructor, derived from above.
    public Interpolator(Interpolator frm) {
        this(frm, 0);
    }

    // Copy and add constructor. Since it uses the copy-and-allocate-extra
    // constructor, it is faster than adding after a copy.
    public Interpolator(Interpolator frm, double xn, double yn) throws DoubleInterpolationException {
        this(frm, 1);
        Add(xn,yn);
    }

    // Add a point to what needs to be interpolated.
    // Not threadsafe because it plays with pts.
    public synchronized void Add(double xn, double yn) throws DoubleInterpolationException {
        int i;

        // Not enough room: allocate more
        if (pts == alloced)
            ReAlloc();

        // Add the point
        x[pts]      = xn;
        ddt[pts]    = new double[pts+1];
        ddt[pts][0] = yn;

        // Construct the DDT
        for (i=1;i<=pts;i++)
            if (xn == x[pts-i])
                throw new DoubleInterpolationException();
            else
                ddt[pts][i] = (ddt[pts][i-1]-ddt[pts-1][i-1])/(xn-x[pts-i]);

        // Class methods would assume that a[pts-1] is valid, so
        // we only increment this after the invariant is satisfied.
        pts++;
    }

    // Nested multiplication evaluation
    public double Eval(double xval) {
        if (pts == 0) return 0;

        double d = getA(pts-1);
        int    i;

        for (i=pts-2;i>=0;i--)
            d = d*(xval-x[i]) + getA(i);

        return d;
    }


    // Handy-dandy function for doing the homework
    public double[] MultiEval(double xval[]) {
        double y[] = new double[xval.length];
        int    i;

        for (i=0;i<xval.length;i++) {
            y[i] = Eval(xval[i]);
            System.out.println(double_to_string(y[i]));
        }

        return y;
    }

    // Another homework-helper. Gives m points in [a,b) and the m+1st
    // point given is b. Wants m > 0.
    public double[] MiddlePoints(double a, double b, int m) {
        if (m <= 0) return null;

        double x[] = new double[m+1];
        double d = (b-a)/((double)m);
        int    i;

        x[0] = a;
        x[m] = b;

        for (i=1;i<m;i++)
            x[i] = x[i-1] + d;

        return x;
    }

    // Reset the Interpolator. Not threadsafe, but it runs REALLY fast
    // so we don't care too much.
    public synchronized void Reset() {
        pts = 0;
    }

    // I presume the user of getX|A will be careful about checking
    // i against size(). Design decision to make things run faster.
    public double getX(int i) {return x[i];}
    public double getY(int i) {return ddt[i][0];}
    public double getA(int i) {return ddt[i][i];}
    public int    size()   {return pts; }

    // This is a local function... Java doesn't like scientific format.
    // I wrote this little helper to do Scientific I/O
    public String double_to_string(double d) {
        // The main portion of the code croaks if d is zero.
        // Detect that special case now.
        if (d == 0)
            return "0";

        // isFlipped: The code also doesn't like negative numbers
        //            Set this flag if it is negative
        // absolute:  Well, we need to print the absolute value of the
        //            number after we print the sign.
        // logOfTen:  This is useful for some stuff we do later...
        // log_10:    The base ten logarithm of the number: for the e field
        // actpow:    The "actual power", that goes to the e field
        // absActPow: The absolute value of the actual power.
        // mantissa:  The mantissa. Divide absolute by 10^actpow
        // myStuff:   For putting partial solutions into
        // head:      Character to write

        boolean isFlipped = (d < 0);
        double  absolute  = Math.abs(d);
        double  logOfTen  = Math.log(10.0);
        double  log_10    = Math.log(absolute)/logOfTen;
        int     actpow    = (int)Math.floor(log_10);
        int     absActPow = Math.abs(actpow);
        double  mantissa  = absolute/Math.pow(10,actpow);
        char    myStuff[] = new char[50];
        int     head      = 0;
        int     thisDigit;

        // OK we were negative, let's output the minus sign
        if (isFlipped)
            myStuff[head++] = '-';

        // Find and output the leading digit
        thisDigit       = (int) Math.floor(mantissa);
        myStuff[head++] = (char)(((char)thisDigit) + '0');
        myStuff[head++] = '.';

        // For each digit we're precise to, find the nth digit and print it
        for (int i=0;i<15;i++) {
            mantissa        -= thisDigit;
            mantissa        *= 10;
            thisDigit        = (int) Math.floor(mantissa);
            myStuff[head++]  = (char)(((char)thisDigit) + '0');
        }

        // Time to round up...
        if (mantissa - thisDigit > 0.5) {
            while (myStuff[head-1] == '9')
                head--;

            if (myStuff[head-1] == '.') {
                head--;
                if (myStuff[head-1] == '9') {
                    myStuff[head-1] = '1';
                    actpow++;
                    absActPow = Math.abs(actpow);
                } else {
                    myStuff[head-1] = (char)(myStuff[head-1] + 1);
                }
            } else {
                myStuff[head-1] = (char)(myStuff[head-1] + 1);
            }
        } else {
            // Purge trailing zeros
            while (myStuff[head-1] == '0')
                head--;

            // If all that's left is the decimal point, purge that too
            if (myStuff[head-1] == '.')
                head--;
        }

        if (actpow != 0) {
            // Do an exponent field
            myStuff[head++] = 'e';

            if (actpow < 0)
                myStuff[head++] = '-';

            if (absActPow >= 1000)
                myStuff[head++] = (char) ((absActPow/1000)     + '0');

            if (absActPow >= 100)
                myStuff[head++] = (char) (((absActPow/100)%10) + '0');

            if (absActPow >= 10)
                myStuff[head++] = (char) (((absActPow/10)%10)  + '0');

            myStuff[head++] = (char) ((absActPow%100) + '0');
        }

        // Finally, done!
        return new String(myStuff,0,head);
    }

    // Reallocate storage: need more. Synchronized for paranoid reasons.
    protected synchronized void ReAlloc() {
        double tx[]     = new double[alloced];
        double tddt[][] = new double[alloced][];

        // Backup x, reallocate x, and restore the backup
        System.arraycopy(x, 0, tx, 0, alloced);
        x = new double[alloced*2];
        System.arraycopy(tx, 0, x, 0, alloced);

        // Backup ddt, reallocate ddt, and restore the backup
        System.arraycopy(ddt, 0, tddt, 0, alloced);
        ddt = new double[alloced*2][];
        System.arraycopy(tddt, 0, ddt, 0, alloced);

        alloced *= 2;
    }

    protected double x[];
    protected double ddt[][];
    protected int pts;     // How many points are we interpolating
    protected int alloced; // How much room have we allocated
}