Model EquationsIn general, when we develop a numerical scheme to approximate a PDE, weanticipate that our scheme will not allow us to solve the PDE exactly: rather,it will produce a solution that contains some error. In particular, the localtruncation error of a given method is merely a measure of how well the truesolution of the difference equation satisfies our numerical method. An interestingquestion to ask, then, is the following:Is there a PDE to which our numerical approximation Qniis actuallythe exact solution?This question may be difficult to answer, but we should believe that the followingis somewhat easier:Can we at least find an equation that is better satisfied by Qnithanthe original PDE we were attempting to solve?If we can find such an equation, we can often learn a great deal about thenumerical method used to generate it, since it is usually much easier to studythe solutions of PDEs than those of finite difference formulas.In fact, using a Taylor series expansion, we can find a PDE which satisfiesthe difference equation exactly, but it will have infinitely many terms. The ideais to truncate this series at some point, yielding a PDE that is simple enough tostudy while simultaneously giving a good indication of the behavior of Qni. Oneinteresting fact is that, if the method is accurate to order s, the new equation(which we call the model equation) is generally a modification of the originalPDE with new terms of order s.Example Consider the the first order upwind method for the one dimensionaladvection equation φt+ uφx= 0 in the case u > 0.Qn+1i− Qni∆t+ uQni− Qni−1∆x= 0Qn+1i= Qni− u∆t∆xQni− Qni−1(1)We can insert a function v(x, t) into the numerical method (much in the sameway we insert the true solution q(x, t) when determining the local truncationerror) in order to find a differential equation that is satisied by v. Note thatv is a function that agrees with Qniexactly at the grid points, and thus v(x, t)satisfies (1) exactly:v(x, t + ∆t) = v(x, t) −u∆t∆x[v(x, t) − v(x −∆x, t)] .Now, if we Taylor expand about (x, t) and simplify, we getvt+12∆t vtt+16(∆t)2vttt+ ···+uvx−12∆x vxx+16(∆x)2vxxx+ ···= 0.1This can be rewritten asvt+ uvx=12(u∆x vxx− ∆t vtt) −16u(∆x)2vxxx+ (∆t)2vttt + ··· .This resulting equation is precisely the PDE that v satisfies. If we assume∆t/∆x is fixed, then the terms on the right hand side are O(∆t), O(∆t2), etc.,so for small ∆t we can truncate the series to obtain a PDE that is well satisfiedby the Qni. In particular, if we drop all the terms on the right, we recover theoriginal advection equation. Since this is equivalent to dropping terms of O(∆t),we expect that Qnisatisfies this equation to O(∆t), which we know to be correctsince this upwind method is first order accurate.If, instead, we keep the O(∆t) terms, we get:vt+ uvx=12(u∆x vxx− ∆t vtt) . (2)This involves second derivatives in both x and t, but we can derive a slightlydifferent model equation with the same accuracy by differentiating (2) withrespect to t to obtainvtt= −uvxt+12(u∆x vxxt− ∆t vttt)and with respect to x to obtainvtx= −uvxx+12(u∆x vxxx− ∆t vttx) .Combining these equations (by reording the partials) gives usvtt= u2vxx+ O(∆t).Combining this equation with (2) givesvt+ uvx=12(u∆x vxx− u2∆t vxx) + O(∆t2).Since we have already dropped O(∆t2) terms, we may do so here to obtainvt+ uvx=12u∆x(1 − ν)vxxwhere ν = u∆t/∆x is the Courant number.We have now transformed our model equation into a more familiar advection-diffusion equation, and the grid function Qnican be viewed as giving a secondorder accurate approximation to the true solution of this equation. The fact thatthe model equation for the upwind method is an advection-diffusion equationexplains a great deal about how the numerical solution behaves. Solutions tothe advection-diffusion equation translate at the proper speed u, but becomesmeared out over time.2If we examine the diffusion coefficient in our equation, we note that it van-ishes in the special case u∆t = ∆x. In this case, the exact solution to theadvection equation is recovered by the upwind method. Also, we note that thediffusion coefficient is positive only if 0 < u∆t/∆x < 1. This is precisely thestability limit of the upwind method! If it is violated, the diffusion coefficientin the model equation is negative, giving an ill-posed backward heat equation.Dissipation v. DispersionAs we have already seen, dissipation is essentially a kind of energy loss. Addingdissipation to the advection equation essentially says that the change in φ overtime results mostly from the bulk motion of the fluid flow, but not entirely.Dissipation has the net effect of making a wave form decay over time (things getsmeared out) and can thus be useful in damping spurious oscillations. Formallywe say that a one-step scheme has dissipation of order 2r if there exists apositive constant c independent of ∆t and ∆x, such that|g(∆xξ)| ≤ 1 −c(sin12∆xξ)2r.If dissipation causes a wave to decay over time, dispersion is a phenomena thatleads to the gradual separation of a waveform into a trail of oscilations. Werecall (via a few judicious applications of the Fourier inversion formula) that wecan write the solution of the one-way wave equation asu(t, x) =1√2πZ∞−∞eiωxe−iωatˆu0(ω)dω.From this, we can conclude that the Fourier transform of the solution satisfiesˆu(t + ∆t, ω) = e−iωa∆tˆu(t, ω).Recall too that, when we consider a one-step finite difference scheme, we haveseen thatˆvn+1= g(∆xξ)ˆvn.By comparing these two equations, we see that we can expect that g(∆xξ) willbe a good approximation to e−iωa∆x. In particular, we can write:g(∆xξ) = |g(∆xξ)|e−iξα(∆xξ)∆t.The quantity α(∆xξ) is called the phase speed and is the speed at which wavesof frequency ξ are propogated by the finite difference scheme. If α(∆xξ) wereequal to a for all ξ, then waves would propogate with the correct speed. How-ever, in practice, this is almost never the case. This is the precise definition ofdispersion: the finite different scheme propogates waves of different frequencieswith different speeds. To precisely quantify the error generated by dispersion,it is oftentimes useful to examine the phase error of a finite difference scheme,which is given by a − α(∆xξ).3Kinematic Description of Rigid Body