New User Special Price Expires in

Let's log you in.

Sign in with Facebook


Don't have a StudySoup account? Create one here!


Create a StudySoup account

Be part of our community, it's free to join!

Sign up with Facebook


Create your account
By creating an account you agree to StudySoup's terms and conditions and privacy policy

Already have a StudySoup account? Login here


by: Mossie Pfannerstill V


Mossie Pfannerstill V
GPA 3.59

David George

Almost Ready


These notes were just uploaded, and will be ready to view shortly.

Purchase these notes here, or revisit this page.

Either way, we'll remind you when they're ready :)

Preview These Notes for FREE

Get a free preview of these Notes, just enter your email below.

Unlock Preview
Unlock Preview

Preview these materials now for free

Why put in your email? Get access to more of this material and other relevant free materials for your school

View Preview

About this Document

David George
Class Notes
25 ?




Popular in Course

Popular in Applied Mathematics

This 12 page Class Notes was uploaded by Mossie Pfannerstill V on Wednesday September 9, 2015. The Class Notes belongs to AMATH 584 at University of Washington taught by David George in Fall. Since its upload, it has received 12 views. For similar materials see /class/192272/amath-584-university-of-washington in Applied Mathematics at University of Washington.

Similar to AMATH 584 at UW

Popular in Applied Mathematics




Report this Material


What is Karma?


Karma is the currency of StudySoup.

You can buy or earn more Karma at anytime and redeem it for class notes, study guides, flashcards, and more!

Date Created: 09/09/15
ODE Lecture Notes Amath 584 David George University of Washington Fall 2007 1 Introduction In these notes we will look at one type of numerical method nite difference methods for one problemea rst order initial value problem of the form 275 fy7t7 1a 24050 2407 1b where in general y E C is a vector and fyt E C is a vector valued possibly nonlinear functionl The goal will be to nd a numerical approximation to the solution yt over some xed interval 25 t0T as a discrete set of values Y E C i 0N where Y is a pointwise approximation to This means that we have broken the domain 25 250 T into a set of points 250251 tN with tN T and we nd an approximation to the continuous solution at those discrete points see Figure 1 Note that the system 1 is more general that it might rst appear since you can turn a system of ODEs with higher order derivatives into a rst order system by adding components to the solution vector For instance consider the following example of a set of two second order ODEs m t hx z mz 2a z t gz z mz 2b where z and z are scalar functions with initial conditions z0 a 20 b m 0 c and z 0 c By letting 241 95 242 7 Z 9 ya i m 7 3 244 2 we can derive a rst order system the governs the vector yt 2 fy7 4 1We will assume that f satis es properties such that our solution is unique and multiply di erentiable etc etc but not concern ourselves with such issues here yt Figure 1 Numerical approximation for yt is the set of values Yi i O7 7N7 where Y is the discrete approximation to 3027 for i O7 N where 95 242 2 7 ya hm7z7m7Z hy 7 5 9907271 72 924 with initial conditions a 1 gm 7 G lt6 1 So by being able to solve 17 we really can solve any higher order nonlinear system Additionally7 we can eliminate the direct dependence of y t on 257 by letting t be a component of the solution vector y E C That is7 we add another component to y which is t The governing ODE for that last component is just ym y t 1 This is not necessary to use the numerical methods we ll discuss and is in fact rarely done numerically However it will save a lot of typing by letting us focus on the system where the direct dependence of f on t has been removed Of course in reality the system still depends on time but that dependence occurs through the in uence of the last component of 2 Problem Stability Just as we were interested in the condition of a problem in linear algebra which is not a property of the numerical method we will be interested in problem stability77 in ODEs Problem stability is not a property of the numerical method but a property of the ODEs themselves The accuracy of our numerical method will depend on problem stability just as in many linear algebra problems our relative error was related to the condition of a matrix Problem stability indicates how a perturbation to the initial problem 7 changes the solution That is suppose that 7 is changed by adding a small perturbation e to the initial data 205 y 8a 24050 yo 6 8b If we let 3725 be the solution to the perturbed initial value problem 8 and let gt be the solution to the original problem 7 then if MW 7 WM 9 grows in time the problem is said to be unstable and if it shrinks as time proceeds the problem is stable The norm can be any norm due to norm equivalence Note that problem refers not only to the governing ODE but also the initial condition Some ODEs can be stable for some initial condition and unstable for others A real scalar problem is stable at a given state yt if 1 lt 07 10 at that state Similar expression can be derived if y is more general If y E C then only the real part of y need be considered assuming 1 is analytical and needs to be less than 0 for stability If y E C and if the Jacobian matrix y is diagonalizable the real part of all eigenvalues should be considered and need to be less than zero for problem stability For problems where the Jacobian is defective stability can get much more complicated and we will totally ignore such casesiexcept maybe on the exam just kidding The reason why this implies stability is explained heuristically below If we let 2t 3725 7 yt then 2t obeys the governing lVP Zt Hi H137 11a 2t0 6 11b If the norm of the solution to 11 grows in time the original problem 7 is unstable if it damps to zero in time the original problem is stable and if it remains constant the original problem is sometimes said to be neutrally stable77 In order to analyze 11 for nonlinear functions we will do a form of linear stability analysis The idea is as follows We know yt0 it is just the initial condition of the original problem We want to consider arbitrary yet small values for e We then do a Taylor expansion of the right hand side of the ODE in 11 about the state g which gives 1 7 Hz Ham 7 11 022 mm 022 12 Again in the case of systems where y is a vector the term 1quot refers to the Jacobian matrix in the Taylor expansion of 1 Now if e is small enough we can ignore 022 initially which implies that 2t is initially governed by something close to 2 05 HEM 13a 2t0 6 13b This is a linear system for z resembling zt z 14a 2to 20 6 1410 the solution to which is of course 2t zoeMt tO lf 9 lt 0 then 2t decays to zero meaning that y approaches 3 and 022 terms in 12 can continue to be neglected Now y is really a function of time in general since y does not need to be stationary So problem stability can change as the solution evolves if the sign of the real part of y changes For problems where y is constant such as if it is a critical point of f stability doesn t change as the solution evolves If y E C we might write 13 as zt A2 15a 2t0 20 6 15b where A E mem is the m by m Jacobian y The solution to 15 is 2t eAtzo 16 where 5quot E mem is the matrix exponential de ned by the series eAtIAtw m 17 If A is diagonalizable Ak RAR Uk RAR 1RAR 1RAR 1 RAkR l 18 since all of the R lR terms are the identity We can therefore write the series for the matrix exponential as A2252 Aktk 5ARIAtTHTHR 17 19a 22 kk RileAtRIAtA2t 1AT55M 19b The signi cance of this is that e is diagonal since all the terms in its series are diagonal the products of diagonal matrices are still diagonal That is eAt t 20 where 7 is the jth eigenvalue of A So we have 2t eAtzo R712t RileAtzo wt RileAtRwo 21a wt eAtwo 2110 where we R lzo In other words by changing to the basis of eigenvectors we see that each component w7 of the solution obeys the scalar equation Mt e vtwg 22 Therefore if the real part of all of the eigenvalues of 1quot 37 have negative real parts each component of wt decays so 2t decays and 37 approaches 3 When looking at a general nonlinear system then it is the eigenvalues of the Jacobian at any state y that determine problem stability When doing linear stability analysis terms such as small when referring to e depend on the relative size of rst and higher order terms in the Taylor expansion 12 which depends on the speci c problem at hand since it depends on the size of f g in relation to y etc etc Don t be too concerned about whether this analysis is always valid It depends on the problem For some problems linear stability analysis is sufficient and others it may not be A good example is when the real part of all eigenvalues is zero Stability then depends on higher order terms in the Taylor expansion For our purposes we are content to examine how our numerical methods will behave on stable systems This section therefore is really just an explanation of why later we will focus on the test problem wt u ut0 u If we consider how a numerical method treats 23 it gives us an indication of how it will treat more general problems where the true solution is perturbed from gt to 3725 at any time 3 Finite Difference Methods Finite difference methods provide formulas to generate a numerical sequence U07 U17 U27 7 UN7 24 as an approximation to the solution ut of the ODE Each Un is an approximation to utn and t0 tN to T the time interval we want the solution over There are various ways to derive a given nite difference method for a set of ODEs Since 25 is a rst order system one way to derive a given method is to use quadratureiformulas for approximating to integrals That is note that for any At tAt tAt ttAt utAtiutt uTdTt fuTdTt Fde 26 We can derive various nite difference methods by approximating in various ways the last integral in 26 For instance if we look at two points tn1 and tn we have ulttn1gte um fut dt 3 Amman lt27 where Atn tag 7 tn Of course the smaller Atn is the closer the quadrature approximations are to the integrals Unless otherwise stated assume that each Atn is the same so we will just call it At However this is not strictly necessary in any of the methods we will look at This quadrature approximation gives rise to a numerical method known as Forward Euler Un1 Un AtfUn7 28a U0 ut0 28b Note that we just start with U0 and then build the sequence U1 U2 UN Since nding Un only requires evaluating fUn the method is known as explicit Another quadrature approximation gives rise to Backward Euler Un1 Un AtfUn17 29a U0 Ultt0 Note that now we need to solve 29a for Un1 at each time step Unless u is linear this may require a root solving iteration The solution Un is the root of the equation CUW E UM1 7 U 7 AtfUn1 30 If we use Newton iteration this requires nding Gq which we can only do if we can nd fq For systems this requires nding the Jacobian of Cq analytically on paper Be sure and realize that Un is a constant when nding the Jacobian If the function f is too complicated to be able to do this some other form of root nding would be needed such as secant iteration For Newton iteration you would iterate by solving G WWW Fqk7 31a qk1 qW 5 31b until some tolerance is reached so that qk is close enough to a root The test might be HGqk lt e for some small epsilon such as emachme Then set Un1 gm Note that you need some q 0 to start the iteration A good choice is q0 Un since we presume that Un will typically not be too far from the root Un especially if At is small Schemes that require solving implicit functions at each time step due to the fUn1 such as Backward Euler are called implicit lmplicit schemes therefore come with an added cost however as we will see when we talk about stability they convey some bene t and are in fact imperative for some ODEs Another implicit scheme is the Trapezoid Method UM U MU fUn17 3 U0 uto 32b The Trapezoid Method 32 would also need Newton iteration 31 The methods we have discussed so far are all Z step methods since Un1 only depends directly on U More generally we might de ne Linear Multistep Methods of the form 2 O iUnH Z W Un t 33 j0 7390 where 047 and B are constants that de ne this general r step method Even the 1 step methods we ve looked at above are Linear Multistep Methods with r 1 Note that if B 0 the method is explicit otherwise it is implicit An example of an explicit 2 step method r2 is the Midpoint Method Un1 UH 2AtfU 34 You might notice that we need to nd U1 some other way before we can start using 34 Typically some other 1 step method would be used to nd U1 and then then 34 is used for U2 U3 UN More generally an r step method must nd the rst U0 U1 UT1 starting values using other methods However since only a small number of steps are used to generate these values it doesn t change the properties of the method overall 31 Global Error The global error of our numerical method is simply the distance between the numerical solution and the actual solution in some norm It is de ned at each point tn n O 1 N That is en M 7 mm 35gt 32 Truncation Error Ultimately we are interested in the global error however it can be dif cult to ascertain the global error directly since of course we don t know the true solution generally The truncation error is another type of measure of the methods accuracy Note that our numerical solution sequence Un0 is the solution to the difference equation 2 O iUnH Z Bj UnH 0 36 j0 7390 However if we plug the true solution into 36 as 7 7 Zajutnj Z j u m t 37 j0 j0 it will not be equal to zero That is if we turn the true solution into a discrete set of points ut0ut1 utN by evaluating the true solution at the points t0t1 tN it won t satisfy the same difference equation 36 that the numerical solution satis es The truncation error is a measure of how far the difference equation is from the true ODE This is most easily explained by example Consider Forward Euler 28 To nd the truncation error we arrange 28 so that it resembles the actual ODE not 26 That is we rearrange Forward Euler as Un1 7 Un At 7 039 38 The numerical solution is the solution to this difference equation Now the truncation error is de ned by utn1 Win At 7 mm m 39 The truncation error tells us how good or bad our difference equation approximates the actual ODE Now we need to nd a more informative value for Tn The trick is to use Taylor expansions of ut about the point tn Note that ulttn1gte um At 2 u tAt u tAt 0At3 40 Because of the presence of the At term on the right hand side of 45 the left hand side is known as a rst order approximation to u tn Now we know from the ODE that futn u tn So plugging that and 45 into 39 gives u tnAt u tnAt2 2 3 OAt3 41 Tn In fact we often write that 42 since that indicates the lowest order error Since this involves At it is known as a rst order method If the lowest order term was a constant times At2 we would say that it is a second order method and so on Therefore higher order methods are more accurate since the truncation error is smaller as At gets small As one more example we ll look at the Trapezoid Method We write the method as W 7 g fUn fUn1 0 43 The truncation error is de ned by e mm 4441 lt44 We can again use 45 for the left hand side of 44 and again use the fact that futn u tn So it remains to nd futn1 The trick is to note that futn1 u tn1 and then do a Taylor expansion of u tn1 about tn This gives tn Atz 4441 mm w 44mm t OW 45 Now when we plug these expressions into 44 we get 7 tn AtZ T 0At3 46 So we see that the Trapezoid Method is second order 33 Onestep Error The truncation error is related to the onestep error The onestep error is actually easier to understand but is not as useful theoretically as the truncation error The one step error can always be found by the truncation error and some authors refer to the onestep error as the truncation error The onestep error n is the error that is made in one step of the numerical method That is if we consider Un Un71 U0 to be exact the one step error is the difference Un1 7 utn1 The norm of this is not the global error because in reality Un Un1 U0 are not exact As an example for Forward Euler the one step error is r utn17 15 7 At 1014 47 In this case the onestep error is just At times the truncation error More generally the one step error is CAtTn for some constant C for the method A naive analysis might suggest that the global error at tN 5N HUN 7 utNH is the sum of all of the onestep errors on the order of O nN O nTAt 0739 However this is only true if the one step errors do not accumulate adversely This involves the numerical stability of the method 3 4 Consistency Consistency of a numerical method merely means that the truncation error goes to zero as At a 0 So a method is said to be consistent if m OAtk k 1 48 Alternatively the onestep error must satisfy En OAt2 or better The onestep error would always be OAt or better since step size is only At 3 5 Convergence Convergence of a method means that as At a 0 UN a utN More precisely if we let to O we say that a method on a problem 25 is convergent if 1 U 7 t 0 49 AtaoflglAtT N u We say that a method is convergent if 49 is satis ed for all T and all problems 25 where u guarantees a unique differentiable solution for all nite T 4 Numerical Stability Numerical stability addresses the question of how numerical errors accumulate There are different types of numerical stability We will only look at one type in any detail known as absolute stability However here are the main types of stability as a summary 41 Types of Numerical Stability for ODEs 411 Zero Stability Zero stability addresses whether convergence is achieved Therefore it is only concerned with the stability of the method in the limit as At a 0 In fact we can say that for all Linear Multistep Methods 33 Consistency 1 Zero Stability ltgt Convergence 50 Furthermore the global error will typically be of the same order as the truncation error We always hope to use a convergent method so we always expect zero stability For nonlinear PDEs this might even be a lofty goall Determining zero stability of a method will be explained when we look at a more general type of stability known as absolute stability 412 Absolute Stability Absolute stability is a more general type of stability theory used to investigate how error accumu lates even for nite At not necessarily in the limit At a O This is often much more practical and useful since in practice we always use some nite At in applications Absolute stability is more relevant particularly in the eld of scienti c computation where scientists simulate real world77 problems Unfortunately in practice often it is erroneously assumed that convergence is all that matters As we will see a zero stable convergent method can be unstable for all nite At You saw this in your homeworkiEuler s method is a convergent method but is unstable for the pen dulum problem We will look at absolute stability in detail in Section 5 Absolute stability theory is usually con ned to stable problems For unstable problems absolute stability of a numerical method will not tell you much about global error 2From ODE theory must be Lipschitz continuous for such a result 10 413 Relative Stability One might ask as evidenced in class what if I m dealing with an unstable problem Relative stability tells you how the relative error in the solution grows or decays even if the solution is unstable Additionally one might hope that even for a stable problem the error should decay faster than the true deviation y 7 37 which itself decays for a stable problem Relative stability can be a complicated subject If you are interested you might do a Google search on order stars and relative stability We obviously won t address this fascinating topic No fascinating wasn t meant as sarcasm 5 Absolute Stability As mentioned absolute stability tells how error decays for a nite timestep At It is usually only informative for stable or neutrally stable problems We consider the scalar test problem ut Au 51 Remember that we are thinking of this test problem as the governing equation for the difference between the true solution and a perturbed solution Since we are considering linear methods the numerical solution will treat this deviation similarly to the way that it treats 51 Conveniently then we will isolate ourselves to 51 The idea is best introduced by example To adopt standard notation in the literature let s denote At k For Forward Euler the numerical solution satis es the difference equation Un1 Un kfUn Un kAUn 52 Solving for Un1 gives Un1 1 kUn 53 lterating this implies that the numerical solution at any tn is Un1 1 kA Uo 54 Now we look at values of k such that the solution either grows or decays We allow to be complex Now if 11 kAl gt 1 the modulus of Un E C grows as 71 increases It decays if l1kl lt 1 and remains constant if 1 1 ml 1 We can graph these regions of the complex plane 2 Ak The stable region 11 kM lt 1 is separated from the unstable region 1 1 ml gt 1 by the unit circle 11 kM 1 centered at k 71 and intersecting the origin This unit circle is known as the stability region for Forward Euler For a given this means that k must be chosen such that k is in the stability region if we expect the solution to 51 to decay This is somewhat nonintuitiveifor very stable problems where 9 ltlt 0 we would need a very small step size k to put k in the stable region guaranteeing a decaying numerical solution Otherwise the numerical solution would grow in magnitude even for a true solution that is rapidly decaying to zero Forward Euler is therefore usually not a good choice for problems with eigenvalues large in magnitude with negative real parts Backward Euler has a very different stability region Looking at the test problem again we have Un1 Un kfUn1 Un kAUnH 55 Solving for Un1 gives U U 56 n1 7 n lterating this implies that the numerical solution at any tn is Un1 U0 57 This means that 03le 1 separates the stable from the unstable region The stable region occurs where 1 7 ml gt 1 The boundary is the unit circle in the complex plane 2 k centered at k 1 The stable region is therefore the exterior of this unit circle Therefore the entire left half plane and more is the stable region for Backward Euler Therefore for a stable problem we can take whatever step size we wish with Backward Euler and the numerical solution will still decay to zero Therefore our step size is not limited by stability but only by accuracy concerns due to the truncation error This is true more generally of implicit methodsithey tend to have stability regions that include the entire left half plane That is the justi cation for the added computational cost of needing to do Newton iteration Interestingly since Backward Euler has a stability region the contains much of the right half plane even for unstable problems where 9 gt O a large step size with Backward Euler will produce a decaying numerical solution It is easy to show that the Trapezoid Method has a stability region that is exactly the left half plane with the boundary being the imaginary axis It is possible but somewhat more dif cult to show that the Midpoint method has a stability region that is only the interval 72 on the imaginary axis The explicit Midpoint Method is therefore the cheapest choice for problems that are neutrally stable with 9 O ie imaginary eigenvalues Of course k must be less in magnitude than 1ll for such problems so that k is in the stability region keeping the solution from growing in magnitude Note that the solution will not decay in such cases because we are always on the boundary of the stability region The same is true for the Trapezoid Method since the imaginary axis is the boundary Note also that for neutrally stable problems with imaginary eigenvalues Forward Euler will always produce a growing numerical solution with nite k and Backward Euler will always produce a decaying numerical solution with nite k since k is on the imaginary axisithe exterior and interior of the stability regions of those methods respectively However both methods will be convergent as k a 0 504 Relation to Zero Stability For Linear Multistep Methods they are zero stable if the origin of the complex plane is on the boundary of their absolute stability regions All of the methods discussed are therefore zero stable


Buy Material

Are you sure you want to buy this material for

25 Karma

Buy Material

BOOM! Enjoy Your Free Notes!

We've added these Notes to your profile, click here to view them now.


You're already Subscribed!

Looks like you've already subscribed to StudySoup, you won't need to purchase another subscription to get this material. To access this material simply click 'View Full Document'

Why people love StudySoup

Bentley McCaw University of Florida

"I was shooting for a perfect 4.0 GPA this semester. Having StudySoup as a study aid was critical to helping me achieve my goal...and I nailed it!"

Janice Dongeun University of Washington

"I used the money I made selling my notes & study guides to pay for spring break in Olympia, Washington...which was Sweet!"

Jim McGreen Ohio University

"Knowing I can count on the Elite Notetaker in my class allows me to focus on what the professor is saying instead of just scribbling notes the whole time and falling behind."


"Their 'Elite Notetakers' are making over $1,200/month in sales by creating high quality content that helps their classmates in a time of need."

Become an Elite Notetaker and start selling your notes online!

Refund Policy


All subscriptions to StudySoup are paid in full at the time of subscribing. To change your credit card information or to cancel your subscription, go to "Edit Settings". All credit card information will be available there. If you should decide to cancel your subscription, it will continue to be valid until the next payment period, as all payments for the current period were made in advance. For special circumstances, please email


StudySoup has more than 1 million course-specific study resources to help students study smarter. If you’re having trouble finding what you’re looking for, our customer support team can help you find what you need! Feel free to contact them here:

Recurring Subscriptions: If you have canceled your recurring subscription on the day of renewal and have not downloaded any documents, you may request a refund by submitting an email to

Satisfaction Guarantee: If you’re not satisfied with your subscription, you can contact us for further help. Contact must be made within 3 business days of your subscription purchase and your refund request will be subject for review.

Please Note: Refunds can never be provided more than 30 days after the initial purchase date regardless of your activity on the site.