您的位置：首页 > 编程语言 > Java开发

Java in Science: Data Interpolation and Extrapolation Using Numerical Methods of Polynomial Fittings

2004-08-20 21:26 781 查看

Java in Science: Data Interpolation and Extrapolation Using Numerical Methods of Polynomial Fittings, Part 1
By Sione Palu

The value of a function g(x) at a specific set of points x1,x2,x3,...,xn is sometimes known, but there is no analytical expression for g(x) that enables a calculation of its value at an arbitrary point. An estimation of g(x) for an arbitrary x in a way that drawing a smooth curve through (xi,g(xi)) can be achieved by numerical analysis.

As an example of looking at a sales profit (in millions of dollars) P = [3.4, 5.1, 7.6, 9] for every two month period T = [2, 4, 6, 8] from a manufacturing company, and when one is trying to estimate P where T = 5(inside T-boundary) or T = 10(outside T-boundary) it would be inaccurate without data smoothing techniques of numerical analysis. If the desired x is bounded between [x1,...,xn] that is x is in between the largest and the smallest of the xi's, this is called interpolation and if x is outside the boundaries then it is called extrapolation.

Polynomials

Interpolation and extrapolation techniques must model the function, within and beyond the boundaries, [x1,...,xn] by some functional form. The functional form that is going to be discussed in this tutorial is known as polynomial fittings. Polynomials of order n have the general form of:

y = anxn + a(n-1)x(n-1) + ... + a2x2 + a1x + a0

n is an integern >= 0.an is a Real number.an is the coefficient of the nth term.Example of a 4th-order polynomial: y = -2x4 + 5x3 + 3x2 - 3x + 7 , thus polynomial coefficients are: a4 = -2, a3 = 5, a2 = 3, a1 = -3, a0 = 7

In this tutorial, we will define and derive mathematical operations mainly in matrix algebra for how to fit a polynomial of order n into a set of known or observed data pairs , [xi,yi]. A polynomial of any order is specified and then the coefficients are returned as an array of doubles (Java data-type).

Linear (straight-line) fitting is a polynomial of order 1 (also known as a first-order polynomial, and it has wide applications in modeling from simulations, graphics to financial forecast . Before we delve in to the Java coding , I want to introduce the concept of linear algebra and matrix operations so the general reader will understand the code with minimum difficulties.

Matrix Algebra

CAD (computer-aided design) and CAM (computer-aided manufacturing) software have revolutionized the automotive industry, for instance. Today computer graphics forms the heart, and linear algebra the soul of modern car design. Before a new car is built, engineers design and construct a mathematical car -- a wire-frame model that exists only in computer memory and on graphics display terminals. The car is also mathematically tested (to simulate impact crash-force, engine heat conduction, metal vibration-fatigue frequency, etc.) before even a prototype is built. Linear algebra is indeed the foundation of any graphics software because all the manipulations of screen images are accomplished via algebraic techniques.

If A is an m x n matrix -- that is, a matrix with m rows and n columns -- then the scalar entry in the ith row and the jth column of A is denoted by aij and is called the (i,j) -- entry of A. The class Matrix.java has an internal array storage as a private instance variable named A and it has an accessor method to read its value.

Jama and Jamlab

Jama is a basic linear algebra package for Java. The classes in this package enable the construction and manipulation of real, dense matrices. Jama is sufficient enough to provide functionality for routine problems and packaged in a way that is natural and understandable to non-experts. It is intended to serve as the standard matrix class for Java and will be proposed as such to the Java Grande Forum, which represents those from academia and corporations drafting proposals to Sun on Java APIs for scientific and engineering applications. Jama was developed by the National Institute of Standards and Technology (NIST) and The MathWorks Inc.. It is available for free download.

The Jama package has six classes:

Matrix
CholeskyDecomposition
LUDecomposition
QRDecomposition
SingularValueDecomposition
EigenvalueDecomposition

The Matrix class provides the fundamental operations of numerical linear algebra. Various constructors create matrices from two-dimensional arrays of double-precision floating point numbers. Various gets and sets provide access to submatrices and matrix elements. The basic arithmetic operations include matrix addition and multiplication, matrix norms and selected element-by-element array operations. In this tutorial , we discuss matrix operations using Jama.

Jamlab (Java Matrix Laboratory) is my own developed Java matrix package, which is based on the functionality of the popular scientific and engineering language called MatLab from MathWorks. Jamlab is my own attempt to write or translate functions written in MatLab to a Java version. This package contain classes that are subclassed from class Matrix of the package Jama. There are currently six classes in Jamlab:

Jdatafun
Jelfun
Jelmat
RowColumnIndex
Polyfit
Polyval

Jamlab is not yet complete -- a full version will be available for free download in the future. The reader is encouraged to download Jama and Jamlab prior to continuing on to the next part of the series.

Download Jama and Jamlab and [i]Jamlab documentation.

References

Java for Engineers and Scientists by Stephen J. Chapman, Prentice Hall, 1999.
Introductory Java for Scientists and Engineers by Richard J. Davies, Addison-Wesley Pub. Co., 1999.
Applied Numerical Analysis (Sixth Edition) by Curtis F. Gerald and Patrick O. Wheatly, Addison-Wesley Pub. Co., 1999.
Linear Algebra and Its Applications (Second Edition), by David C. Lay, Addison-Wesley Pub. Co.
Numerical Recipes in Fortran 77, the Art of Scientific Computing (Volume 1) by William H. Press, Saul A. Teukolsky, William T. Vetterling, and Brian P. Flannery, Cambridge University Press, 1997.
Mastering MATLAB 5: A Comprehensive Tutorial and Reference by Duane Hanselman and Bruce Littlefield, Prentice Hall, 1997.
Advanced Mathematics and Mechanics Applications using MATLAB by Louis H. Turcotte and Howard B. Wilson, CRC Press, 1998.

Java as a Scientific Programming Language (Part 1): More Issues for Scientific Programming in Java
Scientific Computing in Java (Part 2): Writing Scientific Programs in Java
Using Java in Scientific Research: An Introduction

About the Author

double[][] arrayA = {{4.,0.},{5.,3.},{-2.,1.}};
//Instantiate matrix object A
Matrix A = new Matrix(arrayA);
double[][] arrayB = {{4.,0.},{5.,4.},{-2.,1.}};
//Instantiate matrix object B
Matrix B = new Matrix(arrayB);
double[][] arrayC = {{1.,7.},{7.,4.}};
//Instantiate matrix object C
Matrix C = Matrix(arrayC);
double[][] arrayD = {{1.,0.},{0.,1.}};
//Instantiate matrix object D
Matrix D = new Matrix(arrayD);

The additions (sums) and subtractions (difference) of matrices is done element-wise. The sum of A+B would result in a matrix the same size as A or B with all entries as the sum of aij and bij , for example, the resultant entry of (1,1) = 7, which is a11+b11 = 3+4 = 7. The sum of two different size (dimemsion) matrices is undefined such as A+C, because size(A) = 3 by 2 and size(C) = 2 by 2. In multiplying a matrix by a scaler results in a matrix with all the entries scaled by this scaler-value. 3*B is a multiplication of B by a scaler value of 3 in which all entries are now 3bij.

Matrix A_plus_B = A.plus(B) ;
//returned matrix has all its entries scaled by 3.
Matrix matrix_3_times_B = B.times(3.);
//This following line will throw an IllegalArgumentException because matrix dimensions do 
//not agree, matrix A is a 4 by 2 but C is a 2 by 2 size matrix. 
Matrix A_plus_C = A.plus(C);

Matrix Multiplication

Under matrix multiplication, such as A*C, the inner dimensions have to be the same for it to be defined, regardless of the outer dimensions. In this case size(A) = [3 x 2] and size(C) = [2 x 2]; therefore, matrix multiplication is defined because they have the same inner dimensions of two, [3 x 2][2 x 2]. The dimension of the result of multiplying two matrices would be the outer dimensions of the respective matrices. A*C would result in a matrix with size(A*C) = [3 x 2], that is [3 x 2][2 x 2]. Matrix multiplication results in all entries(i,j) equals the summation of element-wise multiplication of each row of A to each column of C. Multiplication of A *B is undefined, because the inner dimensions do not agree: A* B is [3 x 2]*[3 x 2]. Multiplication of matrices is not always commutative (order of operation), B *C is defined, but reversing the order as C*B is undefined but C* D = D*C. Matrix multiplication is directly in contrast with its elementary counterpart, 3 x 4 = 4 x 3 = 12. Matrix D is a special type called Identity, and this is the equivalent in elementary arithmetic of the number 1, which is the multiplication identity, 6 x 1 = 6. Identity is a square-matrix (number of rows = columns) with any dimension and all its elements at the diagonal, from top-left to bottom-right are ones and zeros everywhere. Any matrix that is multiplied to the Identity is always that matrix, examples are: A*D = A, B* D = B and C*D = C.

Matrix A_times_C = A.times(C) ;	
//The following line will throw an IllegalArgumentException because inner matrix 
//dimensions do not agree.
Matrix A_times_B = A.times(B);

Transpose of a Matrix

Given an m x n matrix E, the transpose of E is the n x m matrix denoted by ET, whose columns are formed from the corresponding rows of E. General matrix transposition operations:

1. (F + G)T = FT + GT

Matrix result1 = F.plus(G).transpose() ;						
Matrix result2 = F.transpose().plus(G.transpose());
//Matrix result1 and result2 are equals (same value in all entries)

2. (FT)T = F

Matrix result3 = F.transpose().transpose() ;										
//Matrix result3 and F are equals (same value in all entries)

3. (F * G)T = GT * FT

Matrix result4 = F.times(G).transpose() ;
Matrix result5 = G.transpose().times(F.transpose()) ;							
//Matrix result4 and result5 are equals (same value in all entries)

4. For any scalar r, (r * F)T = r * FT

Matrix result6 = F.times(r).transpose() ;
Matrix result7 = F.transpose().times(r) ;									
//Matrix result6 and result7 are equals (same value in all entries)

Inverse of a Matrix

If K is an [n x n] matrix, sometimes there is another [n x n] matrix L such that:
K * L = I and L * K = I
where I is an [n x n] Identity matrix

We say K is invertible and L is an inverse of K. The inverse of K is denoted by K-1.

K * K-1 = I, where I - is the Identity matrix
K-1 * K = I

Matrix inverseMatrix = K.inverse() ;												
//Matrix inverseMatrix is an inverse of matrix K

QR Matrix Factorization

Factorizing a matrix follows a similar concept in elementary arithmetics. The number 30 has prime factors of = 2 x 3 x 5. If P is an [m x n] matrix with linearly independent columns, then it can be factored as:

P = Q * R:
where Q is an [m x m] whose columns form an orthonormal basis and R is an [m x n] upper triangular invertible matrix with positive entries on its diagonal.

double[][] arrayP = {{1.,0.},{1.,1.},{1.,1.},{1.,1.}};
Matrix P = new Matrix(arrayP);

//Create a QR factorisation object from matrix P.
QRDecomposition QR = new QRDecomposition(P);
Matrix Q = QR.getQ();
Matrix R = QR.getR(); 
//Now that we have two matrix factors for matrix P which:  Q*R = P

Elementary Matrix Function and Special Matrices

repmat (Matrix mat, int row, int col) - static method from JElmat class of Jamlab package
Repeat a matrix in a tiling manner. If matrix A is to be tiled [3 by 2], then A would be tiled in a 3-rows by 2-columns pattern.

double[][] arrayA = {{7,5},{2,3}}; 
Matrix A = new Matrix(arrayA);
//Create matrix object B by tiling matrix A in a 3 by 2 pattern.
Matrix B = JElmat.repmat(A,3,2);

vander (double[] a,int numCol) - static method from JElmat.
Create a Vandermonde matrix in which the second to last column are the elements of array a[]. The matrix has numCol number of columns.

//Create Vandermonde-matrix
double[] arrayA = {5,-1,7,6};
//vander is also overloaded with one input argument a double[] array, where columns = rows
Matrix A = JElmat.vander(arrayA);
double[] arrayB = {5,-1,7,6,-3,4};
//Create a Vandermonde-matrix with only four columns.
Matrix B = JElmat.vander(arrayB,4);

sum (Matrix matrix) - static method from JDatafun class.
sum all the elements of the column of a matrix and return a single row matrix.

double arrayA = {{19,18,16,18},{5,15,9,15},{12,9,12,4},{10,0,16,8}};
Matrix A = new Matrix(arrayA);
//Sum all the columns of matrix A.
Matrix B = JDatafun.sum(A);

linespace (double leftBound, double rightBound, int nPoints) - static method from JElmat class
Linear spaced matrix row vector. Generates nPoints between leftBound and rightBound.

//row matrix A, starting at -8.	and ending at 12. with eleven points equally spaced.
//internal array A = (-8    -6    -4    -2     0     2     4     6     8    10    12)
Matrix A = JElmat.linspace(-8.,12.,11);

reshape (Matrix A, int M, int N) - static method from JElmat class.
returns the M-by-N matrix whose elements are taken columnwise from A.

double[][] arrayA = {{1.,2.,3.,4.},{5.,6.,7.,8.},{9.,10.,11.,12.}}; 
Matrix A = new Matrix(arrayA); 
//reshape matrix A from 3 x 4 dimension into 6 x 2
Matrix B = JElmat.reshape(A,6,2);

zeros (int m, int n) - static method from JElmat class.
returns an m x n matrix of zeros.

ones (int m, int n) - static method from JElmat class.
returns an m x n matrix of ones.

//Create a matrix of ones with 3 x 4 size.
Matrix X = JElmat.ones(3,4);
//Create a matrix of zeros with 2 x 2 size.
Matrix Y = JElmat.zeros(2,2);

We will continue with polynomial fittings and Java code in Part 3.

Download Jama and Jamlab and Jamlab documentation.

References

About the Author

Java Codes

Polyfit Class
This class has one constructor that instantiates a Polyfit object. The constructor takes two input argument arrays with an integer to specify the polynomial order to fit data into. All instance variables are private and all instance methods are accessor. The code is descriptive as much as possible for the non-expert reader.

package jamlab;

import Jama.*;
import Jama.util.*;

public class Polyfit
{
//instance variable for upper triangular matrix
private Matrix R;

//array for polynomial coefficients
private double[] polynomialCoefficients ;

//degree of freedom
private int degreeOfFreedom;

//norm - largest singular value.
private double norm;

//y-intercept : f(x)-> x=0, f(x=0)
private double yIntercept ;

//store polynomial coefficients in matrix
private Matrix polyCoeffMatrix ;

//Contructor with arrays of data pairs , _x & _y with approximating order
public Polyfit(double[] _x , double[] _y , int order) throws Exception{
int xLength = _x.length,
yLength = _y.length;
double tempValue = 0.0;

//store array _y as a matrix
Matrix y2D = new Matrix(JElmat.convertTo2D(_y));

//transpose matrix y2D
Matrix yTranspose = y2D.transpose();

//throw Exception if data pairs are not the same size.
if(xLength != yLength){
throw new Exception(" Polyfit :- The lengths of the 2-input array parameters must be equal.");
}

//throw Exception if only one data pair.
if(xLength < 2){
throw new Exception(" Polyfit :- There must be at least 2 data points for polynomial fitting.");
}

//throw Exception if polynomial order is negative.
if(order < 0){
throw new Exception(" Polyfit :- The polynomial fitting order cannot be less than zero.");
}

//throw Exception if polynomial order exceeds or equal the number of data-pairs.
if(order >= xLength){
throw new Exception(" Polyfit :- The polynomial order = "+order+" , must be less than the number of data points = "+xLength);
}

Matrix tempMatrix = null;
try{
//Create Vandermonde matrix
tempMatrix = JElmat.vander(_x,order+1);
}
catch(Exception ex){
ex.printStackTrace();
}

//Factorise matrix tempMatrix by QR decomposition
QRDecomposition qr = new QRDecomposition(tempMatrix);

//Initialize Orthonormal basis matrix Q from object qr
Matrix Q = qr.getQ();

//Create Upper-triangular invertible matrix
R = qr.getR();

//result = R-1 * ( QT * _yT )

Matrix result = R.inverse().times(Q.transpose().times(yTranspose));

//return a reference to the internal array transpose of the result matrix
double[][] arrayResult = result.transpose().getArray();

//polynomial coefficients is the first row of the internl array.
polynomialCoefficients = arrayResult[0];

//degree of freedom
degreeOfFreedom = yLength - (order+1);

Matrix r = yTranspose.minus(tempMatrix.times(result));

//largest singular value of matrix r
norm = r.norm2();

//Create matrix of the coefficients.
polyCoeffMatrix = new Matrix(JElmat.convertTo2D(polynomialCoefficients));

//Polynomial y-axis intercept
yIntercept = polynomialCoefficients[polynomialCoefficients.length-1];

}//end constructor

//return polynomial coefficients array
public double[] getPolynomialCoefficients(){
return polynomialCoefficients ;
}

//return upper-triangular invertible matrix
public Matrix getR(){
return R;
}

//return the degree of freedom
public int getDegreeOfFreedom(){
return degreeOfFreedom;
}

//return norm
public double getNorm(){
return norm;
}

//return y-intercept
public double getYIntercept(){
return yIntercept;
}

//return polynomial coefficients matrix
public Matrix getPolyCoeffMatrix(){
return polyCoeffMatrix ;
}

}//----- End Class Definition -------

Polyval Class
package jamlab;

import Jama.*;
import Jama.util.*;

public class Polyval
{
//Interpolated and extrapolated values corresponding to domain _x .
private double[] Yout;

//Interpolated and extrapolated values in a matrix object .
private Matrix youtMatrix ;

//Error estimates in an array.
private double[] ErrorBounds;

//Error estimates in a matrix
private Matrix ErrorBoundsMatrix ;

//Constructor with an array of values for interpolation and extrapolation based on coefficients form polyfit object
public Polyval( double[] _x , Polyfit polyfit){

//Retrieve polynomial coefficients from polyfit object
double[] polyCoeff = polyfit.getPolynomialCoefficients();

//Length of the coefficients array
int coeffLen = polyCoeff.length;

//polyfit object degree of freedom
int degree = polyfit.getDegreeOfFreedom();

//polyfit object norm
double normr = polyfit.getNorm();

//Retrieve polyfit upper triangular matrix instance
Matrix R = polyfit.getR();

//Convert the domain _x into a matrix
Matrix _xMatrix = new Matrix(JElmat.convertTo2D(_x));

//Initialize interpolated and extrapolated output to a matrix of zeros of size 1 x length(_x)
youtMatrix = JElmat.zeros(1,_xMatrix.getColumnDimension());

Matrix pC;

//Loop to calculate estimate output extrapolated or interpolated values
for(int i=0 ; i

//youtMatrix = (youtMatrix .* _xMatrix) + pC

youtMatrix.arrayTimesEquals(_xMatrix).plusEquals(pC);
}

//Order of polynomial
int order = coeffLen;

//Declare Vandermonde matrix;
Matrix V = null;

//Declare error-bound matrix
Matrix errorDelta = null;

//Vandermonde matrix with number of columns = order
try{ V = JElmat.vander(_x,order); }
catch(Exception ex){}

//Rinv = R-1

Matrix Rinv = R.inverse();

//E = V * R-1

Matrix E = V.times(Rinv);

//Declare a matrix of ones.
Matrix matOnes = null;

//Declare a matrix for evaluation of error.
Matrix e = null;

//Horizontal line (Order = 0)
if(order == 0){

//matrix of ones with size = [E.getRowDimension(),E.getColumnDimension()]
matOnes = JElmat.ones(E.getRowDimension(),E.getColumnDimension());

//e = sqrt(matOnes + E.*E)

e = JElfun.sqrt( matOnes.plus( E.arrayTimes(E)));
}
else{ //Higher order, 1,2,3,...

//SumT = (sum((E.*E)T))T

Matrix SumT = JDatafun.sum(E.arrayTimes(E).transpose()).transpose();

//row size
int rowDim = SumT.getRowDimension();

//column size
int colDim = SumT.getColumnDimension();

//matrix of ones , same size as SumT
matOnes = JElmat.ones(rowDim,colDim);

//e = sqrt(matOnes + SumT)

e = JElfun.sqrt(matOnes.plus(SumT));
}

if(degree==1){//When degree of freedom = 1, error-bounds would be infinite.
double[][] inf = new double[1][1];
inf[0][0] = Double.POSITIVE_INFINITY ;
try {

//error matrix of size =[e.getRowDimension(),e.getColumnDimension()] with all entries is an infinity.
errorDelta = new Matrix( JElmat.repmat(inf,e.getRowDimension(),e.getColumnDimension()) );
}
catch(Exception ex) {
ex.printStackTrace();
}
}
else{ // Higher degree of freedom
errorDelta = e.times(normr/Math.sqrt(degree));
}

double[][] tempDelta = {{0.},{0.}} ; //initialize

try{

//reshape matrix errorDelta into tempDelta array with size = [1 , _x.length]
tempDelta = JElmat.reshape(errorDelta,1,_x.length).getArrayCopy() ;
}
catch(Exception ex){ ex.printStackTrace(); }

//Assign ErrorBounds first row of tempDelta
ErrorBounds = tempDelta[0] ;

//ErrorBounds in a matrix
ErrorBoundsMatrix = new Matrix(JElmat.convertTo2D(ErrorBounds));

//get copy of the output array from output matrix
double[][] yM = youtMatrix.getArrayCopy();

//Assign output to first row of the array
Yout = yM[0];
}//End constructor

//return output matrix
public Matrix getYoutMatrix(){
return youtMatrix ;
}

//return output array
public double[] getYout(){
return Yout;
}

//return output error-bounds array
public double[] getErrorBounds(){
return ErrorBounds;
}

//return output error-bounds matrix
public Matrix getErrorBoundsMatrix(){
return ErrorBoundsMatrix;
}

}// --- End Class Definition ---

We will conclude by running the complete program in Part 4.

Download Jama and Jamlab and Jamlab documentation.

References

About the Author

Sione Palu is a Java developer at Datacom in Auckland, New Zealand, currently involved in a Web application development project. Palu graduated from the University of Auckland, New Zealand, double majoring in mathematics and computer science. He has a personal interest in applying Java and mathematics in the fields of mathematical modeling and simulations, expert systems, neural and soft computation, wavelets, digital signal processing, and control systems.

Previous article: Java in Science: Data Inter- and Extrapolation Using Numerical Methods of Polynomial Fittings, Part 2
Next article: Java in Science: Data Inter- and Extrapolation Using Numerical Methods of Polynomial Fittings, Part 4

Java in Science: Data Inter- and Extrapolation Using Numerical Methods of Polynomial Fittings, Part 4
By Sione Palu

Review Part 3

Running the Complete Program
Here is the complete program. Please review the lessons of Parts 1 through 3 before using it.

import Jama.*;
import jamlab.*;

public Class Main
{
public static void main(String[] args) throws Exception
{
//actual array data pairs (xVal, yVal)
double[] xVal = {-8.,-6.,-4.,-2.,0.,2.,4.,6.,8.,10.,12.};
double[] yVal = {-4.,-7.,-2.,1.,2.,1.5,4.,5.,4.9,6.7,9.};

//the total number of polynomials to be fitted to the data pairs
int MAXIMUM_ORDER = 8;

//store the polynomial coefficients
double[][] Coef = new double[MAXIMUM_ORDER][];

//upper error-bounds
double[][] U = new double[MAXIMUM_ORDER][];

//lower error-bounds
double[][] L = new double[MAXIMUM_ORDER][];

//evaluation of best fit polynomial
double[][] yEval = new double[MAXIMUM_ORDER][];

//order of polynomial
int ORDER = 0;

//number of points for smoothing the interpolation and extrapolation of data pairs
int number_of_points = 2000;

//Array of Polyfit
Polyfit[] polyf = new Polyfit[MAXIMUM_ORDER];

//array of Polyval
Polyval[] polyv = new Polyval[MAXIMUM_ORDER];

//array of matrices for error estimation for each polynomial order
Matrix[] error = new Matrix[MAXIMUM_ORDER];

//array of matrices for polynomial evaluation for each polynomial order
Matrix[] yOutMat = new Matrix[MAXIMUM_ORDER];

//array of matrices for upper-bound error estimation for each polynomial order
Matrix[] upper = new Matrix[MAXIMUM_ORDER];

//array of matrices for lower-bound error estimation for each polynomial order
Matrix[] lower = new Matrix[MAXIMUM_ORDER];

try{
for(int i=0; i < MAXIMUM_ORDER ; i++){
ORDER = i ;
//instantiate polyfit objects with different order, starting at order=0
polyf[i] = new Polyfit(xVal,yVal,ORDER);

//array of polynomial coefficients at order = ORDER
Coef[ORDER] = polyf.getPolyCoeffMatrix().getColumnPackedCopy();
}
}
catch(Exception ex)
{ throw new Exception(ex.getMessage());}

Matrix xm = JElmat.linspace(-8,22,number_of_points);
double[] xEval = new JElmat(xm.getArrayCopy()).matrixArrayRow(0);

for(int k=0 ; k < MAXIMUM_ORDER ; k++){
polyv[k] = new Polyval(xEval,polyf[k]);

//matrix estimate of errorBounds
error[k] = polyv[k].getErrorBoundsMatrix();

//matrix evaluation of polynomial fitting
yOutMat[k] = polyv[k].getYoutMatrix();

//matrix estimate of upper-limit errorBounds
upper[k] = yOutMat[k].plus(error[k]);

//matrix estimate of lower-limit errorBounds
lower[k] = yOutMat.minus(error[k]);

//array estimate of upper-limit errorBounds
U[k] = new JElmat( upper[k].getArrayCopy()).matrixArrayRow(0);

//array estimate of lower-limit errorBounds
L[k] = new JElmat( lower[k].getArrayCopy()).matrixArrayRow(0);

//array evaluation of polynomial fitting
yEval[k] = polyv[k].getYout();
}

//Readers who use Java2D graphics for plotting, can implement the following loop
// to plot the output evaluation and error-bounds (upper and lower).
for(int j=0 ; j < MAXIMUM_ORDER ; j++) {
// (xVal,yVal) -> plot actual data pairs
// (xEval,yEval[j]) -> plot evaluation data pairs of best fit polynomial with order = j
// (xEval,L[j]) -> plot lower error-bounds for best fit polynomial of order = j
// (xEval,U[j]) -> plot upper error-bounds for best fit polynomial of order = j
}

}//-- end main --

}//--- End Class Definition ---

Interpolations
Intermediate values between data-pairs can be estimated for any polynomial order by evaluating the specific polynomial in the domain boundry with a large number of points. The data pairs used in this tutorial were :
double[] xVal = {-8.,-6.,-4.,-2.,0.,2.,4.,6.,8.,10.,12.};
double[] yVal = {-4.,-7.,-2.,1.,2.,1.5,4.,5.,4.9,6.7,9.};
There were two thousand points between the minimum , -8. and the maximum, 12, to make the estimation as smooth as possible.

The following graphs show different plottings for different order of best fit polynomials.
[i]blue line - actual data pairs
red line - best fit polynomial
black line - upper bound error estimate
yellow line - lower bound error estimate

Y = 1.918182

Y = 0.684091*X + 0.55

Y = 0.010111*X2 + 0.724534*X + 0.913986

Y = 0.000874*X3 - 0.015355*X2 + 0.672786*X + 1.031469

Y = 0.000754*X4 - 0.005157*X3 - 0.072654*X2 + 0.950233*X + 1.61049

Y = 0.000131*X5 + 0.002061*X4 + 0.006171*X3 - 0.161537*X2 + 0.739704*X + 2.27972

Y = 0.000033*X6 - 0.000525*X5 - 0.000927*X4 + 0.040596*X3 - 0.092734*X2 + 0.172263*X + 2.004319

Y = 0.000004*X7 + 0.000088*X6 - 0.000177*X5 - 0.006612*X4 + 0.033296*X3 + 0.049072*X2 + 0.205033*X + 1.506746

Extrapolations

Extending the estimation of the best fit polynomials to extrapolate outside the data pairs domain gives wildly osscilated error-bounds. Error-bounds increase as the estimation moves away from the minimum = -8 or maximum = 12 point of the data pairs domain. Error-bounds would have reasonable values around the vicinity of the minimum and the maximum.

The following graphs plot the extrapolation beyond the data pairs domain extending from minimum = -8 to maximum = 22:

Applications

Polynomial fittings is one of the numerical techniques that engineers and scientists use most often in estimation, simulation, and even model prediction. There are other numerical methods for estimation, such as splines for example. The choice depends on the domain of interest for the engineer. A NASA scientist might prefer the cubic spline as an estimator, while a financial engineering economist in the Tokyo stock market prefer to use a polynomial fitting of ORDER = 1, (straight line) for future outcome predictions.

Final Word

We have introduced the concepts of matrix algebra in this tutorial and how it is applied in numerical techniques of polynomial fittings. Linear algebra is a cornerstone of today's scientific software. There are other branches of engineering mathematics, such as discrete mathematics, differential equations, multi-variable calculus, statistics and probabilities, complex variables, and operations research, etc., which are also the foundations of scientific software. I will try to visit some of these in future tutorials.

Book References

Java for Engineers and Scientists by Stephen J. Chapman, publisher : PRENTICE HALL
Introductory Java for Scientists and Engineers by Richard J. Davies, publisher : ADDISON-WESLEY
Applied Numerical Analysis (Sixth Edition), by Curtis F. Gerald and Patrick O. Wheatly, publisher : ADDISON-WESLEY
Linear Algebra and its Applications (Second Edition), by David C. Lay, publisher : ADDISON-WESLEY
Numerical Recipes in Fortran 77, the Art of Scientific Computing (Volume 1), by William H. Press , Saul A. Teukolsky , William T. Vetterling and Brian P. Flannery, publisher : CAMBRIDGE UNIVERSITY PRESS
Mastering MATLAB 5, by Duane Hanselman and Bruce Littlefield, publisher : PRENTICE HALL
Advanced Mathematics and Mechanics Applications using MATLAB by Howard B. Wilson , publisher : CRC Press

Download Jama and Jamlab and Jamlab documentation.

References

About the Author

内容来自用户分享和网络整理，不保证内容的准确性，如有侵权内容，可联系管理员处理

标签：

相关文章推荐

新的分享

章节导航

Java in Science: Data Interpolation and Extrapolation Using Numerical Methods of Polynomial Fittings

Polynomials

Matrix Algebra

Jama and Jamlab

References

Related Articles

About the Author

Sums and Scalar Multiples

Matrix Multiplication

Transpose of a Matrix

Inverse of a Matrix

QR Matrix Factorization

Elementary Matrix Function and Special Matrices

References

Related Articles

About the Author

Polynomial Fittings

Java Codes

References

Related Articles

About the Author

Running the Complete Program

Interpolations

Extrapolations

Applications

Final Word

Book References

References

Related Articles

About the Author