NUMERICAL METHODS - CHEAT SHEET Compiled from lecture notes of Mrs.V.Valliammal by Sarath Chandar A P 1 SOLUTION OF EQUATIONS AND EIGENVALUE PROBLEMS 1.1 Solution of equation Fixed point iteration: x=g(x) method Newton’s method 1.2 Solution of linear system Gaussian elimination method Gauss-Jordan method Iterative method Gauss-Seidel method 1.3 Inverse of a matrix by Gauss Jordan method 1.4 Eigen value of a matrix by power method 1.5 Eigen value of a matrix by Jacobi method for symmetric matrix 2 INTERPOLATION AND APPROXIMATION Interpolation The process of computing the value of a function inside the given range is called interpolation. Extrapolation The process of computing the value of a function outside the given range is called extrapolation. 2.1 Newton’s forward and backward diﬀerence formulas E[f(x)] = f(x+h) ∆[f(x)] = f(x+h) - f(x) E=1+∆ Here ∆ is called Newton’s forward operator and E is called shift operator. 2 Newton’s backward diﬀerence operator is =1-E Newton’s forward interpolation formula is u(u−1) 2 u(u−1)(u−2) 3 y = y0 + u∆y0 + 2! ∆ y0 + 3! ∆ y0 + ... where u = (x - x0 ) / h The above formula is used mainly for interpolating the values of y near the begin- ning of a set of tabular values. Newton’s backward interpolation formula is u(u+1) 2 u(u+1)(u+2) 3 y = yn + u yn + 2! yn + 3! yn + ... where u = (x - xn ) / h The above formula is used mainly for interpolating the values of y near the end of a set of tabular values. 2.2 Lagrangian Polynomials Newton’s forward and backward interpolation formulae can be used only when the val- ues of dependent variables are equally spaced; where the values of independent variables are not equally spaced, then we can use following Lagrangian’s Interpolation formula and inverse interpolation formula. Lagrangian interpolation formula for unequal spaces is (x−x1 )(x−x2 )...(x−xn ) y = f(x) = (x0 −x1 )(x0 −x2 )...(x0 −xn ) y0 + (x−x0 )(x−x2 )(x−x3 )...(x−xn ) (x1 −x0 )(x1 −x2 )(x1 −x3 )...(x1 −xn ) y1 + (x−x0 )(x−x1 )(x−x3 )...(x−xn ) (x2 −x0 )(x2 −x1 )(x2 −x3 )...(x2 −xn ) y2 + (x−x0 )(x−x1 )(x−x2 )...(x−xn−1 ) ... + (xn −x0 )(xn −x1 )(xn −x2 )...(xn −xn−1 ) yn Lagrangian inverse interpolation formula for unequal spaces is (y−y1 )(y−y2 )...(y−yn ) x = f(y) = (y0 −y1 )(y0 −y2 )...(y0 −yn ) x0 + (y−y0 )(y−y2 )(y−y3 )...(y−yn ) (y1 −y0 )(y1 −y2 )(y1 −y3 )...(y1 −yn ) x1 + (y−y0 )(y−y1 )(y−y3 )...(y−yn ) (y2 −y0 )(y2 −y1 )(y2 −y3 )...(y2 −yn ) x2 + (y−y0 )(y−y1 )(y−y2 )...(y−yn−1 ) ... + (yn −y0 )(yn −y1 )(yn −y2 )...(yn −yn−1 ) xn 3 2.3 Divided diﬀerences If the values of x are given at unequal intervals then we use divided diﬀerences. Divided diﬀerences take into consideration the change of the values of the function f(x) and also the changes in the values of the arguments of x. Suppose the function y = f(x) assume the values of f(x0 ), f(x1 ), f(x2 ),....f(xn ) re- spectively, where the intervals x1 − x0 , x2 − x1 ,..... xn − xn−1 need not be equal. First divided diﬀerence of f(x) f (x1 )−f (x0 ) f(x0 , x1 ) = x1 −x0 similarly, f (x2 )−f (x1 ) f(x1 , x2 ) = x2 −x1 f (xn )−f (xn−1 ) f(xn−1 , xn ) = xn −xn−1 Second divided diﬀerence is f (x1 ,x2 )−f (x0 ,x1 ) f(x0 , x1 , x2 ) = x2 −x0 Third divided diﬀerence is f (x1 ,x2 ,x3 )−f (x0 ,x1 ,x2 ) f(x0 , x1 , x2 , x3 ) = x3 −x0 The above steps are called divided diﬀerences of order 1,2,3 respectively. Property: The divided diﬀerences are symmetrical in all the arguments. i.e. the value of any diﬀerence is independent of order of arguments. Newton’s divided diﬀerence formula is f(x) = f(x0 ) + (x − x0 )f (x0 , x1 ) + (x − x0 )(x − x1 )f (x0 , x1 , x2 ) + ... + (x − x0 )(x − x1 )(x − xn−1 )f (x0 , x1 , ..xn ) 2.4 Interpolating with a cubic spline Suppose n+1 data points are given we have to ﬁnd the value of y corresponding to x where xi x xi+1 where i = 0,1,2,.. n+1. In this case we have to use cubic spline method. 1. f(x) is a cubic polynomial 2. f (x) is second degree polynomial. 4 3. f (x) is a ﬁrst degree polynomial. 4. f(x), f (x), f (x) are continuous at each point (xi , yi ) where i = 0,1,2,..n. If the above all conditions are satisﬁed, then f(x) follows cubic spline. cubic spline formula The cubic spline in xi−1 x xi is h2 Mi−1 h2 Mi f(x) = 1 6h −x (xi − x)3 Mi−1 + (x − xi−1 )3 Mi + xih yi−1 − 6 + x−xi−1 yi − h 6 where h = xi − xi−1 for all i 6 Mi−1 + 4Mi + Mi+1 = h2 [yi−1 − 2yi + yi+1 ] ∀i here i = 1,2,3,.. n-1 M0 = 0 Mn = 0 f (xi ) = Mi 3 NUMERICAL DIFFERENTIATION AND INTEGRATION Numerical diﬀerentiation is a process of calculating the derivative of the given function by means of a table of given values of that function. In this process, we have to calculate dy d2 y dx , dx2 at a given point is called numerical diﬀerentiation. In this case, we have to use Newton’s forward and backward diﬀerence formula to compute the derivatives. In the previous units, we are ﬁnding the polynomial curves y = f(x) passing through n+1 points. But in this unit we have to calculate the derivative of such curves at a par- ticular point say Xk . If the values of X are equally spaced, we get the interpolating polynomial using Newton’s formula. If the derivative required at a point nearer to starting value in table, use newton’s forward method and for value at end of the table use newton’s backward method. In the case of unequal intervals, we can use, Newton’s divided diﬀerence formula or Lagrangian interpolation formula. 3.1 Diﬀerentiation using interpolation formulae Newton’s forward diﬀerence formula is u(u−1) 2 u(u−1)(u−2) 3 y = y0 + u∆y0 + 2! ∆ y0 + 3! ∆ y0 + ... 5 where u = (x - x0 ) / h dy 1 dx at(x = x0 ) = h 1 1 ∆y0 − 2 ∆2 y0 + 3 ∆3 y0 − 1 ∆4 y0 + ... 4 d2 y 1 11 4 13 5 dx2 at(x = x0 ) = h2 ∆2 y 0 − ∆3 y 0 + 12 ∆ y0 − 16 ∆ y0 + ... d3 y 1 dx3 at(x = x0 ) = h3 ∆3 y0 − 3 ∆4 y0 + ... 2 Newton’s backward diﬀerence formula is u(u+1) 2 u(u+1)(u+2) 3 y = yn + u yn + 2! yn + 3! yn + ... where u = (x - xn ) / h dy 1 1 2 1 3 1 4 dx at(x = xn ) = h yn + 2 yn + 3 yn + 4 yn + ... d2 y 1 2 3 11 4 dx2 at(x = xn ) = h2 yn + yn + 12 yn + ... d3 y 1 3 4 dx3 at(x = xn ) = h3 2 yn + ... 3.2 Numerical integration by trapezoidal and Simpson 1/3 and 3/8 rules The term numerical integration is the numerical evaluation of a deﬁnite integral say b a f(x) dx where a and b are given constants and f(x) be a given function. When we apply numerical integration to any function of two independent variables, using double integration is called mechanical cubature. Trapezoidal rule The follwoing rule is said to be trapezoidal rule. h I= 2 [(y0 + yn ) + 2(y1 + y2 + ... + yn−1 )] h is called length of the interval. U h= no.ofL−LL intervals Simpson’s 1/3 rule is h I= 3 [(y0 + yn ) + 2(y2 + y4 + ...) + 4(y1 + y3 + ...)] Note: In simpson’s 1/3 rule, y(x) is a polynomial of degree 2. To apply this rule the number of intervals should be even or the number of ordinates must be odd. In simpson’s 1/3 rule, the truncation error is of the order h4 . i.e. (b−a)h4 M |E| 180 iv iv iv where M = max(y0 , y1 , y2 , ..) 6 Simpson’s 3/8 rule is 3h I= 8 [(y0 + yn ) + 3(y1 + y2 + y4 + y5 + ... + yn−3 ) + 2(y3 + y6 + ...)] Note: In simpson’s 3/8 rule, y(x) is a polynomial of degree 3. To apply this rule the number of intervals must be a multiple of 3. 3.3 Romberg’s method b Let I = a f(x) dx Using trapezoidal or simpson’s rule, let us ﬁnd the value of the given deﬁnite inte- gral by taking the diﬀerent values of h or ∆x. Let us apply trapezoidal rule for several times by the values as h, h/2, h/4, h/8, h/16... In this type of problems, let I2 −I1 I = I2 + 3 where I1 = the value of the deﬁnite integral, dividing h into two parts. Similarly I2 = the value of the deﬁnite integral, dividing h into four parts. Similarly we can calculate the values of I3 −I2 I = I3 + 3 This method can continue till we get two successive values are equal. 3.4 Two and Three point Gaussian quadrature formulae Two point formula In this case change the given interval a to b into -1 to 1 using the following procedure. Transformation formula is a+b b−a x= 2 + 2 t where b is the upper limit and a is the lower limit. 1 √ √ −1 f(t) dt = f(1/ 3) + f(-1/ 3) Three point formula In this case change the given interval a to b into -1 to 1 using the following procedure. Transformation formula is a+b b−a x= 2 + 2 t 7 where b is the upper limit and a is the lower limit. 1 −1 f(t) dt = A1 f(t1 ) + A2 f(t2 ) + A3 f(t3 ) where A1 = A3 = 0.5555 A2 = 0.88888 t1 = −0.7745 t2 = 0 t3 = +0.7745 3.5 Double integral using trapezoidal and Simpson’s rules trapezoidal rule is I = hk [sumof thef ourcornervalues(circle)+ 4 2(sum of values available in rectangular box) + 4(sum of remaining values) ] simpson rule is I = hk [sumof thef ourcornervalues(circle)+ 9 4(sum of values available in rectangular box) + 16(sum of remaining values) ] 4 INITIAL VALUE PROBLEMS FOR ODE There are several methods for solving diﬀerential equations numerically. 1. Taylor’s series method 2. Euler’s method 3. Rungi-kutta method 4. Milne’s predictor and corrector method 5. Adam’s predictor and corrector method 4.1 Single step methods Taylor’s series method Taylor’s series about x = x0 is given by (x−x0 )2 (x−x0 )3 (x−x0 )4 y(x) = y0 + (x - x0 ) y0 + 2! y0 + 3! y0 + 4! y0 v where n0 , y0 denotes the initial values of x and y Euler method for ﬁrst order equation eulers method formula or eulers algorithm. yn+1 = yn + h f(xn , yn ) where n=0,1,2... modiﬁed eulers method formula : yn+1 = yn + h f xn + h ; yn + h f (xn , yn ) 2 2 8 Fourth order Runge-Kutta method for solving ﬁrst and second order equa- tions in this method , let us assume that dy dx = f(x,y) be given diﬀ eqn, the above eqn can be solved under the condition y(xo )= y0 let the increment h = length of the interval between two values fourth order rangi kutta method K1 = h f(x0 , y0 ) k1 K2 = h f(x0 + h , y0 + 2 2 ) k2 K3 = h f(x0 + h , y0 + 2 2 ) K4 = h f(x0 + h, y0 + k3 ) (K1 +2K2 +2K3 +K4 ) ∆y = 6 y1 = y0 + ∆y x1 = x0 + h The increments in y for the second interval is computed in a similar manner by the following formulae. K1 = h f(x1 , y1 ) k1 K2 = h f(x1 + h , y1 + 2 2 ) k2 K3 = h f(x1 + h , y1 + 2 2 ) K4 = h f(x1 + h, y1 + k3 ) (K1 +2K2 +2K3 +K4 ) ∆y = 6 y2 = y1 + ∆y x2 = x1 + h 4.2 Multistep methods Milne’s predictor and corrector methods 1. To use Milne’s predictor and corrector formula, we need atleast 4 values. 9 2. If the initial values are not given, we can obtain these values either by using Taylor’s series method or by R-K method. 3. We can apply milne’s corrector formula only after applying milne’s predictor for- mula. Milne’s predictor formula is 4h yn+1,p = yn−3 + 3 2yn−2 − yn−1 + 2yn Milne’s corrector formula is h yn+1,c = yn−1 + 3 yn−1 + 4yn + yn+1 Adam’s predictor and corrector methods Adam’s predictor formula is h yn+1,p = yn + 24 55yn − 59yn−1 + 37yn−2 − 9yn−3 Adam’s corrector formula is h yn+1,c = yn + 24 9yn+1 + 19yn − 5yn−1 + yn−2 5 BOUNDARY VALUE PROBLEMS IN ODE AND PDE to be included later.