NARRATOR: The following content is provided by MIT
OpenCourseWare under a Creative Commons license.
Additional information about our license and MIT
OpenCourseWare in general is available at ocw.mit.edu.
PROFESSOR: So, ready for the second lecture about -- second
to last lecture about the linear one way wave equation,
so the model for hyperbolic equation.
That word hyperbolic tells me that signals are traveling
with finite speed.
In this case, the speed is c and in this case, the signals
are only going to the last.
The two way wave equation will travel left and right, but
again, finite speed.
I emphasize finite speed, because for the heat equation,
the speed is infinite.
As soon as you change an initial condition, at any
later time, you immediately feel the effect --
probably a small effect way, way out.
So diffusion travels with infinite speed.
Waves travel with finite speed.
The different sort of is that the true --
the true growth factor g -- you remember, the true growth
factor was, for the wave equation.
That This is the g for the -- that multiplies the
exponential e to the ikx and the difference equation -- and
g, the true g, is e to the ikct.
So that --
what I was going to say is, that has absolute value 1.
Again, a difference from the heat equation.
For the heat equation with the fusion, g will be
smaller than 1.
You're losing energy.
For the wave equation, the energy is constant and the
wave is just moving.
So it's a little more delicate, because we're trying
to -- this g in the difference equation has to stay -- still
be somehow close to that one, which is right on the unit
circle and yet, we got to stay -- we can't go outside the
unit circle.
So here are two methods.
We had one method last time that worked.
That was one sided differences and of course, you had to take
the difference on the good side, the upwind side and that
was first order.
The Courant-Freidrichs-Lewy condition was the same as the
von Neumann stability condition and it was r
between 0 and 1.
r is this ratio, c delta t over delta x, it comes in.
It comes into the wave [? limit. ?]
So now we really have to think about better methods, other
methods and here are two: Lax-Freidrichs, which will
only be first order.
So when I -- these are the questions I always ask, and
then finally, experiment.
I mean, answering these questions will certainly be a
guide to what we'll see experimentally, but there's
more to learn from actually doing the numbers.
The accuracy, Lax-Freidrichs, will be only first order,
whereas Lax-Wendroff will be second order.
So with respect to accuracy, that's certainly better.
Can we see that maybe without too much calculation?
The notes would have the more details.
Can you see -- here's Lax-Wendroff, Lax-Freidrichs.
So what's the idea of Lax-Freidrichs?
You remember the total failure, the complete and
absolute failure when we used a time difference there and a
space difference there?
You remember that that gave the time difference that
produced a 1 ng and this space difference produced an
imaginary part ng, so we had g as 1 plus something imaginary
and we were outside the circle automatically.
Freidrichs thought of a way to get around that.
Instead of using ujn here, instead of the backward
difference, he replaced ujn by the average of that value and
that value, so that both sides now involved just two numbers
-- uj and uu, this one and this one.
This one is not used now.
Again, for a still first order, you won't be surprised
that this approximation, this has got first order error.
But stability now is back with us and easy to check because
now each new value depends on what old values -- just these
two, and with what coefficients?
So let me just check stability for Lax-Freidrichs.
So uj plus 1 n, the new value, is some multiple of this guy,
the one on the left and some multiple of this one, the one
on the right.
What are those multiples?
Let's see.
I'm going to move things over to the right hand side of the
equation, so I think there's a half and then I think there's
an r over --
minus, I guess.
This is the one that comes in with a minus sign --
it's the one on the left and the one on the
right, uj plus 1 n.
There's a half plus r over 2.
Simple idea.
Simple idea and it has one important advantage over
upwind, one sided that we saw last time.
It's not a great method.
It's only the first order here.
If I --
and you see when it's going to be stable?
The coefficients add to 1.
Our coefficients always add to 1, right?
Because if we have a constant state, this next step should
be that same state.
It's not going to change, but instead of having the typical
1 in the g, I have -- it's split into a half, so what is
g for this guy?
g for this guy is this number: 1 minus r over 2 times e to
the minus ik delta x plus this number times e to the
plus ik delta x.
You're seeing now the pattern.
So I've sort of jumped to g without plugging in the e to
the ikx , taking a step and then factoring
out e to the ikx.
I did it without writing down the e to the ikx.
This is what factors out.
Maybe a 1 plus r would be more correct.
So civility requires this complex number, which I could
write in different ways.
The real part there would be a cosine.
This would actually be the cosine plus r times the sine
-- plus ir times the sine.
This is really --
I could rewrite that as cos k delta x plus ar
sine k delta x.
I don't know which way you --
here, this shows you the real and imaginary parts.
That makes it nice.
This shows you the two coefficients.
That makes it nice.
Can you see stability -- so what is
Courant-Freidrichs-Lewy?
The characteristic discussion, the argument based on
characteristics, says what?
Stability is going mean that delta t --
stability is going to happen when r is below 1.
It's just like upwind, except that r could also be coming --
c could have the opposite sine.
Like Freidrichs is prepared for c to be positive or
negative, because it's not purely
upwind or purely downwind.
It's both winds.
You see that if r is anywhere --
the condition is going to be that minus 1 is -- but r is
between minus 1 and 1 --
see what's nice about Freidrichs in that case.
This number is positive.
This number is positive, because the 1 is
bigger than the r.
So each step --
and they add to 1.
So each step takes the new value as a combination of the
two old values and of course, the magnitude can't
be bigger than 1.
So under this condition, I have positive coefficients.
That's positive and that's positive.
They multiply things of absolute value 1, so they
can't -- the sum can't get beyond 1 -- or you might like
it here, that magnitude squared is cosined squared
plus r squared sine squared, and if r is below 1 where
smaller then cosine squared plus sine squared, which is 1.
So if I draw the --
just for a moment, if I draw the g in the complex plane --
so there's the unit circle. k equals 0 is always there.
That's at 1, then these points --
I think probably those points fall on an elipse here, where
we go up as high as r or magnitude of r and we go left,
right as far as --
this is k delta x equal pi as far as minus 1.
So it's stable if this is less than 1.
Nothing hard.
In fact, easy.
Two positive numbers that add to 1.
Now I've added a problem in the problem set.
Not requiring it to it be turned in, but just to
mention it here --
that whenever the coefficients are all positive, the accuracy
could only be first order, so we can't get too far with
positive coefficients.
First order accuracy is the best possible, except in the
extreme case when we're right on the characteristic.
So if r was exactly 1 --
you remember that special case when r is exactly 1 and we're
following exactly the line that carries the true
information --
r is 1.
That's gone and this is just --
so if r is 1, then uj plus 1 n is u -- is this, is that
number and we're just carrying that information exactly.
I have to say, again, that the idea of trying to follow that
characteristic is a totally natural one and I'm not
writing down here a message that tries to do that, but it
wouldn't be a bad idea in one dimension.
You couldn't do it, really, well in higher dimensions
where these ideas do continue to apply.
So that's Lex-Freidrichs.
That's Lex-Freidrichs.
So now let me turn to Lax-Wendroff.
These are natural ideas, but Lax and Wendroff did it first.
They said, if I use these three values, all three, then
certainly I can get the accuracy up one more because I
have one more coefficient to choose and to see what they
did, the easiest way is to write it as a time difference
on the left to see how they match.
All they want to do with match that second term, that delta x
squared term in the Taylor series.
Because this is going to give us a d by dt, this will give
us a d by dx centered, so we get extra accuracy there and
the matching that we need is this term.
Do you recognize that?
What does that term approximate?
Hope I've got it right.
That term approximates uxf, right?
So this is -- they've thrown in a term that approximates
uxx and its purpose is to catch the next order of
accuracy and hopefully not lose stability, hopefully
remain stable.
So what Lax-Wendroff are doing is staying even
closer to this curve.
I guess I have forgotten exactly what
their curve will be.
If I write down g of, maybe you could ask
MATLAB to do that.
Write down the g 4 -- this for Lax-Wendroff.
It's got that extra term and plot absolute value of g.
I guess you might plot it along with the
Lax-Freidrichs out.
So what am I guessing?
I'm pretty sure that the Lax-Wendroff one will hug the
unit circle more closely, but it'll stay inside.
So can I just --
Lax-Wendroff is g closer to 1, at least for a small k, at
least for a small k.
That's really where the accuracy is checked.
How close does these Taylor series match?
Lax-Wendroff match better and they stay inside.
That takes a little calculation.
Took me awhile yesterday to rewriter the --
write out the g, square the real part, square the
imaginary part, add them up and this is what I got.
So for Lax-Wendroff, I got g squared is 1 minus r squared
minus r to the 4th times 1 minus
cosine k delta x squared.
I put in the notes eventually, much algebra,
but it does the job.
Because if r squared --
so what's the stability condition?
It's again this same guy, the same condition -- r
squared below 1.
Again, we could be --
we could have a left going wave and a right going wave.
Now you'll say, what's the good of that?
Here we just have one wave, but on Friday, we'll have a
wave going both ways.
So Lax-Freidrichs and Lax-Wendroff can
both deal with that.
In fact, their stability conditions are the same and in
fact, you may be wondering why would anybody use
Lax-Freidrichs when Lax-Wendroff gives higher
accuracy at no real significant increase in cost?
So Lax-Wendroff is a favorite.
Don't forget, still experiment has to be done to see what it
really does in practice, because right now,
everything's looking in Lax-Wendroff's
favor with one exception.
Lax-Freidrichs, which only used these two coefficients
and had them positive.
This is positive.
This is positive, where Lax-Wendroff is going to have
a negative.
Let's see.
Where will our negative --
I think there'll be a negative coefficient in here somewhere.
Probably one of the three outside guys.
I know that we couldn't have three positive coefficients
and get up to the second order actors.
Now where do we pay the price for that?
Where do we pay the price --
I should finish here.
Suppose r squared's bigger than 1.
If you accept this calculation, which has a nice
result, but took awhile --
if r squared's below 1, then what do I know about g?
Suppose r squared's bigger than 1.
What can we say from this formula?
Does that hit you right away?
g --
if r squared --
this is certainly --
absolute g.
I've taken the actual g and squared its real parts,
squared its imaginary part, did all the algebra I could
think of and got this answer.
For the magnitude squared -- so the magnitude squared is
never going to go below 0.
I don't have to think of that.
It's just I want to be sure it's below 1.
So all I really want to do is be sure that I'm subtracting
something positive.
Am I subtracting something positive?
That's certainly positive.
Is this positive?
r squared'll be bigger than r to the fourth.
If r squared's less than 1, then r to the
fourth is even smaller.
So this is positive, this is positive.
I'm subtracting something so g is certainly
not bigger than 1.
So that formula really makes it evident and I guess also
evident that we're getting some high accuracy.
You see, what's the size of this when k is small?
When k is small, if I'm taking the cosine of a small angle --
what am I looking at there?
If this is a small angle theta, I'm looking at 1 minus
cosign theta -- that's of what order?
Theta squared.
Then that's squared again, so that I'm up to,
theta to the fourth.
That's why I'm getting off to a good start.
So stability we've got, two sided we've got, experiment --
can I just try to draw a little graph of what I think
the two methods will do for our model problem?
Like this model problem of a wave moving along?
I want to try to draw it at, let's say, t
equal to 1, for example.
So at t equal 1, the wave should be here.
I'll draw --
you might say this is totally unsatisfactory, and it is,
because I'm just drawing by hand what you can see.
Maybe you should.
This would be a perfect experiment and I really could
use for the notes --
the notes could use a proper graph.
I think what you see --
let's see.
At t equal 1, both methods are going to be 0 way out here and
they'll be 1 way out here.
Why's that?
Because the signal speed is finite, right?
This way out here, it hasn't heard about
this wave of water.
Way out here, it doesn't know that it's not a constant.
The question is, and I guess the thing has
moved, so maybe --
can I try roughly?
So the true wall of water would be that.
I think that of what Lax-Wendroff --
Lax-Freidrichs --
so Lax-Freidrichs with positive coefficients will
show a smooth profile, but the perfect signal
will get spread out.
The perfect signal will get spread out.
Why is that?
Let me just draw a spread out signal.
So apologies to do this by hand.
I think it's probably going to be pretty standard, but it's
going to be spread out.
So that, I think, would be Lax-Freidrichs.
That's very rough, right?
But I really --
hearty thanks for anybody who gives me a code that does it
and graphs it.
But the main point is, we have a lot of smearing here.
The accuracy is not going to be terrific.
But there isn't any bumpiness.
The reason there's no bumpiness is that we're
Lax-Freidrichs with --
has these positive coefficients, so that every
value always has to stay between 0 and 1.
The starting values were between 0 and 1.
At every step, the values are between 0 and 1, because we're
always taking a combination of a positive multiple of one
value and another positive multiple of the other value
and the two coefficients adding to one.
We're taking 30% of one and 70% of another and we never --
we're always in between the 0 and the 1.
What about Lax-Wendroff?
What are you expecting for Lax-Wendroff?
Expecting better accuracy, expecting --
if I were on the ball, I'd have another color.
Lax-Wendroff will be a much better profile, sharper
profile, but you know what the price is going to be?
There'll be some oscillation.
There'll be some oscillations behind --
I think behind -- oh gosh, I should know.
Forgive me for --
I can fix these if you'll tell me how on the next --
because I'm not sure that I should know whether they're at
both or whether the oscillations
happens at both ends.
I doubt it.
No, I doubt it.
I don't think it will.
Why do I not think it well?
Because in Lax-Wendroff, the three coefficients --
one of them will be negative and two will be positive.
They'll add to 1, of course.
So the oscillations, which are coming from sometimes a
negative coefficient, sometimes a positive, will all
be on one side.
So only one of those oscillations [UNINTELLIGIBLE]
So this is my picture of --
can I --
of Lax-Wendroff with oscillations
behind the wave front.
Alright.
If you've done any computing, you, when you see these
oscillations, you immediately think of somebody's name,
which is Gibbs.
Of course, Gibbs -- you remember the Gibbs phenomenon
that comes when you're approximating a square wave,
like there by pure Fourier series?
So the Gibbs, the classical Gibbs phenomenon is, take the
4 a series for a square wave while it's smooth and then
truncate that series, cut that series off at 100 terms or
1,000 terms.
Gibbs noticed -- and it's always a problem --
that these oscillations appear right near the front.
So far away, no problem, but here we really see it only on
one side, because we're not doing Gibbs Fourier series
experiment.
We're doing Lax-Wendroff wave equation.
When I say we, I really hope you'll do it and new lines of
MATLAB will be greatly appreciated.
Maybe Mr. Cho can --
you could email to him, for example, a code and the plot.
So what does that tell us?
Tells us Lax-Wendroff has got a lot of potential, but this
oscillation isn't good.
How are you going to get rid of it?
How do you get rid of bumpiness?
You introduce diffusion.
So diffusion, heat equation stuff smooth things out.
Wave equations don't smooth things out.
This wall of water in the wave equation stayed a sharp wall.
But we'll see very soon in the heat equation, diffusion,
viscosity, all those physical things connect numerically to
a smoother solution, smoother profile.
So that is really a key idea, is --
and we'll see it, but it will make the problem not linear.
If we want to stay with a linear equation and increase
the accuracy, we're going to get oscillation.
But if we locate the oscillation, introduce some
diffusion in the right place, we can reduce the oscillation.
So that's where --
this is what you really have to do if you want good output
from numerical methods.
Even on this simple model problem, you'll have to have
this idea of numerical viscosity.
Artificial viscosity, numerical diffusion popped in
at the right spot.
So we'll see that.
Maybe the right place to see it will be for when the
equation itself is already nonlinear.
THat's where it has developed a lot.
So we'll -- when we get to nonlinear wave equations,
which isn't far off, that's what we'll see.
So now I guess I have one more topic in this lecture.
that is convergence.
How do you show that --
we know, we sort of know that we need stability.
So I guess I'm going to have to come to the question of
convergence.
Alright.
Of course, I have to mention --
more than mention --
the key result in the whole theory, the connection --
so we have methods that are first order or second order --
stable or unstable, and then the key idea is this Lax --
again, he enters --
Lax equivalency.
So this is a little bit of pure math, comes into this and
it gives exactly the result we would expect.
Stability implies convergence and
convergence implies stability.
So that's the two ideas or equivalent, provided we're
talking about a decent problem, for a consistent
approximation to well-posed problems.
I've got four words on the board, four
keywords words here.
Stability, we kind of know what that means.
Convergence.
These are the new words and I'll just --
and actually -- so consistent means that the method is at
least first order and that the difference equation approaches
the differential equation at a single step.
Well-posed means that the differential equations is OK.
So, for example, what would be a non well-posed problem?
When we get to the heat equation, if you ran the heat
equation backward in time, that would give us an equation
dudt as minus uxx.
That would be not well-posed.
If you tried to solve Laplace's equation equation --
say, utt plus uxx equals 0 from initial values.
That would not be well-posed.
But we're dealing with well-posed problems.
The wave equation, which has the right side and the heat
equation with the right side.
So our problems are going to be well-posed.
You can safely think we're only talking about consistent
approximations to well-posed problems.
Here is the main point.
Stability holds if and only if the convergence holds.
So there's actually two things to show there.
How does stability give convergence?
Why does convergence requires stability?
So that instability could not [UNINTELLIGIBLE]
So that's sort of maybe the fundamental theorem of
numerical analysis.
The fact that stability, a balance --
we better say what stability is and what is convergence.
So let me capture these x in any quality, in mathematical
terms rather than just words.
So those are the ideas -- maybe I'm going to have to
write a little above here, but --
so what does --
what's our methods?
So I need some notation and then I'll lift the board and
we'll see what each thing means.
So our method, we're allowing any difference method that
moves forward in time.
There's time.
Let me use r for the true operator that takes the
function, follows the differential equation
for the next time.
So I'm going to say that true u, which is always a small u
at delta t is sum r, something r times u of t, where the
approximate guy is sum s times u of t.
So r is the true solution.
So well-posed.
So of course to take --
if this is n delta t and we've taken r n times, so that the
true solution is r to the nth times the initial value and
the approximate guy is sn times the initial value and
we'll -- we might as well start out with those.
That's what r and s are.
So now we got three ideas here.
What does well-posed --
well-posed means that the true thing is bound.
That r to the nth times u is less than some
constant times u.
These inequalities are sort of part of the mathematics of
partial differential equations.
I think you've got to take them as not --
I mean, take them as sort of a reasonable statement that the
-- this is energy this u typically will measure the
energy and in our wave equations, c 1, the constant c
1 will be 1.
The energy doesn't change.
In another problem, we could allow the possibility of c 1
being a little bigger than 1, but think of c 1 as 1.
Then stability will be abound on s to the n, say, c 2.
And again, in our examples, c 2 has been 1.
This is like -- c 2 is like --
I mean, s to the nth on each exponential is g to the nth
and g is below 1.
So often c 1 and c 2 will be 1.
Then finally, what is this ideas of consistency?
And this means --
consistent means that --
so what's consistency?
Consistent is comparing the differential equation with the
difference equation.
It's comparing r with s and at a single step,
they should be close.
So at a single step, r minus su --
that's the --
starting with the same u, I'm just taking one step in the
differential equation and comparing it with one step in
the difference equation and I would want that to be c 3,
maybe delta t to the p plus 1.
Was that -- is that --
times some size of u, I think.
Order delta t to the p plus 1.
Like delta t squared.
This would be the usual situation,
that when it's stable.
So at the moment then, I want to say, if these are all true,
then I get convergence because what I'm now going to do is
show that stability implies convergence.
It's what we did for ordinary differential equations.
This was the local error.
This is the local error that I called de, the
discreditization error.
These are saying -- so this is saying that the --
here's really what happens.
So I want to look at the error up here after n steps.
So the error after n steps --
typical term is, I follow the differential equation.
That's u.
Then I have this little error de, which I make at that step,
the difference between r and s and then I follow the
difference equation up to the remaining time.
What I want to show is that I can track the
error and bounce it.
I don't know if you remember the discussion for ordinary
differential equations, but it went the same way.
We had --
we looked at the new error de that was produced at
each time step k.
When we took steps of -- when we got up to this point with
the true solution, we plugged the true solution into the
difference equation, it had a little error and whatever
error there was, that got thrown in with other errors
and propogated up to any later time s.
Am I bringing back what we did there or just -- you can just
see that's what we're going to have here.
Now how will I express that?
Somehow I have to estimate the difference between r to the
nth and s to the nth.
That's really my problem.
Actually, it's a nice problem.
Suppose I have control of powers of r and s, and I know
that r minus s is small somehow.
How do I show that r to the nth minus s
to the nth is small?
That's my job.
This is the error when I applied to a starting value u.
This is the difference between little u and big u.
So how do I -- what I want to do is split it out into this
situation where I just have 1 r minus s.
I think this'll work.
So suppose I have --
I think I could start with r to the n minus 1 times r minus
s, and then in r to the n minus 2 times an r minus s
times an s.
Maybe I don't -- sorry, this guy is going to carry and I'd
rather do -- sorry.
Got to do it right because I want the rs to come first.
I want maybe r minus s times s to the n minus 1 and then --
sorry.
I want rs on the right.
Good.
Finally got it.
So that gives me the r of the nth, but I've subtracted it
off so I need an sr minus sr to the n minus 2 and so on.
Sorry for this bit of algebra in the middle of practical
things, but do you see that this is the case where I
followed the differential equation up and at the last
step I made an error?
This is the case where I followed the differential
equation almost to the top, I made an error and then I had
one more step to go with the difference.
The s doing this part, so a typical term will be a power
of s, s to some power, r minus s, r to some power.
I'm breaking the error at n steps into n different terms,
each from the error at one step.
This is the error at the last step.
This is the error at the next to last step and so on.
We've got it.
Because I just estimate every term, every term has a power
of r, which gives me a c 1, a power of s which gives me a c
2 and the discreditization error, the one step actual
error, r minus s, which gives me a c 3 and this delta t to
the p plus 1.
You see that I'm just piling up n of these, so I'm getting
n of those, delta t to the p plus 1, so the final result
will be, I'll have n of these delta t to the p plus 1s and
that will be t and delta t, delta t to the p.
So that'll be the error.
Since I went to all the trouble of defining a c 1, c
2, and c 3, those three numbers.
So the error --
there'll be a c 1 times c 2 times c 3 constant
there in the bound.
Have a look at the notes for that and let me just concludes
with Lax's insight on that other part.
The fact that convergence require stability.
So you see, what I've shown is that if we're stable, if I
bound these, bound these, bound these, then I have a
bound on the error.
That's convergence.
The other way is, if it's unstable, how do you know it
couldn't converge?
How do we know that
convergence requires stability?
It's a little delicate, so I just -- it's the last 30
seconds of the lecture, just to mention this point, because
it's easy to miss a key point here.
We've tracked a different --
now I'm supposing that the method is unstable and I want
to show that it will fail, right?
That's my job now.
That's Lax's -- other half of Lax's theorem: If it's
unstable, the method will fail.
There'll be some initial value for which it blows up.
Now you might say, that initial value is
just e to the iks.
Buw we didn't follow this point, but if I take it e to
the iks where g is bigger than 1 --
of course, that is getting amplified a bit, but for that
particular e to the ikx, when delta t decreases --
so if I have a --
you see, this is a k delta x.
What I'm trying to say is that for a fixed k, even an
unstable method can work.
For a fixed frequency k, as delta t and delta x
gets smaller --
these guys, whether they're outside or inside, will be
coming down toward 1 and if you're patient, you could see
it actually worked.
So it's somehow the initial function that doesn't work,
that blows up, is a combination of frequencies.
A single frequency, actually, you move down here and it's
not catastrophic in the end, but a combination of
frequencies, so that as one frequency moves down, another
one is moving into this position, another
one into this --
there is a typical initial function and it is disaster.
So that's what Lax had to catch.
The notes will use the single word, uniform boundedness.
Anyway, this is one point where pure math --
a kind of neat result in pure math says that if stability
fails, then there is some single starting value on which
it will fail.
That starting value won't be a pure e to the ikx, because it,
as we said, actually in the end gets to be OK, but there
will be a combination of frequencies that'll fail and
so convergence will fail.
So that's the key theorem.
I probably won't use the word theorem again, but this is one
time when it's justified and when a little sort of abstract
functional analysis gives the conclusion that convergence
requires stability.
So you'll see that in the notes and next time we move on
to the true wave equation, second order.
Thanks.