INTRODUCTION TO NUMERICAL ANALYSIS

1 Cho, Hyoung KyuDepartment of Nuclear EngineeringSeoul National UniversityCho, Hyoung KyuDepartment of Nuclear EngineeringSeoul National UniversityINTRODUCTION TO NUMERICAL ANALYSIS10. NUMERICAL Euler's Modified Euler's Midpoint Runge Kutta Multistep Predictor Corrector System of First Order Ordinary Differential Solving a Higher Order Initial Value Use of MATLAB Built In Functions for Solving Initial Value Local Truncation Error in Second Order Range Kutta Step Size for Desired Stiff Ordinary Differential Background Ordinary differential equation A differential equation that has one independent variable A first order ODE The first derivative of the dependent variable with respect to the independent variable Example Rates of water inflow and outflow The time rate of change of the mass in the tank Equation for the rate of height Background Time dependent problem Independent variable: time Dependent variable.

2 Water level To obtain a specific solution, a first order ODE must have an initial condition or constraint that specifies the value of the dependent variable at a particular value of the independent variable. In typical time dependent problems Initial condition Initial value problem (IVP) Background First order ODE statement General form Ex) Flow lines Analytical solution In many situations an analytical solution is not possible! NUMERICAL solution of a first order ODE A set of discrete points that approximate the function y(x) Domain of the solution: , N Background Overview of NUMERICAL methods used/or solving a first order ODE Start from the initial value Then, estimate the value at a second nearby point third point .. Single step and multistep approach Single step approach: Multistep approach: , , , Explicit and implicit approach Right hand side in explicit method: known values Right hand side in implicit method: unknown value In general, non linear equation non linear equation solution method !

3 Implicit methods provide improved accuracy over explicit methods , but require more effort at each Background Errors in NUMERICAL solution of ODEs Round off errors Truncation errors NUMERICAL solution of a differential equation calculated in increments (steps) Local truncation error: in a single step Propagated, or accumulated, truncation error Accumulation of local truncation errors from previous steps Single step explicit methods Euler s explicit method: slope at , Modified Euler s explicit method: average slope at , and , Midpoint method: slope at /2 Runge Kutta methods A weighted average of estimates of the slope of at several Euler s Method Euler s method Simplest technique for solving a first order ODE Explicit or implicit Euler's Explicit Method The error in this method depends on the value of and is smaller for smaller h.

4 Derivation NUMERICAL integration or finite difference approximation of the derivativestep size is exaggerated !(rectangle method)(forward Euler method) Euler s Method Example 10 1: Solving a first order ODE using Euler's explicit Euler s Method Example 10 1: Solving a first order ODE using Euler's explicit Solving Example 8-1clear alla=0; b= ; h= ; yINI = 3;[x, y] = ylabel('y')function[x, y] = odeEULER(ODE,a,b,h,yINI)x(1) = a; y(1) = yINI;N = (b-a)/h;fori = 1:Nx(i+1) = x(i) + h;y(i+1) = y(i) + ODE(x(i),y(i))*h;endfunctiondydx = Chap8 Exmp1 ODE(x,y)dydx = *y + 7*exp( *x); Euler s Method ANALYSIS of truncation error in Euler's explicit method Local truncation error Global truncation error Taylor series expansion at position 1 NUMERICAL solution Local truncation error Total error The difference between the NUMERICAL solutionand the true Euler s Method ANALYSIS of truncation error in Euler's explicit method Global truncation error Truncation error is propagated or accumulated !

5 Mean value Euler s Method ANALYSIS of truncation error in Euler's explicit method Global truncation error Suppose Propagation of error At the first point to the second point to the third point At the fourth point At the point Suppose Euler s Method ANALYSIS of truncation error in Euler's explicit method Global truncation error Suppose Difficult to determine the order of magnitude directly Possible to determine the bound Euler s Method Euler's implicit method In general, this equation is non linear! Must be solved with a NUMERICAL solution method In the derivation Backward difference formula for the derivative backward Euler method The local and global truncation errors Same as those in the explicit Euler s Method Example 10 2: Solving a first order ODE using Euler's implicit methodTo solve the non linear equationusing Newton methodIteration Euler s Method Example 10 2: Solving a first order ODE using Euler's implicit method% Solving First Order ODE with Euler's implicit alla = 0; b = ; h = ;N = (b - a)/h;n(1) = 2000; t(1) = a;fori=1:Nt(i+1) = t(i) + h;x = n(i);% Newton's method starts.

6 Forj = 1:20num = x + *x^(3/2)*h - *n(1)..*(1-exp(-3*t(i+1)))*h - n(i);denom = 1 + * *x^(1/2)*h;xnew = x - num/denom;ifabs((xnew - x)/x) < = xnew;endendifj == 20fprintf(' NUMERICAL solution could no be calculated at t = %g s', t(i))breakend% Newton's method ends. n(i+1) = xnew;endplot(t,n)axis([0 0 2000]), xlabel('t (s)'), ylabel('n') Modified Euler s Method Modified Euler s Method Main assumption in Explicit method Constant derivative (slope) between , and , Equal to the derivative at point , Modified Euler method To include the effect of slope changes within the subinterval Uses the average of the slope at points , and , Slope at the beginning: at point , Using Euler s explicit method Estimating slope at the end of interval: at point , main source of Modified Euler s Method Modified Euler s Method Better estimation of with new slope Also derived by integrating ODE using the trapezoidal method Modified Euler s Method Example 10 3: Solving a first order ODE using the modified Euler [x, y] = odeModEuler(ODE,a,b,h,yINI)x(1) = a.

7 Y(1) = yINI;N = (b-a)/h;fori = 1:Nx(i+1) = x(i) + h;SlopeEu = ODE(x(i),y(i));yEu = y(i) + SlopeEu*h;SlopeEnd = ODE(x(i+1),yEu);y(i+1) = y(i) + (SlopeEu+SlopeEnd)*h/2;endExplicit EulerModified Midpoint Method Midpoint method Another modification of Euler s explicit method Slope at the middle point of the interval Two slope at the midpoint , NUMERICAL solution : Runge Kutta methods Runge Kutta methods A family of single step, explicit, NUMERICAL techniques for a first order ODE Slope: obtained by considering the slope at several points within the subinterval Different orders of Runge Kutta method the number of points within the subinterval Second order RK: two points Third order RK: three points Fourth order RK (classical RK): four points Global truncation error Second order RK: second order accurate globally.

8 Third order accurate locally More accurate than simple Euler s explicit method However, require several evaluations of function for the derivative Runge Kutta methods Second order Runge Kutta methods General form The values of these constants vary with the specific second order method. Modified Euler method and the midpoint method Two versions of a second order RK method Modified Euler method: , , 1, 1 Midpoint method: 0, 1, , Runge Kutta methods Second order Runge Kutta methods Heun s method , , , Truncation error in second order RK methods Local truncation error: Global truncation error: A larger step size can be used for the same accuracy ! However, the function , is calculated Runge Kutta methods Second order Runge Kutta methods and Talyor series expansion First and second derivatives: given by the differential equation Then, General form of Runge Kutta methods Second order Runge Kutta methods and Talyor series expansion Two different equations for Three equations with four unknownsModified Euler: , , 1, 1 Midpoint method: 0, 1, , Heun s: , , , Runge Kutta methods Example 10 4.

9 Solving by hand a first order ODE using the second order RungeKutta Runge Kutta methods Third order RK methods General form , , , , , , , Four term Talyor series expansion Classical third order RK method , , , , 1, , 1, Runge Kutta methods Third order RK methods Truncation error Local truncation error: Global truncation error: Other third order RK Runge Kutta methods Fourth order RK methods General form 13 constants: , , , , , , , , , , , , Classical fourth order RK Runge Kutta methods Fourth order RK Runge Kutta methods Fourth order RK methods Truncation error Local truncation error: Global truncation error: Runge Kutta methods Example 10 5: Solving by hand a first order ODE using the fourth order Runge Kutta method Example 10 6: A user defined function for solving a first order ODE using the fourth order Runge Kutta Runge Kutta methods Example 10 6: A user defined function for solving a first order ODE using the fourth order Runge Kutta alla=0; b= ;h= ; yIni=3;[x,y] = )error=yExact-yfunction[x, y].

10 = odeRK4(ODE,a,b,h,yIni)x(1) = a; y(1) = yIni;n = (b-a)/h;fori = 1:nx(i+1) = x(i) + h;K1 = ODE(x(i),y(i));xhalf = x(i) + h/2;yK1 = y(i) + K1*h/2;K2 = ODE(xhalf,yK1);yK2 = y(i) + K2*h/2;K3 = ODE(xhalf,yK2);yK3 = y(i) + K3*h;K4 = ODE(x(i+1),yK3);y(i+1) = y(i) + ..(K1 + 2*K2 + 2*K3 + K4)*h/6;endfunctiondydx = Chap8 Exmp6 ODE(x,y)dydx = *y + 7*exp( *x); Multistep methods Single step method Use only the value of and at the previous point Multi step method Two or more previous points Explicit and implicit methods Explicit: The first few points can be determined by single step methods or by multistep methods that use fewer prior points. Implicit method: appears on both sides. non linear equation solution method! Adams Bashforth Method Adams Moulton Multistep methods Adams Bashforth Method Explicit multistep method for solving a first order ODE Second order formula uses , and , Third order formula uses , , , , , Integration of ODE Integration is carried out with a polynomial that interpolates the value of , at , and at previous Multistep methods Adams Bashforth Method Second order AB method First order polynomial Integration Third order AB method Fourth order AB methodEuler s explicit Multistep methods Adams Moulton Method Implicit multistep method for solving first order ODEs Second order formula Implicit form of the modified Euler method Third order and fourth formulas Two ways to use AM method If they are used by themselves, they have to be solved numerically.