WEBVTT
00:00:00.000 --> 00:00:01.960
The following
content is provided
00:00:01.960 --> 00:00:06.310
by MIT OpenCourseWare under
a Creative Commons license.
00:00:06.310 --> 00:00:08.690
Additional information
about our license
00:00:08.690 --> 00:00:10.490
and MIT OpenCourseWare
in general
00:00:10.490 --> 00:00:11.800
is available at ocw.mit.edu.
00:00:15.122 --> 00:00:16.580
PROFESSOR: All
right, so I'll start
00:00:16.580 --> 00:00:23.270
with this review board, which
shows the conservation law.
00:00:23.270 --> 00:00:27.140
Scalar, one unknown,
one dimension,
00:00:27.140 --> 00:00:34.310
x, one space variable, and
the flux function is f of u.
00:00:34.310 --> 00:00:37.100
That's the differential form.
00:00:37.100 --> 00:00:42.900
We have seen that at some
points that differential form
00:00:42.900 --> 00:00:48.400
produces a discontinuous
solution, a discontinuity.
00:00:48.400 --> 00:00:51.650
In that case we
are better off with
00:00:51.650 --> 00:00:54.720
the integrated integral form.
00:00:54.720 --> 00:01:00.610
So I just integrate that one
from a point x_left to a point
00:01:00.610 --> 00:01:04.900
x_right, from an a to
b, and that gives me-- I
00:01:04.900 --> 00:01:07.620
have then the integral.
00:01:07.620 --> 00:01:11.070
So if u was density,
for example,
00:01:11.070 --> 00:01:13.630
the integral of density
will be the mass.
00:01:13.630 --> 00:01:17.530
That's the quantity
that I'm seeing here,
00:01:17.530 --> 00:01:20.210
so this is the rate
of change of mass
00:01:20.210 --> 00:01:29.310
in the interval comes from flux
going out minus flux coming in.
00:01:29.310 --> 00:01:34.530
So that's the conservation
law in its pure statement.
00:01:34.530 --> 00:01:40.880
And we learned that
solution for this--
00:01:40.880 --> 00:01:43.570
while we're looking at
the differential problem,
00:01:43.570 --> 00:01:46.170
we can follow these
characteristics
00:01:46.170 --> 00:01:49.490
from the initial time zero.
00:01:49.490 --> 00:01:55.960
So at time zero we look
at a point x and the value
00:01:55.960 --> 00:02:01.090
at that point just
travels-- is sort of safely
00:02:01.090 --> 00:02:07.510
carried up the characteristic
line, which is that equals x_0.
00:02:07.510 --> 00:02:14.970
And the problem arose when
characteristics either collided
00:02:14.970 --> 00:02:17.690
into each other or separated.
00:02:21.340 --> 00:02:23.480
Those are the two
difficulties and they
00:02:23.480 --> 00:02:28.480
have different resolutions.
00:02:28.480 --> 00:02:32.370
So this is what we met
in the traffic example
00:02:32.370 --> 00:02:35.580
when a light goes red.
00:02:35.580 --> 00:02:41.160
When a light goes red cars
are coming up to the light,
00:02:41.160 --> 00:02:46.950
in front of the light
the density is low.
00:02:46.950 --> 00:02:49.410
The cars have to brake quickly.
00:02:49.410 --> 00:02:52.600
So something happens--
that's a question mark there.
00:02:52.600 --> 00:02:55.060
What do you do when
characteristics collide?
00:02:55.060 --> 00:02:59.430
When the value here cannot
be both the initial value
00:02:59.430 --> 00:03:02.380
that's coming up that
characteristic and the initial
00:03:02.380 --> 00:03:05.720
value that's coming up there
because they're different?
00:03:05.720 --> 00:03:08.180
So we need a rule
for what to do there.
00:03:08.180 --> 00:03:10.780
And we also need a
rule for what to do
00:03:10.780 --> 00:03:13.600
in case of an empty space.
00:03:13.600 --> 00:03:16.920
When no characteristics--
here, two characteristics
00:03:16.920 --> 00:03:20.160
are hitting, here none
are touching there.
00:03:20.160 --> 00:03:23.280
Characteristics are separating
in this case and the question
00:03:23.280 --> 00:03:25.540
is, what solution do
you fill in the middle
00:03:25.540 --> 00:03:27.580
when you haven't got
a characteristic?
00:03:27.580 --> 00:03:32.030
It's not coming from the start.
00:03:32.030 --> 00:03:35.590
OK, what solution to choose?
00:03:35.590 --> 00:03:38.210
The answers will be different.
00:03:38.210 --> 00:03:44.490
And let me take the first one,
the collision question, first.
00:03:44.490 --> 00:03:52.850
So what will happen
is a shock will form.
00:03:52.850 --> 00:03:59.330
There'll be a shock line so
that these characteristics carry
00:03:59.330 --> 00:04:04.500
information into the shock, by
the way, that's a big deal--
00:04:04.500 --> 00:04:07.710
up until they hit the shock,
and these characteristics
00:04:07.710 --> 00:04:08.970
carry their information.
00:04:08.970 --> 00:04:14.110
And across the shock we
have different values,
00:04:14.110 --> 00:04:16.960
u on the left and
u on the right.
00:04:16.960 --> 00:04:23.490
So that's one possible
form of the solution
00:04:23.490 --> 00:04:26.660
because that will solve
the equation in this region
00:04:26.660 --> 00:04:28.760
where the characteristics
from the right
00:04:28.760 --> 00:04:30.920
are following the formula.
00:04:30.920 --> 00:04:35.060
It'll solve it from the left,
so everything comes down
00:04:35.060 --> 00:04:39.650
to the jump condition
at the shock.
00:04:39.650 --> 00:04:41.260
What's the shock speed?
00:04:41.260 --> 00:04:44.230
Where does this shock go?
00:04:44.230 --> 00:04:47.500
It's to find that
shock in the xt-plane
00:04:47.500 --> 00:04:52.970
that I have to use
this integral form.
00:04:52.970 --> 00:04:54.730
I must use the integral
form because I've
00:04:54.730 --> 00:04:57.880
got a discontinuity and
the differential form just
00:04:57.880 --> 00:04:59.280
doesn't make sense.
00:04:59.280 --> 00:05:04.790
OK, I wrote down the key idea
here or the conclusion, but not
00:05:04.790 --> 00:05:06.360
the reasoning.
00:05:06.360 --> 00:05:08.460
I'm looking for
this jump condition
00:05:08.460 --> 00:05:13.370
at the shock, which is
going to locate the shock.
00:05:13.370 --> 00:05:16.230
And the key quantity
that I have to find
00:05:16.230 --> 00:05:19.040
is the shock speed,
which will tell me
00:05:19.040 --> 00:05:24.520
where that shock is going,
how it's moving up in time.
00:05:24.520 --> 00:05:30.010
So time is going that way, x is
going this way in this picture.
00:05:30.010 --> 00:05:33.400
So this has a speed, s of t.
00:05:36.790 --> 00:05:38.330
And here's the conclusion.
00:05:38.330 --> 00:05:41.270
There's the jump condition.
00:05:41.270 --> 00:05:43.730
So what does it mean?
00:05:43.730 --> 00:05:48.190
This notation, this bracket
notation is the jump.
00:05:51.640 --> 00:05:53.350
The shock equation is beautiful.
00:05:53.350 --> 00:05:57.230
It says that s, the shock
speed, times the jump
00:05:57.230 --> 00:06:00.510
in u is the jump in f of u.
00:06:00.510 --> 00:06:03.640
And let me take
this example that
00:06:03.640 --> 00:06:07.950
is sort of the Burgers'
equation example.
00:06:07.950 --> 00:06:11.360
So this is the example that
comes from Burgers' equation.
00:06:18.260 --> 00:06:20.730
That's the example
where the flux
00:06:20.730 --> 00:06:27.420
is this nice parabolic
function u squared over 2.
00:06:27.420 --> 00:06:32.520
Then the change in flux, the
jump in the flux function
00:06:32.520 --> 00:06:36.630
is u right squared minus
u left squared over 2
00:06:36.630 --> 00:06:41.330
and I divide by the jump in u,
which is u_right minus u_left
00:06:41.330 --> 00:06:44.600
and that cancels
beautifully and leaves me
00:06:44.600 --> 00:06:48.680
with the 1/2 of
u_right plus u_left.
00:06:48.680 --> 00:06:54.270
It tells me that the shock speed
is exactly halfway between.
00:06:54.270 --> 00:06:58.820
It's exactly the average,
in this Burger's example,
00:06:58.820 --> 00:07:02.420
of the characteristic
speed from the right
00:07:02.420 --> 00:07:04.720
and the characteristics
speed from the left.
00:07:04.720 --> 00:07:07.770
And notice, it's in between.
00:07:07.770 --> 00:07:09.320
That's crucial.
00:07:09.320 --> 00:07:12.130
That tells us that
this figure is good.
00:07:12.130 --> 00:07:17.010
It tells us that characteristics
are entering the shock
00:07:17.010 --> 00:07:18.440
from both sides.
00:07:18.440 --> 00:07:20.550
Characteristic lines
are entering the shock.
00:07:20.550 --> 00:07:25.460
The shock is somehow, in
this sort of nonlinear case,
00:07:25.460 --> 00:07:31.070
it's somehow kind of cutting
out bits of the solution,
00:07:31.070 --> 00:07:33.590
reducing the energy,
reducing the entropy,
00:07:33.590 --> 00:07:37.830
reducing the total variation
as we'll see it in a minute.
00:07:37.830 --> 00:07:41.740
So that's the answer here and
I guess I have to justify,
00:07:41.740 --> 00:07:44.880
where did this come from?
00:07:44.880 --> 00:07:49.600
Well, it came from
this integral form.
00:07:52.210 --> 00:07:55.750
Suppose I integrate
from a little bit
00:07:55.750 --> 00:07:58.117
to the left of the
shock to a little bit
00:07:58.117 --> 00:07:59.200
to the right of the shock.
00:07:59.200 --> 00:08:02.630
Of course, that's where I
want to integrate in order
00:08:02.630 --> 00:08:04.540
to capture the shock.
00:08:04.540 --> 00:08:07.700
OK, so if I just
take that equation
00:08:07.700 --> 00:08:12.510
and write it over here, and
I'll imagine-- because I'm not
00:08:12.510 --> 00:08:15.410
going to be-- x_left and x_right
are going to be pretty near
00:08:15.410 --> 00:08:19.590
the shock-- u up to the
shock will pretty much
00:08:19.590 --> 00:08:21.630
have a constant value u_left.
00:08:21.630 --> 00:08:26.960
So approximately, that
is the shock position.
00:08:26.960 --> 00:08:30.600
Let me say the shock
is at-- capital
00:08:30.600 --> 00:08:32.100
X is the position of the shock.
00:08:32.100 --> 00:08:42.460
So I have x minus x_left times
u on the left of the shock.
00:08:45.790 --> 00:08:47.570
I'm doing that integral.
00:08:47.570 --> 00:08:51.450
Oh, I need d by dt of all this.
00:08:51.450 --> 00:08:52.200
So d by dt.
00:08:55.000 --> 00:08:57.540
It's d by dt of the
integral of u dx.
00:08:57.540 --> 00:09:01.520
Let me just put here
what I'm computing.
00:09:01.520 --> 00:09:05.840
I'm computing that
from x_left to x_right.
00:09:05.840 --> 00:09:08.560
Well, I'm going up
to x, up to capital
00:09:08.560 --> 00:09:10.880
X where the shock
is, and then onward.
00:09:10.880 --> 00:09:15.260
So here I went up to x
using the value u left.
00:09:15.260 --> 00:09:18.750
And now will be the
onward part, which
00:09:18.750 --> 00:09:25.050
will be x_right minus x, the
half of the interval that's
00:09:25.050 --> 00:09:27.250
on the right of the
shock, times u right.
00:09:30.160 --> 00:09:38.430
And then copying this term is
exactly our jump in f of u.
00:09:38.430 --> 00:09:46.530
Plus the jump in f of u is zero.
00:09:49.610 --> 00:09:53.020
Those words, total variation
are a reminder to myself
00:09:53.020 --> 00:09:54.270
for what comes next.
00:09:54.270 --> 00:09:58.330
OK, so is that OK?
00:09:58.330 --> 00:10:03.440
That will be approximately
zero because u wasn't exactly
00:10:03.440 --> 00:10:06.700
constant on the left and right.
00:10:06.700 --> 00:10:09.340
It didn't stay at the
value u_L but we're
00:10:09.340 --> 00:10:11.680
squeezing the
interval down where
00:10:11.680 --> 00:10:13.560
it's essentially constant.
00:10:13.560 --> 00:10:20.850
OK, so now, just look
and see what we've got.
00:10:20.850 --> 00:10:30.990
dx/dt is the shock speed, so
that's s, multiplying u_L,
00:10:30.990 --> 00:10:34.680
and here I have-- the
derivative of that is zero,
00:10:34.680 --> 00:10:37.230
the derivative of x_R is zero.
00:10:37.230 --> 00:10:39.750
Here I have a minus,
so it's minus u_R.
00:10:42.540 --> 00:10:50.320
dx/dt times that number plus the
jump in f of u, left to right,
00:10:50.320 --> 00:10:52.250
is 0.
00:10:52.250 --> 00:10:54.560
And that's our equation.
00:10:57.570 --> 00:11:00.790
If I put this on
the opposite side
00:11:00.790 --> 00:11:02.710
it becomes u_right minus u_left.
00:11:02.710 --> 00:11:04.330
Shall I do that?
00:11:04.330 --> 00:11:09.360
So I'll keep this on this
side, got a good sign,
00:11:09.360 --> 00:11:11.820
but I'll move this
onto the other side
00:11:11.820 --> 00:11:16.030
so it will be the shock
speed, d capital X dt times
00:11:16.030 --> 00:11:22.110
u_R minus u_L, which
is the jump in u.
00:11:22.110 --> 00:11:23.790
So that's great.
00:11:23.790 --> 00:11:28.420
We now have a way to deal
with shocks that collide.
00:11:31.740 --> 00:11:40.690
So that was a derivation of
this jump condition, which first
00:11:40.690 --> 00:11:45.490
had a name, two authors Rankine,
Hugoniot in gas dynamics,
00:11:45.490 --> 00:11:48.860
but then it was seen
as the right condition
00:11:48.860 --> 00:11:52.420
in all conservation laws.
00:11:52.420 --> 00:11:53.170
OK.
00:11:53.170 --> 00:11:56.630
And now I'm ready to
comment on the second issue.
00:11:56.630 --> 00:12:00.350
What do you do when
characteristics separate?
00:12:00.350 --> 00:12:03.220
Well, let me say
what you don't do.
00:12:03.220 --> 00:12:08.560
You could, without
wrecking the equation,
00:12:08.560 --> 00:12:12.440
you could try to
put a shock in there
00:12:12.440 --> 00:12:15.680
and give this the
value u on the left,
00:12:15.680 --> 00:12:18.250
even though a
characteristic isn't
00:12:18.250 --> 00:12:22.450
telling us what's going on
in this space or this one.
00:12:22.450 --> 00:12:24.810
You could try to put
a shock in there.
00:12:24.810 --> 00:12:26.900
You could maybe satisfy
the jump condition
00:12:26.900 --> 00:12:31.740
if you put the right shock in,
but the characteristics would
00:12:31.740 --> 00:12:36.150
be coming out and that's wrong.
00:12:36.150 --> 00:12:40.400
And it's called the
entropy condition, which
00:12:40.400 --> 00:12:42.740
says that this can't happen.
00:12:42.740 --> 00:12:46.390
That the derivative of
entropy has a certain sign
00:12:46.390 --> 00:12:48.860
and it rules this out.
00:12:48.860 --> 00:12:52.230
It rules out this possibility
of characteristics
00:12:52.230 --> 00:12:55.330
emerging from the shock and
we need a totally different
00:12:55.330 --> 00:12:56.170
solution.
00:12:56.170 --> 00:12:59.050
And it turns out to
be a fan in there.
00:12:59.050 --> 00:13:04.710
So we put in a fan
of characteristics
00:13:04.710 --> 00:13:09.500
coming from here and that
gives a solution that actually
00:13:09.500 --> 00:13:14.490
will depend on x over t.
00:13:14.490 --> 00:13:23.870
This solution will be a
simple-- a simple ratio x over t
00:13:23.870 --> 00:13:24.820
enters that.
00:13:28.500 --> 00:13:30.910
In other words, that gives
us a smooth solution,
00:13:30.910 --> 00:13:34.410
and I guess, thinking about
the traffic example, that's
00:13:34.410 --> 00:13:35.640
exactly what happens.
00:13:35.640 --> 00:13:37.650
So this is the shock part.
00:13:37.650 --> 00:13:42.850
Now for the fan
part, what happens
00:13:42.850 --> 00:13:48.740
when we have a bunch of cars
sitting there, ready to go?
00:13:48.740 --> 00:13:56.430
The light goes green, the
first cars start moving,
00:13:56.430 --> 00:14:02.150
cars behind them see the car in
front moving, they start moving
00:14:02.150 --> 00:14:11.140
and so the initial profile
would be, in this case,
00:14:11.140 --> 00:14:20.180
for cars, would be high
density and then zero density.
00:14:20.180 --> 00:14:28.800
So this would be maximum density
sitting at the light and zero
00:14:28.800 --> 00:14:33.180
density in front of the light
and then the light goes green
00:14:33.180 --> 00:14:35.180
now and what happens?
00:14:35.180 --> 00:14:38.920
Well, as we all know those
first cars pull ahead
00:14:38.920 --> 00:14:47.670
and the profile is these cars
haven't yet started to move.
00:14:47.670 --> 00:14:50.340
These have started to
move and there's a profile
00:14:50.340 --> 00:14:56.160
there and then zero density.
00:14:56.160 --> 00:15:00.290
So this represents the
first car, the lead car.
00:15:00.290 --> 00:15:08.900
In front of that car is zero
density and that's the fan.
00:15:08.900 --> 00:15:14.600
That's the x over t expression.
00:15:14.600 --> 00:15:19.880
So a simple linear expression,
which will flatten out;
00:15:19.880 --> 00:15:24.190
as time goes on, these cars
will be going faster and faster
00:15:24.190 --> 00:15:28.440
and then eventually
up to full speed.
00:15:28.440 --> 00:15:30.860
OK, so that's the fan.
00:15:34.950 --> 00:15:41.060
Maybe now is my chance to
speak about this concept
00:15:41.060 --> 00:15:44.350
of the variation in a function.
00:15:44.350 --> 00:15:49.790
The idea will be now, how do
we get the finite difference
00:15:49.790 --> 00:15:54.370
method to do this automatically?
00:15:54.370 --> 00:16:00.180
Well, OK, one way to do it
is put in a little viscosity.
00:16:03.540 --> 00:16:10.180
This rule, the rule that picks
out this shock in one case,
00:16:10.180 --> 00:16:17.810
of colliding
characteristics, and this fan
00:16:17.810 --> 00:16:22.250
in the case of separating
characteristics, one way
00:16:22.250 --> 00:16:25.520
to get to those
two solution would
00:16:25.520 --> 00:16:32.560
be to put in a little
bit of viscosity
00:16:32.560 --> 00:16:35.930
and then let epsilon go to zero.
00:16:35.930 --> 00:16:42.070
And the result would
be these two solutions.
00:16:42.070 --> 00:16:45.950
And in a way, that's what finite
differences is going to do.
00:16:45.950 --> 00:16:47.760
Finite difference
or finite volume
00:16:47.760 --> 00:16:52.530
is going to be a
numerical viscosity
00:16:52.530 --> 00:16:58.400
and the epsilon in
the Burgers' equation
00:16:58.400 --> 00:17:05.844
would become delta t or a delta
x in the numerical method,
00:17:05.844 --> 00:17:07.260
in the finite
difference equation.
00:17:10.810 --> 00:17:12.600
And then what's the concept?
00:17:12.600 --> 00:17:15.740
So we have to find
finite differences.
00:17:15.740 --> 00:17:25.100
And we want those-- the
point is finite differences
00:17:25.100 --> 00:17:27.300
sort of takes more care.
00:17:27.300 --> 00:17:31.310
We can't just create
something that's
00:17:31.310 --> 00:17:35.620
consistent in the smooth case.
00:17:35.620 --> 00:17:40.210
Up till now we've just
had no hesitation.
00:17:40.210 --> 00:17:43.680
We replaced derivatives
by differences
00:17:43.680 --> 00:17:46.270
without too many
worries, really.
00:17:49.040 --> 00:17:51.500
And that was fine as long
as the solution is smooth.
00:17:51.500 --> 00:17:55.630
But now if we're going
to replace derivatives
00:17:55.630 --> 00:17:59.370
by differences, we have
to do it in such a way
00:17:59.370 --> 00:18:05.400
that this entropy
condition, these solutions
00:18:05.400 --> 00:18:07.800
emerge automatically.
00:18:07.800 --> 00:18:15.420
And really, in a way that--
our finite difference
00:18:15.420 --> 00:18:18.780
or our finite volume
methods should somehow
00:18:18.780 --> 00:18:24.100
conserve mass the same way
the true equation does.
00:18:24.100 --> 00:18:27.020
So those are two
different points.
00:18:27.020 --> 00:18:31.720
First point is this
question of total variation.
00:18:31.720 --> 00:18:33.610
What's the total
variation in a function?
00:18:36.240 --> 00:18:41.220
Well, let me draw a function.
00:18:41.220 --> 00:18:44.190
The total variation
is-- well, let me see.
00:18:44.190 --> 00:18:50.480
Total variation in a
function u is going to be,
00:18:50.480 --> 00:18:57.130
by definition, the integral of
the absolute value of du/dx.
00:19:03.180 --> 00:19:09.060
OK now we have to
understand this concept.
00:19:09.060 --> 00:19:21.350
It's beautiful that the analysis
has led to a quantity as simple
00:19:21.350 --> 00:19:25.750
as that, as a guide.
00:19:25.750 --> 00:19:31.250
So when we create finite
difference methods we will ask,
00:19:31.250 --> 00:19:36.260
are they total
variation diminishing?
00:19:36.260 --> 00:19:42.060
So TVD will be total
variation diminishing methods.
00:19:42.060 --> 00:19:46.600
Methods where this quantity
does not increase in time.
00:19:49.180 --> 00:19:52.050
Let me just understand
this quantity first.
00:19:52.050 --> 00:19:53.600
The absolute value is crucial.
00:19:53.600 --> 00:19:57.550
If the absolute
value was not there,
00:19:57.550 --> 00:20:01.770
then I would just have the
integral of the derivative,
00:20:01.770 --> 00:20:06.250
so it would be u at that
point minus u at that.
00:20:06.250 --> 00:20:11.110
But here I'm taking absolute
values, so for a while
00:20:11.110 --> 00:20:14.220
the derivative is
positive, up to there.
00:20:14.220 --> 00:20:16.770
So the integral gives me
this value minus this.
00:20:19.840 --> 00:20:24.160
Then, in this region, the
derivative is negative,
00:20:24.160 --> 00:20:26.560
but I'm taking absolute value.
00:20:26.560 --> 00:20:30.450
Do you see that that
integral will split
00:20:30.450 --> 00:20:32.180
into four separate integrals?
00:20:32.180 --> 00:20:37.170
The total variation will be this
difference plus this difference
00:20:37.170 --> 00:20:39.350
plus this plus this.
00:20:39.350 --> 00:20:42.470
I'll have four--
in this example--
00:20:42.470 --> 00:20:46.960
have four positive
contributions to the variation.
00:20:51.350 --> 00:20:54.100
And that's the key
quantity there.
00:20:54.100 --> 00:20:56.470
And what a shock
will do-- can I maybe
00:20:56.470 --> 00:21:03.630
mention in this figure what
a shock does, can do, would
00:21:03.630 --> 00:21:10.970
be to somehow cut out
part of the variation.
00:21:10.970 --> 00:21:14.940
When a shock occurs it might
happen that these things are
00:21:14.940 --> 00:21:17.410
going into the shock,
these are going
00:21:17.410 --> 00:21:23.310
in from that side, and
the shock kind of jumps
00:21:23.310 --> 00:21:28.710
from this value, u_left,
to this value, u_right,
00:21:28.710 --> 00:21:32.670
following the shock
condition when it emerges.
00:21:32.670 --> 00:21:37.030
And this additional
variation is wiped out.
00:21:37.030 --> 00:21:43.020
You see that the variation
for the solution that
00:21:43.020 --> 00:21:46.810
contains a shock is still--
this part is still there,
00:21:46.810 --> 00:21:50.150
but this part is all
downhill, so it's just this.
00:21:52.710 --> 00:21:54.460
It's that minus that.
00:21:54.460 --> 00:21:58.350
And the middle two
parts have been cut out.
00:21:58.350 --> 00:22:01.840
OK, so that's the total
variation in a function.
00:22:01.840 --> 00:22:06.560
Let me just take the
most important example
00:22:06.560 --> 00:22:09.274
that's not total variation
diminishing, which
00:22:09.274 --> 00:22:09.940
is Lax-Wendroff.
00:22:12.720 --> 00:22:16.250
The one-- even in
the linear case.
00:22:21.330 --> 00:22:27.400
The one experiment we've done
with that profile, where,
00:22:27.400 --> 00:22:34.220
you remember, we had a wall of
water and the exact solution
00:22:34.220 --> 00:22:39.820
brought the wall to the left
with velocity c, with speed c,
00:22:39.820 --> 00:22:44.000
but Lax-Wendroff,
which did quite well
00:22:44.000 --> 00:22:53.410
staying near the
discontinuity, near the jump,
00:22:53.410 --> 00:23:00.970
had a little oscillation behind.
00:23:00.970 --> 00:23:05.590
That was essentially our one
objection to Lax-Wendroff,
00:23:05.590 --> 00:23:09.900
was that oscillation showing up.
00:23:09.900 --> 00:23:14.970
And so that's telling us
that the variation-- see,
00:23:14.970 --> 00:23:20.550
the variation in the true
solution is just that step.
00:23:20.550 --> 00:23:24.520
The variation in the
Lax-Wendroff answer
00:23:24.520 --> 00:23:28.540
is this bigger
step plus this one
00:23:28.540 --> 00:23:31.750
plus this one plus
this, larger than 1.
00:23:36.010 --> 00:23:45.000
Roughly the match is that
monotone methods-- so here's
00:23:45.000 --> 00:23:46.750
the general picture.
00:23:46.750 --> 00:23:57.730
We can get monotone, or I'll
say, all positive coefficients,
00:23:57.730 --> 00:24:02.330
and those will be total
variation diminishing.
00:24:02.330 --> 00:24:04.460
Great.
00:24:04.460 --> 00:24:07.990
Why don't I draw--
recall a monotone example
00:24:07.990 --> 00:24:12.650
was Lax-Friedrichs or upwind.
00:24:12.650 --> 00:24:17.320
Upwind will be our model, so
let me draw the upwind one.
00:24:17.320 --> 00:24:24.830
The upwind one didn't increase
the variation, it stayed at 1.
00:24:24.830 --> 00:24:31.630
We had no oscillations,
but it was low accuracy.
00:24:31.630 --> 00:24:32.860
So first order.
00:24:38.440 --> 00:24:45.355
We can't do better than
first-order accuracy, the order
00:24:45.355 --> 00:24:50.860
of accuracy I'm measuring in
the smooth part of the solution
00:24:50.860 --> 00:24:54.950
because that's where I can take
Taylor series and match terms
00:24:54.950 --> 00:24:58.730
and see which is the first
one that doesn't match.
00:24:58.730 --> 00:25:05.950
And Lax-Friedrichs and
upwind was first order.
00:25:05.950 --> 00:25:22.540
Then our problem is that if I
allow a negative coefficient
00:25:22.540 --> 00:25:32.080
then variation increases,
oscillations will come in,
00:25:32.080 --> 00:25:37.730
but I can get the accuracy
up to second order or higher.
00:25:41.480 --> 00:25:46.130
And Lax-Wendroff is the
example, the key example
00:25:46.130 --> 00:25:47.960
there that I'll continue with.
00:25:47.960 --> 00:25:53.780
OK, that's the idea
of total variation.
00:25:53.780 --> 00:25:59.650
It just gives us
a quantity to work
00:25:59.650 --> 00:26:07.040
with when we see oscillation.
00:26:07.040 --> 00:26:10.770
And of course, I'm
interested in oscillations,
00:26:10.770 --> 00:26:14.920
so this is a linear case
with a discontinuity.
00:26:14.920 --> 00:26:18.540
I don't think it's quite
fair to call it a shock.
00:26:18.540 --> 00:26:21.500
At least, the
shocks that we meet
00:26:21.500 --> 00:26:25.230
in conservation laws in
the nonlinear problems,
00:26:25.230 --> 00:26:28.130
they're created--
the ones we saw here
00:26:28.130 --> 00:26:34.970
are created for nonlinear
reasons, colliding
00:26:34.970 --> 00:26:38.960
characteristics, whereas the
linear case has characteristics
00:26:38.960 --> 00:26:39.830
all parallel.
00:26:42.930 --> 00:26:46.240
So here, we're in
the lecture, we're
00:26:46.240 --> 00:26:52.460
really into the second lecture
into the scientific computing
00:26:52.460 --> 00:26:54.010
part of our problem.
00:26:54.010 --> 00:26:58.930
So that's the analysis
part, the PDE part.
00:26:58.930 --> 00:27:08.380
Now comes the experience of
quite a lot of people working
00:27:08.380 --> 00:27:17.640
hard to identify a way to get
second-order accuracy where
00:27:17.640 --> 00:27:25.220
it's smooth and monotonicity
or something like it
00:27:25.220 --> 00:27:28.090
where the shocks are.
00:27:28.090 --> 00:27:32.640
So that's essentially the goal.
00:27:32.640 --> 00:27:35.850
Our method is going to
turn out to be nonlinear.
00:27:35.850 --> 00:27:42.510
It's going to have some--
it'll be called a flux limiter.
00:27:42.510 --> 00:27:47.680
The key quantity
will be the flux.
00:27:47.680 --> 00:27:56.940
And for a smooth function
when u_right and u_left are
00:27:56.940 --> 00:28:02.000
very close, the flux is small,
we don't have a problem.
00:28:02.000 --> 00:28:06.250
But at a shock where there's
a significant jump in u
00:28:06.250 --> 00:28:14.250
and especially in f
of u, the flux limiter
00:28:14.250 --> 00:28:20.160
is going to come in,
into the discrete method.
00:28:20.160 --> 00:28:23.400
Let me begin discrete
methods with a new word,
00:28:23.400 --> 00:28:24.520
finite volumes.
00:28:27.650 --> 00:28:34.320
And with a new
interpretation of capital U.
00:28:34.320 --> 00:28:40.950
So the discrete value,
capital U_(j, n), which,
00:28:40.950 --> 00:28:44.300
up to now we've been thinking
of as that's an approximation
00:28:44.300 --> 00:28:48.860
to the true solution,
at j delta x, n delta t.
00:28:48.860 --> 00:28:52.290
Now, because we're dealing
with discontinuities,
00:28:52.290 --> 00:28:54.060
we better integrate.
00:28:54.060 --> 00:28:58.230
We integrate over
a delta x interval.
00:28:58.230 --> 00:29:01.130
We take the average over
a little delta x interval
00:29:01.130 --> 00:29:04.120
instead of hitting
just the point.
00:29:04.120 --> 00:29:08.270
So the delta x interval from
j minus 1/2 to j plus 1/2,
00:29:08.270 --> 00:29:13.340
we take the average of
u at the time n delta t
00:29:13.340 --> 00:29:23.670
and that's what finite volume--
and let's think of as capital
00:29:23.670 --> 00:29:27.360
U, the discrete approximation.
00:29:27.360 --> 00:29:31.940
The appropriateness of that
is in the integral form
00:29:31.940 --> 00:29:34.260
of the conservation law.
00:29:34.260 --> 00:29:36.330
So the integral form
of a conservation law
00:29:36.330 --> 00:29:49.670
leads us to that natural
discretization, d by dt of U,
00:29:49.670 --> 00:29:52.060
of the average.
00:29:52.060 --> 00:29:55.750
You're keeping your eye
on the integral form
00:29:55.750 --> 00:29:56.960
of the conservation law.
00:29:56.960 --> 00:29:59.720
I'll just come back to it again.
00:29:59.720 --> 00:30:04.280
The derivative of the
mass, the integral
00:30:04.280 --> 00:30:10.430
of u-- so in the discrete case,
the difference of capital U--
00:30:10.430 --> 00:30:15.560
plus the integrated
term of flux at one end
00:30:15.560 --> 00:30:17.530
of the interval minus a
flux at the other end.
00:30:17.530 --> 00:30:22.410
And notice that we are
getting for the flux, half--
00:30:22.410 --> 00:30:32.320
so a staggered grid is
presenting itself naturally.
00:30:32.320 --> 00:30:35.780
And so that
f_(j+1/2), the fluxes,
00:30:35.780 --> 00:30:41.130
which is the critical thing,
is on the j plus 1/2 grid.
00:30:41.130 --> 00:30:45.650
All right, now this
could be exactly defined
00:30:45.650 --> 00:30:52.200
as the flux at that grid point.
00:30:52.200 --> 00:30:55.990
Capital U squared over 2
in the Burger equation,
00:30:55.990 --> 00:30:58.550
at that grid point.
00:30:58.550 --> 00:31:04.560
That would give us a totally
straightforward, first-order
00:31:04.560 --> 00:31:05.800
numerical method.
00:31:05.800 --> 00:31:08.710
First order, I see
first order here in time
00:31:08.710 --> 00:31:11.550
and it might be better
than that in space
00:31:11.550 --> 00:31:13.250
because in some way
that's centered.
00:31:15.970 --> 00:31:19.450
But if we want higher
accuracy, if we
00:31:19.450 --> 00:31:21.700
want to include
Lax-Wendroff, for example,
00:31:21.700 --> 00:31:26.530
that includes
several values of u,
00:31:26.530 --> 00:31:31.510
we go to a more
general flux function.
00:31:34.670 --> 00:31:39.610
The whole job here is create
a smart flux function.
00:31:39.610 --> 00:31:48.790
So a flux function depends not
only on u in the center, u_j,
00:31:48.790 --> 00:31:52.100
at time n-- this
is all at time n--
00:31:52.100 --> 00:31:56.350
it can also depend on
values to the left of it
00:31:56.350 --> 00:31:57.400
and values to the right.
00:32:00.350 --> 00:32:02.970
It may use those
values, it may not,
00:32:02.970 --> 00:32:05.870
but we allow that
possibility then,
00:32:05.870 --> 00:32:08.900
that the flux-- so this
is the numerical method
00:32:08.900 --> 00:32:10.440
that you aim for.
00:32:10.440 --> 00:32:16.210
Is to choose a flux
function and of course,
00:32:16.210 --> 00:32:19.620
this should be consistent
with the original problem
00:32:19.620 --> 00:32:26.180
and that simply says that if
all the u's were the same,
00:32:26.180 --> 00:32:32.159
that the numerical flux
function would give you the flux
00:32:32.159 --> 00:32:33.450
function in the given equation.
00:32:36.860 --> 00:32:39.750
OK, so this is our
general pattern.
00:32:39.750 --> 00:32:43.780
And now, the question is
really, maybe the first step
00:32:43.780 --> 00:33:00.900
is write a-- I would like to
write upwind and Lax-Wendroff--
00:33:00.900 --> 00:33:02.960
I'd like to find
their flux functions.
00:33:02.960 --> 00:33:06.710
I guess I'm saying
that when we look
00:33:06.710 --> 00:33:14.630
at an equation in
this format, we
00:33:14.630 --> 00:33:17.150
can take methods that
we have already used
00:33:17.150 --> 00:33:22.000
and see, OK, even in the
linear case, the one-way wave
00:33:22.000 --> 00:33:26.850
equation, what are we picking
for the flux function?
00:33:26.850 --> 00:33:32.170
OK, now the one distinction
is: the one-way wave equation,
00:33:32.170 --> 00:33:38.150
we had a constant, and
upwind, for example,
00:33:38.150 --> 00:33:44.070
went in the direction
decided by that speed, c.
00:33:44.070 --> 00:33:49.350
If the speed reversed sign
so that the waves were going
00:33:49.350 --> 00:33:52.390
the other way, the wind
is coming the other way,
00:33:52.390 --> 00:33:55.360
the upwind direction
is the other way,
00:33:55.360 --> 00:33:57.870
so our method would change.
00:34:02.770 --> 00:34:05.450
Let me just write
down the flux function
00:34:05.450 --> 00:34:12.210
that turns out to
come for upwind
00:34:12.210 --> 00:34:14.290
and then for Lax-Wendroff.
00:34:14.290 --> 00:34:18.120
OK, so for upwind the
flux turns out to be--
00:34:18.120 --> 00:34:24.730
F turns out to be 1/2
of-- because it's at this,
00:34:24.730 --> 00:34:26.720
this is f_(j+1/2), really.
00:34:33.460 --> 00:34:34.510
It's at the center point.
00:34:34.510 --> 00:34:44.000
The flux function will be--
I'll put in the a-- a delta t--
00:34:44.000 --> 00:34:52.620
sorry, a over 2 U_(j+1)
plus U_j, now, plus--
00:34:52.620 --> 00:35:05.290
or rather minus actually--
comes absolute value of a over 2
00:35:05.290 --> 00:35:08.710
times u_(j+1) minus u_j.
00:35:12.510 --> 00:35:19.660
If we're patient we can track
down the two cases, a positive
00:35:19.660 --> 00:35:28.770
and a negative, and recognize
that this form switches
00:35:28.770 --> 00:35:35.160
direction correctly on
the difference method.
00:35:35.160 --> 00:35:37.230
So it keeps us upwind,
in other words.
00:35:37.230 --> 00:35:39.750
If the wind changes
direction our flux
00:35:39.750 --> 00:35:44.390
changes appropriately to
follow that direction.
00:35:44.390 --> 00:35:48.720
So that's an OK flux function.
00:35:48.720 --> 00:35:57.480
It's a monotone one
and it's first order.
00:35:57.480 --> 00:36:00.190
Now let me put in Lax-Wendroff.
00:36:00.190 --> 00:36:02.420
What would Lax-Wendroff be?
00:36:02.420 --> 00:36:08.790
Would be the same, the same
upwind flux, if I call that UP
00:36:08.790 --> 00:36:12.310
for upwind and
then a correction,
00:36:12.310 --> 00:36:19.670
a Lax-Wendroff term that gets
the second-order accuracy.
00:36:19.670 --> 00:36:26.730
OK, and that turns out to be--
it turns out that I multiply
00:36:26.730 --> 00:36:32.780
this term by that ratio r.
00:36:32.780 --> 00:36:34.340
By the absolute value of r.
00:36:34.340 --> 00:36:39.700
So it turns out that it's--
the Lax-Wendroff thing,
00:36:39.700 --> 00:36:42.240
to get that extra accuracy,
you might remember,
00:36:42.240 --> 00:36:45.200
had an extra delta
t over delta x.
00:36:45.200 --> 00:36:51.150
And so it's this same
thing times a over 2.
00:36:51.150 --> 00:36:55.730
In the nonlinear case
it'll be f of u_(j+1).
00:36:59.050 --> 00:37:03.200
f at u_(j+1) minus f at u_j.
00:37:10.090 --> 00:37:16.640
Next time, I'll come
back to this point
00:37:16.640 --> 00:37:24.040
and pin down this
more carefully.
00:37:24.040 --> 00:37:26.450
All I want to do
in completing today
00:37:26.450 --> 00:37:35.210
is introduce this idea
of a flux limiter.
00:37:35.210 --> 00:37:42.920
So what people say when
they look at Lax-Wendroff
00:37:42.920 --> 00:37:49.010
is that this term
that's been added on
00:37:49.010 --> 00:37:52.280
is fine in the case
of a smooth solution,
00:37:52.280 --> 00:37:59.320
but in that case of a
jump, like our picture
00:37:59.320 --> 00:38:07.030
here, that term is
increasing variation.
00:38:07.030 --> 00:38:09.650
It's an anti-diffusion term.
00:38:09.650 --> 00:38:14.390
Instead of smoothing things
out the way F upwind does,
00:38:14.390 --> 00:38:21.830
that extra term has the
effect of increasing the jump.
00:38:21.830 --> 00:38:23.740
It's anti-diffusive.
00:38:23.740 --> 00:38:31.500
And therefore, when
it's a significant term,
00:38:31.500 --> 00:38:39.120
the flux limiter
cuts back on it.
00:38:39.120 --> 00:38:41.070
The result will be
a numerical method,
00:38:41.070 --> 00:38:43.580
which is not linear
anymore even if we apply it
00:38:43.580 --> 00:38:45.300
to a linear problem.
00:38:45.300 --> 00:38:49.190
The result won't
be strictly linear
00:38:49.190 --> 00:38:54.400
and it will produce
a profile that
00:38:54.400 --> 00:38:56.960
doesn't have that oscillation.
00:38:56.960 --> 00:39:00.630
So it has what we want.
00:39:00.630 --> 00:39:05.410
It gives us good answers
when the solution is smooth,
00:39:05.410 --> 00:39:08.820
it stays near the shock
the way Lax-Wendroff does,
00:39:08.820 --> 00:39:11.490
but it doesn't oscillate.
00:39:11.490 --> 00:39:15.550
OK, so what's the flux limiter?
00:39:15.550 --> 00:39:20.000
The flux limiter is, I multiply
in here by some function,
00:39:20.000 --> 00:39:26.590
the limiter function, phi_(j+1)
and I multiply in here
00:39:26.590 --> 00:39:28.230
by phi_j.
00:39:28.230 --> 00:39:33.050
In other words, there's
an extra factor,
00:39:33.050 --> 00:39:36.500
this flux limiting factor.
00:39:36.500 --> 00:39:40.090
So the question is,
how to construct that?
00:39:40.090 --> 00:39:43.700
What formula do we use
for the flux limiter?
00:39:43.700 --> 00:39:49.160
So the flux limiter,
after lot of experiments
00:39:49.160 --> 00:40:00.540
is chosen to depend on--
so we want it to be 1.
00:40:00.540 --> 00:40:02.980
So what do we want?
00:40:02.980 --> 00:40:10.380
We want it to be 1 or
close to 1 at smooth part
00:40:10.380 --> 00:40:17.480
and we want it to be greater
than 1 at a discontinuity.
00:40:20.080 --> 00:40:23.030
OK, so our question
is, how do you
00:40:23.030 --> 00:40:25.420
identify that discontinuity?
00:40:25.420 --> 00:40:29.850
Well, the convenient
way has turned out
00:40:29.850 --> 00:40:35.980
to be look at the ratio--
look at differences.
00:40:35.980 --> 00:40:40.000
Look at the ratio of--
so what comes in here
00:40:40.000 --> 00:40:46.160
will be the ratio
r, say, j plus 1/2.
00:40:46.160 --> 00:40:53.610
That'll be the ratio of U,
say j plus 2 minus U_(j+1) 1
00:40:53.610 --> 00:40:58.090
to U_(j+1) minus U_j.
00:40:58.090 --> 00:41:00.910
All at time n.
00:41:00.910 --> 00:41:04.250
That ratio tells us--
that ratio is close to 1,
00:41:04.250 --> 00:41:16.870
if-- so the ratio r is close to
1 for smooth and our function,
00:41:16.870 --> 00:41:23.620
phi of r then, we want our
function phi, when the ratio is
00:41:23.620 --> 00:41:26.630
1, to give the value 1.
00:41:26.630 --> 00:41:29.270
So if I try to draw
a graph of-- let
00:41:29.270 --> 00:41:43.960
me draw a graph of people's flux
functions-- of the phi of r.
00:41:43.960 --> 00:41:47.180
So here's a graph of phi of r.
00:41:47.180 --> 00:41:50.230
First on the graph I'll
put upwind and Lax-Wendroff
00:41:50.230 --> 00:41:51.890
without anything.
00:41:51.890 --> 00:41:59.640
So Lax-Wendroff would take-- phi
Lax-Wendroff is identically 1.
00:41:59.640 --> 00:42:02.660
This is the ratio
r and this is phi.
00:42:06.490 --> 00:42:15.250
So Lax-Wendroff, before fixing,
had no flux limiter function
00:42:15.250 --> 00:42:17.100
and it just had 1.
00:42:17.100 --> 00:42:22.680
And actually, if I can
keep the same formula,
00:42:22.680 --> 00:42:27.400
phi_upwind would just be zero.
00:42:27.400 --> 00:42:34.280
Then it'll turn out that--
let me draw these ratios.
00:42:34.280 --> 00:42:45.750
There is a region here and
here where the flux limiter
00:42:45.750 --> 00:42:47.390
is doing the right thing.
00:42:47.390 --> 00:42:53.710
So in that region it's
acting the right way.
00:42:53.710 --> 00:42:57.350
It's limiting the
anti-diffusive flux.
00:42:57.350 --> 00:42:59.670
It's controlling the variation.
00:42:59.670 --> 00:43:01.710
It's controlling
the oscillation.
00:43:01.710 --> 00:43:11.020
So various people have
proposed-- let's see.
00:43:11.020 --> 00:43:12.760
No, I haven't got this right.
00:43:16.180 --> 00:43:20.330
I think it's-- because I
know that part is wrong
00:43:20.330 --> 00:43:21.290
for Lax-Wendroff.
00:43:21.290 --> 00:43:22.270
So maybe it's there.
00:43:25.200 --> 00:43:28.790
No, I've forgotten--
sorry that I'm
00:43:28.790 --> 00:43:41.210
not-- I've lost this because
I know that one possibility is
00:43:41.210 --> 00:43:44.850
this as the flux limiter.
00:43:44.850 --> 00:43:48.850
One possibility
would be this one
00:43:48.850 --> 00:43:51.600
and then there's
another possibility that
00:43:51.600 --> 00:43:54.440
is OK that goes through here.
00:43:59.950 --> 00:44:06.770
Let me stop for today with
this idea of a flux limiter
00:44:06.770 --> 00:44:12.380
and you see that it's a
task to fix the method
00:44:12.380 --> 00:44:15.650
to do the right
thing at smooth parts
00:44:15.650 --> 00:44:18.090
and also the right thing
at discontinuities.
00:44:18.090 --> 00:44:20.130
OK, I'll pick up
on that next time.
00:44:20.130 --> 00:44:21.380
Thanks.