EFFICIENT SOLUTION FOR CUBIC SPLINE INTEPOIATION USING LU DECOMPOSITION FIELD/BACKGROUND OF THE INVENTION [0001] The present invention discussed here relates to the problem of interpolation using cubic polynomials. The proprietary method and algorithm discussed here for the solution of the coefficients of cubic spline interpolation is unique and more efficient than the known algorithms. The proposed algorithm can be easily modeled into any computer program and used extensively in any field of science. [0002] Scientists, mathematicians and statisticians face a lot of challenges when it comes to solving the interpolation problem more efficiently and accurately. There has always been a tradeoff between efficiency and accuracy even though a variety of methods are available for interpolation like the linear interpolation, polynomial interpolation and spline interpolation. [0002a] The most commonly used method of interpolation is the linear interpolation method. It involves joining the adjacent data points using a single order linear polynomial. This method of interpolation is relatively simpler than other methods but it does not always yield an accurate result particularly when the data itself is non linear. [0002b] Similarly, polynomial interpolation is another method which tries to fit all data points using a single higher order polynomial. Usually, the polynomial connecting the given n points; is of the degree n-1. This type of interpolation has a lot of limitations and is not always feasible for highly non linear data points.
[0002c] Even though not easy to comprehend and implement, spline interpolation is considered the most accurate type of interpolation for both linear as well as non linear data. It involves connecting each of the adjacent data points using a higher degree polynomial also termed as "spline". Cubic spline interpolation is a limiting case of spline interpolation with the use of cubic polynomials (splines) to join the adjacent data points. Since the use of polynomials of a degree higher than 3, results in oscillations around the nodes, cubic polynomials give the best fit for the data. There are multiple algorithms available for solving the cubic spline interpolation problem. Next section briefly discusses these exiting algorithms and gives a detailed explanation of the improvisations done as part of this invention. DETAILED DESCRIPTION Building Equations from data set [0003] Let us consider a given data set of n points (xi, yi), (x 2 , y2)... (xn, yn) . In order to join these n points we need n-i cubic polynomials of the form y = anx3 + bnx2 + cnx + dn for each pair of the adjacent points. As a result, the first challenge is to somehow find the value for all 4(n-1) polynomial coefficients i.e. ai, bi.......an, bn, cn, dn. This can only be achieved if we can derive 4(n-1) equations out of the given set of n data points. [0003a] The first set of 2(n-1) equations are obtained by inputting the values for x and y into the polynomial equation at each of the n nodes. [0003b] Next 2(n-2) equations are derived by equating the first and second derivative of the adjacent polynomials at each of the connecting nodes. This is a fair assumption since the curve at each of the adjacent nodes, needs to connect seamlessly and without any breaks. This can only be achieved if the slope of the curve (first derivative) and slope of its slope (second derivative) equate at each of the nodes. [0003c] The last 2 equations are derived by equating the second derivative of the curve at the end points, to 0. This is a common assumption related to natural cubic spline. PROPOSED SOLUTION AND INNOVATION Existing Algorithms [0004] The second challenge is to solve the system of linear equations derived above since the number of linear equations increases linearly with the number of data points. There are multiple methods available which try to reduce the complexity of the problem by modifying the cubic polynomials itself to the form y = a' (x-xo) 3 + b' (x-xo) 2 + c' (x-xo) + d' where x 0 is the arbitrary reference point of the independent variable. After applying a little mathematics and the same constraints as described in [0003a] [0003b] and [0003c], the system of linear equations is converted to a tri-diagonal system of linear equations which can easily be solved using matrix algebra and standard algorithms. Proposed Algorithm [0005] Instead of tampering with the cubic polynomials as described above and trying to fit the linear equations to a tri diagonal system of equations, the proposed algorithm coverts the linear system of equations directly into a matrix form. The 4(n-1) equations originating from [0003a] [0003b] and [0003c] are converted into a matrix form A X = B as below - [0005a] 4(n-1) * 4(n-1) square matrix containing the values of x termed as A. [0005b] n * 1 column matrix containing all the coefficients termed as X. [0005c] n * 1 column matrix on the right hand side containing the values for y termed as B. [0006] Matrix A needs to be built in a specific way so that it can be efficiently decomposed into an upper and lower triangular matrix using LU decomposition algorithm. This is pertinent to the solution otherwise the matrix either would not decompose or will require enormous computing power to rearrange it in a way that it decomposes. The first row of matrix A is created with the coefficients of the equation which equates the second derivative at first node to 0. Next two rows are formed by inputting the values of x and y at each of the first two nodes. Next two rows are produced by equating the value of first and second derivative of the polynomials at the second node. Steps [0006b] and [0006c] are then followed for each of the nodes. Last row is obtained by equating the second derivative at the last node to 0. [0007] An optimizer algorithm is then run to check that the matrix A has valid values. This is done by verifying that that there are no duplicates, the diagonal elements are non zero and all fist row elements are non zero. In case the algorithm encounters any such anomaly, an error is thrown and if possible matrix rows are rearranged so that it can be decomposed. [0008] Once the optimizer algorithm has run, we can then decompose the matrix into a lower triangular matrix L and an upper triangular matrix U, using any standard algorithm.
[0009] Now the matrix form A X = B has been modified to (LU) X = B which can then be solved using forward backward substitution algorithm and thus yielding the value for the coefficients ai, bi.......an, bn, cn, dn. [0010] These coefficients are then inputted to the piecewise cubic polynomials producing a smooth curve joining all points of the given data set. Any value for y can thus be obtained for any given point x in between the given data set i.e. (x 1 , yi) to (xn, yn) . Example [0011] Let us assume a data set of three points (1, 2), (3, 3) and (4, 2). For these three points, we define two cubic polynomials of the form aix3 + bix2 + cix + di and a 2 x3 + b 2 x2 + c 2 x + d 2 for each interval of points. We now build eight equations using the steps mentioned in [0003a], [0003b] and [0003c]. ai + bi + ci + di = 1 Equation 1 27ai + 9bi + 3ci + di = 3 Equation 2 27a 2 + 9b 2 + 3c 2 + d 2 = 3 Equation 3 64a 2 + 16b 2 + 4c 2 + d 2 = 2 Equation 4 27ai + 6bi + ci- 27a 2 - 6b 2 + c 2 = 0 Equation 5 18ai + 2bi - 18a 2 - 2b 2 = 0 Equation 6 6ai + 2bi = 0 Equation 7 24a 2 + 2b 2 = 0 Equation 8 [0011a] We then create Matrices A, X and B from equations mentioned above and using the process described in [0006].
6 2 0 0 0 000- al- 0 1 1 1 1 0 0 0 0 b1 1 27 9 3 1 0 0 0 0 c1 3 27 6 1 0 -27 -6-1 0 dl 0 18 2 0 0 -18 -2 0 0 al 0 0 0 0 0 27 9 3 1 b2 3 0 0 0 0 64 16 4 1 c2 2 -0 0 0 0 24 2 0 0- -d2- -0 A X B [0011b] Now we run an algorithm as described in [0007] to verify that Matrix A does not have any invalid entries. [0011c] Matrix A is then broken down into a lower triangular matrix L and an upper triangular matrix U 1 0 0 0 0 0 0 0- 6 2 0 0 0 0 0 0 0.16 1 0 0 0 0 0 0 0 0.67 1 1 0 0 0 0 4.5 0 1 0 0 0 0 0 0 0 3 1 0 0 0 0 4.5-4.5 1.83 1 0 0 0 0 0 0 0 2.67 -27 -6 -1 0 3 -6 2 1.5 1 0 0 0 0 0 0 0 -22.5 7 1.5 0 0 0 0 0 1.2 1 0 0 0 0 00 0 0.6 1.2 1 0 0 0 02.84-6.5210 0 000 0 0 7.55 7.52 - 0 0 0 0 1.06 -9.111.23 1- -0 0 0 0 0 0 0 -0.18 L U [0011d] This simpler Matrix form (LU) X = B is solved using forward backward substitution method to obtain the value of polynomial coefficients. ai = -0.16667, bi = 0.5, ci = 1.166667, di = -0.5, a 2 = 0.33333, b 2 = -4, c 2 = 14.66667, d 2 = -14. [0012] These coefficients are inputted into the piecewise cubic polynomials as shown below and can be further used to get any value of y for a given value of x. Between points (1, 2) and (3, 3) -0.16667 x 3 + 0.5 x 2 + 1.166667 x - 0.5 Equation 9 Between points (3, 3) and (4, 2) 0.33333 x 3 - 4 x 2 + 14.66667 x - 14 Equations 10