1 / 30

MATLAB 입문 CHAPTER 8 Numerical Calculus and Differential Equations

MATLAB 입문 CHAPTER 8 Numerical Calculus and Differential Equations. Numerical Methods for Differential Equations. The Euler Method Consider the equation (8.5-1) r(t) : a known function From the definition of derivative,

daire
Download Presentation

MATLAB 입문 CHAPTER 8 Numerical Calculus and Differential Equations

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. MATLAB 입문CHAPTER 8 Numerical Calculus and Differential Equations ACSL, POSTECH

  2. Numerical Methods for Differential Equations • The Euler Method • Consider the equation (8.5-1) r(t) : a known function • From the definition of derivative, ( is small enough) (8.5-2) • Use (8.5-2) to replace (8.5-1) by the following approximation: (8.5-3)

  3. Numerical Methods for Differential Equations • Assume that the right side of (8.5-1) remains constant over the time interval . Then equation (8.5-3) can be written in more convenient form as follows: (8.5-4) , where : step size • The Euler method for the general first-order equation is (8.5-5) • The accuracy of the Euler method can be improved by using a smaller step size. However, very small step sizes require longer runtimes and can result in a large accumulated error because of round-off effects.

  4. Euler method for the free response of dy/dt = –10y, y(0) = 2 r= -10; deltaT = 0.01 ; % step size t=[0:deltaT:0.5]; % our vector of time values y(1)=2; % initial value at t=0 for k=1: length(t)-1 % for every value in t vector y(k+1) = y(k) + r*y(k)*deltaT; end y_true = 2*exp(-10*t); plot(t,y,'o',t,y_true), xlabel('t'), ylabel('y');

  5. Euler method solution for the free response of dy/dt= –10y, y(0) = 2. Figure 8.5–1 8-19

  6. Euler method solution of dy/dt= sin t, y(0) = 0. Figure 8.5–2 8-20 More? See pages 490-492.

  7. Numerical Methods for Differential Equations • The Predictor-Corrector Method • The Euler method can have a serious deficiency in problems where the variables are rapidly changing because the method assumes the variables are constant over time interval . (8.5-7) (8.5-8) • Suppose instead we use the average of the right side of (8.5-7) on the interval . (8.5-9) , where (8.5-10)

  8. Numerical Methods for Differential Equations • Let Euler predictor: (8.5-11) Trapezoidal corrector: (8.5-12) • This algorithm is sometimes called the modified Euler method. However, note that any algorithm can be tried as a predictor or a corrector. • For purposes of comparison with the Runge-Kutta methods, we can express the modifided Euler method as (8.5-13) g1 is deltaY given by f(tk,yk) (8.5-14) g2 is deltaY given by f(tk+1,yk+1) (8.5-15)

  9. Modified Euler solution of dy/dt = –10y, y(0) = 2 r= -10; deltaT = 0.01 ; % same step size as before y(1)=2; t=[0:deltaT:0.5]; % our vector of time values for k=1:length(t)-1 g1 = deltaT*r*y(k); % estimate of deltaY at t(k),y(k) g2 = deltaT*r*(y(k) + g1); % est of deltaY at t(k+1),y(k+1) y(k+1) = y(k) + 0.5*(g1 + g2); end y_true = 2*exp(-10*t); plot(t,y,'o',t,y_true), xlabel('t'), ylabel('y');

  10. Modified Euler solution of dy/dt= –10y, y(0) = 2. Figure 8.5–3 8-21

  11. Modified Euler solution of dy/dt= sin t, y(0) = 0. Figure 8.5–4 8-22 More? See pages 493-496.

  12. Numerical Methods for Differential Equations • Runge-Kutta Methods • The second-order Runge-Kutta methods: (8.5-17) , where : constant weighting factors (8.5-18) (8.5-19) • To duplicate the Taylor series through the h2 term, these coefficients must satisfy the following: (8.5-19) (8.5-19) (8.5-19)

  13. Numerical Methods for Differential Equations • The fourth-order Runge-Kutta methods: (8.5-23) (8.5-24) Apply Simpson’s rule for integration

  14. Numerical Methods for Differential Equations • MATLAB ODE Solvers ode23 and ode45 • MATLAB provides functions, called solvers, that implement Runge-Kutta methods with variable step size. • The ode23 function uses a combination of second- and third-order Runge-Kutta methods, whereas ode45 uses a combination of fourth- and fifth-order methods. <Table 8.5-1> ODE solvers

  15. Stiff? • Stiff Differential Equations • A stiff differential equation is one whose response changes rapidly over a time scale that is short compared to the time scale over which we are interested in the solution. • A small step size is needed to solve for the rapid changes, but many steps are needed to obtain the solution over the longer time interval, and thus a large error might accumulate. • The four solvers specifically designed to handle stiff equations: ode15s (a variable-order method), ode23s (a low-order method), ode23tb (another low-order method), ode23t (a trapezoidal method).

  16. Numerical Methods for Differential Equations • Solver Syntax <Table 8.5-2> Basic syntax of ODE solvers

  17. Lets try these on the problem we've been solving

  18. Numerical Methods for Differential Equations <example>Response of an RC Circuit (RC=0.1s, v(0)=0V, y(0)=2V)

  19. Numerical Methods for Differential Equations • Effect of Step Size • The spacing used by ode23 is smaller than that used by ode45 because ode45 has less truncation error than ode23 and thus can use a larger step size. • ode23 is sometimes more useful for plotting the solution because it often gives a smoother curve. • Numerical Methods and Linear Equations • It is sometimes more convenient to use a numerical method to find the solution. • Examples of such situations are when the forcing function is a complicated function or when the order of the differential equation is higher than two. • Use of Global Parameters • The global x y z command allows all functions and files using that command to share the values of the variables x, y, and z.

  20. Extension to Higher-Order Equations • The Cauchy form or the state-variable form • Consider the second-order equation (8.6-1) (8.6-2) • If , , then The Cauchy form or state-variable form

  21. Extension to Higher-Order Equations <example>Solve (8.6-1) for with the initial conditions . ( Suppose that and use ode45.) %Define the function function xdot=example1(t,x) xdot(1)=x(2); xdot(2)=(1/5) * (sin(t) - 4*x(1) - 7*x(2)); xdot = [xdot(1) ; xdot(2)]; %then, use it: [t, x] = ode45('example1',[0,6],[3,9]); plot(t,x); The time interval of interest The initial condition for the vector x. xdot(1) xdot(2) x(1) x(2)

  22. Extension to Higher-Order Equations • Matrix Methods < a mass and spring with viscous surface friction > (8.6-6) ( Letting , ) (Matrix form) (8.6-7) (Compact form)

  23. The following is extra credit (2 pts)

  24. Extension to Higher-Order Equations • Characteristic Roots from the eig Function • Substituting , • Cancel the terms • Its roots are s = -6.7321 and s = -3.2679. (8.6-8) (8.6-8) A nonzero solution will exist for A1 and A2 if and only if the deteminant is zero

  25. Extension to Higher-Order Equations • MATLAB provides the eig function to compute the characteristic roots. • Its syntax is eig(A) <example> The matrix A for the equations (8.6-8) and (8.6-9) is • To find the time constants, which are the negative reciprocals of the real parts of the roots, you type tau = -1./real (r). The time constants are 0.1485 and 0.3060.

  26. Extension to Higher-Order Equations • Programming Detailed Forcing Functions < An armature-controlled dc motor> ▪ Apply Kirchhoff’s voltage low and Newton’s low (8.6-10) (8.6-11) : motor’s current , : rotational velocity ( matrix form ) Letting L: inductance, R: resistance, I: inertia, : torque constant, : back emf constant, c: viscous damping constant, : applied voltage

More Related