1
00:00:00 --> 00:00:01
2
00:00:01 --> 00:00:02
The following content is
provided under a Creative
3
00:00:02 --> 00:00:03
Commons license.
4
00:00:03 --> 00:00:06
Your support will help MIT
OpenCourseWare continue to
5
00:00:06 --> 00:00:09
offer high-quality educational
resources for free.
6
00:00:09 --> 00:00:12
To make a donation, or to view
additional materials from
7
00:00:12 --> 00:00:16
hundreds of MIT courses visit
MIT OpenCourseWare
8
00:00:16 --> 00:00:19
at ocw.mit.edu.
9
00:00:19 --> 00:00:22
PROFESSOR STRANG: Ready?
10
00:00:22 --> 00:00:29
So let me recap the important
lecture from last time that
11
00:00:29 --> 00:00:33
gave the framework that we'll
see in Chapter two and Chapter
12
00:00:33 --> 00:00:35
three, really a major
part of the course.
13
00:00:35 --> 00:00:39
And our example was a line
of masses and springs.
14
00:00:39 --> 00:00:43
So that's like, the first
example to start with.
15
00:00:43 --> 00:00:46
And you remember, there
were three steps.
16
00:00:46 --> 00:00:50
If we want to know how much
the masses displace we have
17
00:00:50 --> 00:00:54
to get to the springs.
18
00:00:54 --> 00:00:56
Difference in displacements
gave the stretching
19
00:00:56 --> 00:00:58
of the springs.
20
00:00:58 --> 00:01:01
Then comes the material law,
Hooke's Law in this case with
21
00:01:01 --> 00:01:05
a matrix C that led to CAu.
22
00:01:06 --> 00:01:09
And then finally, those spring
forces are balanced by
23
00:01:09 --> 00:01:11
the external forces.
24
00:01:11 --> 00:01:13
That brought in A transpose.
25
00:01:13 --> 00:01:17
So that's the picture
to look for.
26
00:01:17 --> 00:01:27
And I just mention it again
because it's so central.
27
00:01:27 --> 00:01:30
So we'll be doing other
examples of that.
28
00:01:30 --> 00:01:35
But I wanted to stay with
springs and masses for today.
29
00:01:35 --> 00:01:37
To do another type of problem.
30
00:01:37 --> 00:01:40
A problem in which time enters.
31
00:01:40 --> 00:01:44
A problem in which we're not
looking for the steady state.
32
00:01:44 --> 00:01:47
But it's Newton's Law
is going to come in.
33
00:01:47 --> 00:01:49
This is actually Newton's Law.
34
00:01:49 --> 00:01:51
You recognize mass.
35
00:01:51 --> 00:01:54
You recognize acceleration.
36
00:01:54 --> 00:01:58
And the force on
the spring is -Ku.
37
00:01:59 --> 00:02:05
So when - Ku, when the force
came over to that side, that's
38
00:02:05 --> 00:02:08
the equation we're looking at.
39
00:02:08 --> 00:02:13
So it's like time-out for
a discussion of the very
40
00:02:13 --> 00:02:16
important subject of
time derivatives.
41
00:02:16 --> 00:02:24
How do I solve initial
value problems now?
42
00:02:24 --> 00:02:31
So I have some starting time,
I'm given u(0), the starting
43
00:02:31 --> 00:02:37
position of the springs, and
I'm given the velocity at zero.
44
00:02:37 --> 00:02:41
So I'm given two facts at zero.
45
00:02:41 --> 00:02:46
That's very different from
the steady state problem
46
00:02:46 --> 00:02:48
where you're given stuff
around the boundary.
47
00:02:48 --> 00:02:50
Here we're given
stuff at the start.
48
00:02:50 --> 00:02:53
Initial values instead
of boundary values.
49
00:02:53 --> 00:02:55
So what does this mean?
50
00:02:55 --> 00:03:00
It means that maybe I pull the
springs down and I let go.
51
00:03:00 --> 00:03:02
And what do they do?
52
00:03:02 --> 00:03:02
They oscillate.
53
00:03:02 --> 00:03:05
And you will have
studied this equation.
54
00:03:05 --> 00:03:08
I hope you've seen
that equation.
55
00:03:08 --> 00:03:10
Because I think of it
as the fundamental
56
00:03:10 --> 00:03:12
equation of mechanics.
57
00:03:12 --> 00:03:14
In fact, I've even
used those words.
58
00:03:14 --> 00:03:17
The fundamental equation of
mechanical engineering.
59
00:03:17 --> 00:03:24
Well, it's Newton's
Law for a system.
60
00:03:24 --> 00:03:29
So the matrix M now
has two masses.
61
00:03:29 --> 00:03:32
In this case, it would
be a two by two system.
62
00:03:32 --> 00:03:35
The matrix M would
have two masses.
63
00:03:35 --> 00:03:42
Suppose mass one may be like,
nine and mass two is one.
64
00:03:42 --> 00:03:48
So copying what our problem is
we'd have a nine and a one
65
00:03:48 --> 00:03:52
multiplying u_1'' and u_2''.
66
00:03:53 --> 00:03:57
Can I use prime
for a derivative?
67
00:03:57 --> 00:04:00
Rather than dots, because
I'm never too sure if
68
00:04:00 --> 00:04:01
the dot is there or not.
69
00:04:01 --> 00:04:04
So prime I can see.
70
00:04:04 --> 00:04:07
And then plus our Ku.
71
00:04:08 --> 00:04:13
And actually you know what K
will look like for this system.
72
00:04:13 --> 00:04:17
It's fixed-fixed.
73
00:04:17 --> 00:04:20
Well, let me write down
what K would look like.
74
00:04:20 --> 00:04:23
Just because we saw
it at the very end.
75
00:04:23 --> 00:04:25
What it looks like if there
are spring constants,
76
00:04:25 --> 00:04:25
c_1, c_2, and c_3.
77
00:04:27 --> 00:04:33
Can I just remind you what
that will look like?
78
00:04:33 --> 00:04:36
So this will multiply
u and give zero.
79
00:04:36 --> 00:04:43
So the mass matrix is
diagonal in the problem.
80
00:04:43 --> 00:04:47
That's diagonal matrices, we
certainly expect to be easy.
81
00:04:47 --> 00:04:51
Almost as easy as the identity
matrix where we would
82
00:04:51 --> 00:04:53
have two equal masses.
83
00:04:53 --> 00:04:55
Here it's just that little
bit more general with
84
00:04:55 --> 00:04:57
two different masses.
85
00:04:57 --> 00:05:00
And in here, the kind of thing
that we created just at the end
86
00:05:00 --> 00:05:05
of the last lecture was
something like c_1+c_2 on the
87
00:05:05 --> 00:05:11
diagonal and a -c_2 and
a -c_2 and a c_2 + c_3.
88
00:05:11 --> 00:05:16
89
00:05:16 --> 00:05:20
So c_2 for the middle spring is
the one that's really there
90
00:05:20 --> 00:05:23
with its full little element
matrix, c_2, -c_2, -c_2, c_2.
91
00:05:23 --> 00:05:26
92
00:05:26 --> 00:05:32
The c_1 part is only coming
in at the top because it's
93
00:05:32 --> 00:05:34
connected to a fixed support.
94
00:05:34 --> 00:05:37
And c_3 is only coming
in at the bottom.
95
00:05:37 --> 00:05:40
Anyway, that's the matrix K.
96
00:05:40 --> 00:05:42
This is the matrix M.
97
00:05:42 --> 00:05:46
And how do we solve
the problem?
98
00:05:46 --> 00:05:50
You've seen before and now
ADINA Professor Bathe's code or
99
00:05:50 --> 00:05:55
any other code would offer,
like two different ways.
100
00:05:55 --> 00:06:02
One way is because we have
constant coefficients here we
101
00:06:02 --> 00:06:07
can expect a pretty clean
formula for the answer.
102
00:06:07 --> 00:06:12
And because we're growing,
we're oscillating in time,
103
00:06:12 --> 00:06:16
things are happening in time,
we expect eigenvectors,
104
00:06:16 --> 00:06:18
eigenvalues to come in.
105
00:06:18 --> 00:06:19
Into that formula.
106
00:06:19 --> 00:06:22
I think it's worth just
repeating that formula.
107
00:06:22 --> 00:06:24
How do you approach it?
108
00:06:24 --> 00:06:30
And then the second, the really
serious issue is how do you
109
00:06:30 --> 00:06:37
solve equation like this
or even more general,
110
00:06:37 --> 00:06:40
time-dependant problems.
111
00:06:40 --> 00:06:43
I mean this is what finite
element codes are created for.
112
00:06:43 --> 00:06:47
How do you study the
crash of a car?
113
00:06:47 --> 00:06:52
So this isn't quite the
equation for a car crash.
114
00:06:52 --> 00:06:57
What would be different
in a car crash?
115
00:06:57 --> 00:07:02
There's a short section later
in Chapter two called the
116
00:07:02 --> 00:07:05
reality of computational
engineering.
117
00:07:05 --> 00:07:10
And the reality is, cars crash,
people drop their cell phones.
118
00:07:10 --> 00:07:15
And those are very difficult
problems to compute with.
119
00:07:15 --> 00:07:18
And of course, absolutely
impossible to use eigenvectors
120
00:07:18 --> 00:07:22
because, well they're
not linear at all.
121
00:07:22 --> 00:07:23
Right?
122
00:07:23 --> 00:07:27
If you crash two cars, if you
want to study, everything's
123
00:07:27 --> 00:07:30
happening in like, hundredth
of a second or something.
124
00:07:30 --> 00:07:35
But what's happening there is
highly non-linear so it's
125
00:07:35 --> 00:07:37
very much more complicated.
126
00:07:37 --> 00:07:40
And you take thousands
of time steps maybe
127
00:07:40 --> 00:07:43
within the crash time.
128
00:07:43 --> 00:07:52
And you're going to use finite
differences or finite elements.
129
00:07:52 --> 00:07:56
So this is a chance to say
something about this special
130
00:07:56 --> 00:08:02
equation for finite
difference method.
131
00:08:02 --> 00:08:06
It's really 18.086 that takes
up these questions seriously.
132
00:08:06 --> 00:08:08
If I choose a finite
difference method.
133
00:08:08 --> 00:08:11
See, the thing is with finite
differences you've got
134
00:08:11 --> 00:08:14
many, many choices.
135
00:08:14 --> 00:08:17
There's only one way
to go in step one.
136
00:08:17 --> 00:08:21
Eigenvectors, this thing
has got eigenvectors, you
137
00:08:21 --> 00:08:23
find them, you're in.
138
00:08:23 --> 00:08:26
But for the much more general,
typical problems that you'll be
139
00:08:26 --> 00:08:30
solving the rest of your lives,
finite differences or finite
140
00:08:30 --> 00:08:35
elements, you've got
lots of choices.
141
00:08:35 --> 00:08:41
And the issues that come up are
the accuracy of the choice.
142
00:08:41 --> 00:08:44
Center differences often
gives you that extra
143
00:08:44 --> 00:08:45
order of accuracy.
144
00:08:45 --> 00:08:50
Stability is something
we have not seen yet.
145
00:08:50 --> 00:08:52
So what does stability mean?
146
00:08:52 --> 00:08:56
Stability means that as I
follow my difference equation
147
00:08:56 --> 00:09:02
forward in time it stays
near the true equation.
148
00:09:02 --> 00:09:04
You'll see by example.
149
00:09:04 --> 00:09:06
So some methods are more
stable than others,
150
00:09:06 --> 00:09:09
some are completely
unstable and unusable.
151
00:09:09 --> 00:09:16
And then, of course, the other
condition, another desired
152
00:09:16 --> 00:09:19
requirement is speed
of calculation.
153
00:09:19 --> 00:09:21
So these balance each other.
154
00:09:21 --> 00:09:23
It's a beautiful subject.
155
00:09:23 --> 00:09:28
And let me just open the
door to that subject today.
156
00:09:28 --> 00:09:32
May I start with the
eigenvector solution.
157
00:09:32 --> 00:09:37
How would I solve this
equation by eigenvectors?
158
00:09:37 --> 00:09:42
Well again, I'm going to
look for special solutions.
159
00:09:42 --> 00:09:45
So let me write that
equation down again.
160
00:09:45 --> 00:09:45
Mu''+Ku=0.
161
00:09:45 --> 00:09:49
162
00:09:49 --> 00:09:53
I'm taking zero external force.
163
00:09:53 --> 00:09:57
So the springs and masses
are just oscillating.
164
00:09:57 --> 00:10:00
Their total energy won't
change, it's a closed system.
165
00:10:00 --> 00:10:04
The kinetic plus the potential
energy will stay constant
166
00:10:04 --> 00:10:07
and we can show why.
167
00:10:07 --> 00:10:11
There's a lot in this section,
by the way, section 2.2.
168
00:10:11 --> 00:10:13
More than I would be able
to do in the lecture.
169
00:10:13 --> 00:10:17
But let me capture
the key ideas.
170
00:10:17 --> 00:10:21
So one key idea is the
natural idea when we have
171
00:10:21 --> 00:10:24
constant coefficients.
172
00:10:24 --> 00:10:26
Look for special solutions.
173
00:10:26 --> 00:10:32
So let me say look for special
solutions of the form-- well,
174
00:10:32 --> 00:10:36
let's keep in mind always
the simplest model of all.
175
00:10:36 --> 00:10:41
Let me put it over here
because I'll focus on it
176
00:10:41 --> 00:10:43
particularly, especially.
177
00:10:43 --> 00:10:48
The simplest case would be
u'' with a mass normalized
178
00:10:48 --> 00:10:50
to one plus u=0.
179
00:10:52 --> 00:10:54
Just a single equation.
180
00:10:54 --> 00:10:57
So that's just a single spring.
181
00:10:57 --> 00:10:59
Single mass.
182
00:10:59 --> 00:11:02
So the mass is just
oscillating up and down.
183
00:11:02 --> 00:11:08
And we all know that it's going
to oscillate up and down, sines
184
00:11:08 --> 00:11:10
and cosines are going
to come into here.
185
00:11:10 --> 00:11:14
In fact, we know that the
solution to that would be, this
186
00:11:14 --> 00:11:22
could have cosine, this would
have some cos(t)'s and some, u
187
00:11:22 --> 00:11:27
could be a combination
of cos(t) and sin(t).
188
00:11:29 --> 00:11:30
Right?
189
00:11:30 --> 00:11:34
And the A and B would be
used to match the two
190
00:11:34 --> 00:11:35
initial conditions.
191
00:11:35 --> 00:11:38
So that's what we've got
in the scalar case.
192
00:11:38 --> 00:11:44
This is one spring. n is one.
193
00:11:44 --> 00:11:47
We don't have a matrix
here, just a number.
194
00:11:47 --> 00:11:53
But that's our guide
to the matrix case.
195
00:11:53 --> 00:11:54
So now what am I
going to look for?
196
00:11:54 --> 00:11:55
I'm going to look for u(t).
197
00:11:56 --> 00:12:05
I'll look for special solutions
of the form cos, they'll go
198
00:12:05 --> 00:12:06
at different frequencies.
199
00:12:06 --> 00:12:09
So the frequency
has to be found.
200
00:12:09 --> 00:12:12
So that's the time dependence.
201
00:12:12 --> 00:12:17
And then some, if I'm lucky
all the springs will
202
00:12:17 --> 00:12:20
be going together.
203
00:12:20 --> 00:12:25
So this x-- Am I going to
use x for eigenvector?
204
00:12:25 --> 00:12:27
I think so, yes.
205
00:12:27 --> 00:12:30
Let's use x for eigenvector.
206
00:12:30 --> 00:12:37
So this is a vector, this
is a constant vector.
207
00:12:37 --> 00:12:49
And this gives the oscillation.
208
00:12:49 --> 00:12:55
Let me just already
start with that.
209
00:12:55 --> 00:12:57
I've got two unknowns
here, right?
210
00:12:57 --> 00:12:59
Two things that I'm
free to choose.
211
00:12:59 --> 00:13:03
I'm free to choose the
frequency omega and I'm free
212
00:13:03 --> 00:13:05
to choose this vector x.
213
00:13:05 --> 00:13:08
And I've like,
separated out time.
214
00:13:08 --> 00:13:10
The way we've been doing
it with an e^(lambda*).
215
00:13:12 --> 00:13:13
We used an e^(lambda*t).
216
00:13:14 --> 00:13:18
So that was sort of right
for a first order equation.
217
00:13:18 --> 00:13:23
Cosine is right for a second
order equation like this.
218
00:13:23 --> 00:13:26
When is that a solution
to my equation?
219
00:13:26 --> 00:13:29
May I just plug it in?
220
00:13:29 --> 00:13:34
So I'm going to
just plug it in?
221
00:13:34 --> 00:13:36
Substitute is a better
word than plug in, right?
222
00:13:36 --> 00:13:42
So shall I just plug it in?
223
00:13:42 --> 00:13:46
So I get M times that
second derivative.
224
00:13:46 --> 00:13:50
So the second derivative of
the cosine is the cosine with
225
00:13:50 --> 00:13:54
a minus M omega squared.
226
00:13:54 --> 00:13:56
Is that right?
227
00:13:56 --> 00:13:57
Yes?
228
00:13:57 --> 00:13:58
Times what?
229
00:13:58 --> 00:14:01
So that's what came from
the time derivative.
230
00:14:01 --> 00:14:04
Times cos times this,
cos(omega*t)x.
231
00:14:04 --> 00:14:08
232
00:14:08 --> 00:14:11
Two time derivatives brought
down a minus omega squared.
233
00:14:11 --> 00:14:16
And the other part is just
K times u, cos(omega*t)x.
234
00:14:16 --> 00:14:19
235
00:14:19 --> 00:14:24
And that should be zero.
236
00:14:24 --> 00:14:27
And this is supposed to
be a good solution that
237
00:14:27 --> 00:14:33
works for all time.
238
00:14:33 --> 00:14:37
This is what I want to be zero.
239
00:14:37 --> 00:14:39
Do you see what's left?
240
00:14:39 --> 00:14:45
Kx from here.
241
00:14:45 --> 00:14:48
And putting that on the other
side, M(omega squared)x.
242
00:14:48 --> 00:14:53
243
00:14:53 --> 00:14:55
Can I write the omega squared
that way? (omega squared)Mx.
244
00:14:55 --> 00:14:59
245
00:14:59 --> 00:15:04
Then if that is satisfied
then the differential
246
00:15:04 --> 00:15:05
equation is satisfied.
247
00:15:05 --> 00:15:06
I've got a solution.
248
00:15:06 --> 00:15:10
So I'm looking for K and I'm
looking for x and omega.
249
00:15:10 --> 00:15:16
And by the way, I could also
have sin(omega*t) here.
250
00:15:16 --> 00:15:19
As well as cosine just
the way over there.
251
00:15:19 --> 00:15:22
Sines and cosines
both possible.
252
00:15:22 --> 00:15:22
So sin(omega*t)x.
253
00:15:22 --> 00:15:26
254
00:15:26 --> 00:15:28
It would be exactly the same
except that it'd be a
255
00:15:28 --> 00:15:31
sin(omega*t) that I'd
be dividing out.
256
00:15:31 --> 00:15:33
And again, I'd come
back to this.
257
00:15:33 --> 00:15:39
This is key problem.
258
00:15:39 --> 00:15:43
And do you see that somehow
here we have an eigenvector
259
00:15:43 --> 00:15:48
x and an eigenvalue
omega squared?
260
00:15:48 --> 00:15:54
An eigenvalue omega
squared, right.
261
00:15:54 --> 00:15:55
But there is one little twist.
262
00:15:55 --> 00:15:59
It's not quite our standard
eigenvalue problem.
263
00:15:59 --> 00:16:04
What's the extra guy that's
present in that box that
264
00:16:04 --> 00:16:07
we don't usually see?
265
00:16:07 --> 00:16:08
M.
266
00:16:08 --> 00:16:11
M is the new person
there, new thing.
267
00:16:11 --> 00:16:13
The mass matrix.
268
00:16:13 --> 00:16:17
Which, if all masses where
one, the mass matrix
269
00:16:17 --> 00:16:18
would be the identity.
270
00:16:18 --> 00:16:19
We would be back.
271
00:16:19 --> 00:16:21
We would just have our
standard problem.
272
00:16:21 --> 00:16:26
I just have to say a word
about the case when the
273
00:16:26 --> 00:16:29
mass matrix is there.
274
00:16:29 --> 00:16:32
It's still an
eigenvalue problem.
275
00:16:32 --> 00:16:37
If you like, I can bring
M inverse over here.
276
00:16:37 --> 00:16:38
I could write it as (M
inverse*K)x=(omega squared)x.
277
00:16:38 --> 00:16:43
278
00:16:43 --> 00:16:47
So I'm looking for the
eigenvalues of M inverse K.
279
00:16:49 --> 00:16:51
That's really what I'm doing.
280
00:16:51 --> 00:16:53
The eigenvalues and
eigenvectors of M inverse K.
281
00:16:55 --> 00:17:03
But I'm always trying to be
like, aesthetically right.
282
00:17:03 --> 00:17:09
And what's M inverse K, has
just a little something
283
00:17:09 --> 00:17:09
wrong with it.
284
00:17:09 --> 00:17:11
Which is?
285
00:17:11 --> 00:17:14
It's not symmetric.
286
00:17:14 --> 00:17:18
If I take my matrix K, which is
beautifully symmetric, and my
287
00:17:18 --> 00:17:21
M, which is beautifully
symmetric, but when I do M
288
00:17:21 --> 00:17:23
inverse K, do see
what'll happen?
289
00:17:23 --> 00:17:27
The inverse of that
matrix will have a 1/9.
290
00:17:27 --> 00:17:31
When I do M inverse*K, that
row will get divided by nine.
291
00:17:31 --> 00:17:35
And the second row won't change
because I've got a one.
292
00:17:35 --> 00:17:39
So it just like, spoiled
things a little.
293
00:17:39 --> 00:17:41
You can deal with it.
294
00:17:41 --> 00:17:46
But actually, in some way it's
better to keep it right.
295
00:17:46 --> 00:17:50
And MATLAB is totally
okay with that.
296
00:17:50 --> 00:17:53
So MATLAB would use the
command eig, and other
297
00:17:53 --> 00:18:00
systems too, of K, M.
298
00:18:00 --> 00:18:03
Just tell MATLAB
the two matrices.
299
00:18:03 --> 00:18:06
Then it solves that problem.
300
00:18:06 --> 00:18:10
It will print out the, well if
you just typed eig, it would
301
00:18:10 --> 00:18:13
print out the eigenvalues.
302
00:18:13 --> 00:18:16
If you ask it for the
eigenvector there matrix,
303
00:18:16 --> 00:18:18
it'll print those too.
304
00:18:18 --> 00:18:24
So that has the grand name
generalized eigenvalue problem.
305
00:18:24 --> 00:18:27
Generalized because it's
got an M in there.
306
00:18:27 --> 00:18:30
I don't know if you've
met these problems where
307
00:18:30 --> 00:18:32
there could be an M.
308
00:18:32 --> 00:18:35
It's just a natural, and it's
not really a big deal,
309
00:18:35 --> 00:18:39
particularly when M is a
positive diagonal matrix.
310
00:18:39 --> 00:18:40
It's not a big deal at all.
311
00:18:40 --> 00:18:45
It's the same codes essentially
finding the eigenvalues
312
00:18:45 --> 00:18:48
and eigenvectors.
313
00:18:48 --> 00:18:52
So suppose we have them.
314
00:18:52 --> 00:19:01
We expect n positive
eigenvalues, and those'll
315
00:19:01 --> 00:19:04
be the omega squareds.
316
00:19:04 --> 00:19:08
And each one comes
with its eigenvector.
317
00:19:08 --> 00:19:14
And complete set,
everything is good.
318
00:19:14 --> 00:19:17
We're talking about symmetric
positive definite matrices.
319
00:19:17 --> 00:19:21
Notice that K is
positive definite.
320
00:19:21 --> 00:19:25
That was what our whole last
weeks have been about.
321
00:19:25 --> 00:19:28
And M is obviously
positive definite.
322
00:19:28 --> 00:19:34
So it's a complete, it's
a perfect set-up for
323
00:19:34 --> 00:19:37
the eigenvalue problem.
324
00:19:37 --> 00:19:41
I just want to think
through the final step.
325
00:19:41 --> 00:19:46
After we find the eignevalues,
lambda, the omega squareds,
326
00:19:46 --> 00:19:51
then we know the omegas and we
find the eigenvectors, x, how
327
00:19:51 --> 00:19:56
do you write the answer?
328
00:19:56 --> 00:19:59
So what have I done so far?
329
00:19:59 --> 00:20:03
I've looked for some
special solutions.
330
00:20:03 --> 00:20:06
And I will find them.
331
00:20:06 --> 00:20:09
So I've found these omegas.
332
00:20:09 --> 00:20:13
They can go with sines or
cosines, I found the x's.
333
00:20:13 --> 00:20:17
Now what's the
general solution?
334
00:20:17 --> 00:20:21
If I have these solutions to my
great equation there, what's
335
00:20:21 --> 00:20:23
the general solution? u(t).
336
00:20:23 --> 00:20:30
So this would be
the big picture.
337
00:20:30 --> 00:20:35
What can I do in
linear equations?
338
00:20:35 --> 00:20:37
Linear combinations.
339
00:20:37 --> 00:20:40
That's what linear algebra is
all about, linear combinations.
340
00:20:40 --> 00:20:44
So I could take a linear
combination of these cosines,
341
00:20:44 --> 00:20:45
so cos(omega_1*t)x_1.
342
00:20:45 --> 00:20:49
343
00:20:49 --> 00:20:52
So that would be a solution
coming from the first time
344
00:20:52 --> 00:20:56
eigenvector and eigenvalue
or the square root
345
00:20:56 --> 00:20:58
times any constant.
346
00:20:58 --> 00:21:00
And then, of course,
I could have a b_1.
347
00:21:01 --> 00:21:05
I'll use b's for the
sines. (omega_1*t)x_1.
348
00:21:05 --> 00:21:08
349
00:21:08 --> 00:21:10
So that's just like that.
350
00:21:10 --> 00:21:15
But it's used the sines, which
would also solve the equation.
351
00:21:15 --> 00:21:19
And then all the
way down to x_n's.
352
00:21:19 --> 00:21:21
Right?
353
00:21:21 --> 00:21:21
So I've got 2n.
354
00:21:23 --> 00:21:30
2n simple solutions. n cosines
times eigenvectors and n sines
355
00:21:30 --> 00:21:32
times those same eigenvectors.
356
00:21:32 --> 00:21:34
Why do I want 2n?
357
00:21:36 --> 00:21:41
I've got 2n constants then to
choose, the a's and the b's.
358
00:21:41 --> 00:21:44
And how do I choose them?
359
00:21:44 --> 00:21:53
What's the next step?
360
00:21:53 --> 00:21:56
Using eigenvectors, maybe
I just mention again
361
00:21:56 --> 00:22:03
the three steps.
362
00:22:03 --> 00:22:09
The three step method.
363
00:22:09 --> 00:22:19
Step one; find a's and b's.
364
00:22:19 --> 00:22:23
Oh, now you have to
answer my question.
365
00:22:23 --> 00:22:25
Where are the a's and
b's going to come from?
366
00:22:25 --> 00:22:28
What's going to determine
these 2n constants,
367
00:22:28 --> 00:22:30
the a's and the b's?
368
00:22:30 --> 00:22:32
The initial conditions,
of course.
369
00:22:32 --> 00:22:36
We've got two initial
conditions, this is a vector.
370
00:22:36 --> 00:22:37
We're talking about
a system here.
371
00:22:37 --> 00:22:41
I've got position
of both masses.
372
00:22:41 --> 00:22:45
I've got the initial
velocity of both masses.
373
00:22:45 --> 00:22:53
And the a's, they come from,
because the a's go with the
374
00:22:53 --> 00:22:57
cosines and so they're going
to be the ones associated
375
00:22:57 --> 00:23:01
with u(0), will give this.
376
00:23:01 --> 00:23:07
Will give a combination, will
give the a's. u(0) will
377
00:23:07 --> 00:23:08
be a_1*x_1 + a_n*x_n.
378
00:23:08 --> 00:23:15
379
00:23:15 --> 00:23:16
This is at t=0.
380
00:23:16 --> 00:23:19
381
00:23:19 --> 00:23:22
So I'm using the initial
conditions. u(0) is a
382
00:23:22 --> 00:23:26
combination of those
eigenvectors.
383
00:23:26 --> 00:23:32
So this is from the
initial conditions.
384
00:23:32 --> 00:23:36
And the b's of course are
going to come from u'(0),
385
00:23:36 --> 00:23:40
the initial velocity
because they go with sines.
386
00:23:40 --> 00:23:47
So sines start at zero but they
have an initial velocity.
387
00:23:47 --> 00:23:50
So that's step one.
388
00:23:50 --> 00:23:53
That finds the constants.
389
00:23:53 --> 00:23:57
It splits the problem
into these normal modes.
390
00:23:57 --> 00:24:00
It finds how much of
each eigenvector is
391
00:24:00 --> 00:24:04
in there at the Start
392
00:24:04 --> 00:24:05
Step two?
393
00:24:05 --> 00:24:08
So it's really worth seeing
these three steps because they
394
00:24:08 --> 00:24:12
just repeat and repeat for all
applications of eigenvalues.
395
00:24:12 --> 00:24:15
Step two is what?
396
00:24:15 --> 00:24:20
Step two is follow each of
these guys forward in time.
397
00:24:20 --> 00:24:22
Follow them to any time.
398
00:24:22 --> 00:24:23
What does that mean?
399
00:24:23 --> 00:24:28
That means just put in the
cosines and the sines.
400
00:24:28 --> 00:24:39
So the a's go to
a*cos(omega*t).
401
00:24:39 --> 00:24:44
402
00:24:44 --> 00:24:48
So I'm just saying what happens
to those coefficients.
403
00:24:48 --> 00:24:48
Each one.
404
00:24:48 --> 00:24:55
Let's see. a_i, say, goes
to a_i is multiplied
405
00:24:55 --> 00:24:57
by the cosine.
406
00:24:57 --> 00:25:02
And the b's, b_i's go
to b_i*sin(omega_i*t).
407
00:25:02 --> 00:25:08
408
00:25:08 --> 00:25:13
So now we followed the
coefficients forward in time.
409
00:25:13 --> 00:25:15
And step three?
410
00:25:15 --> 00:25:19
The final step is?
411
00:25:19 --> 00:25:21
So here, we're
following each one.
412
00:25:21 --> 00:25:25
So the pattern is
always the same one.
413
00:25:25 --> 00:25:33
Take your starting conditions,
split them up into
414
00:25:33 --> 00:25:34
eigenvectors.
415
00:25:34 --> 00:25:37
Follow each eigenvector.
416
00:25:37 --> 00:25:39
And then step three is?
417
00:25:39 --> 00:25:40
Reassemble.
418
00:25:40 --> 00:25:42
Put them back together.
419
00:25:42 --> 00:25:48
Step three is the solution.
u at that later time t is
420
00:25:48 --> 00:25:49
a_1*cos(omega_1*t)x_1.
421
00:25:54 --> 00:25:55
And b_1*sin(omega_1*t)x_1.
422
00:26:00 --> 00:26:06
Those are the two guys
that are reflecting
423
00:26:06 --> 00:26:08
the fundamental mode.
424
00:26:08 --> 00:26:14
And then we have a_2's and
b_2's multiplied by cosines and
425
00:26:14 --> 00:26:18
sines and times x_2 and so on.
426
00:26:18 --> 00:26:22
Put the pieces together.
427
00:26:22 --> 00:26:25
So it's split the initial
conditions, follow each
428
00:26:25 --> 00:26:31
eigenvector, reassemble.
429
00:26:31 --> 00:26:36
It's the heart of Fourier
methods, it's the heart
430
00:26:36 --> 00:26:38
of using eigenvectors.
431
00:26:38 --> 00:26:43
And there's a little
code in the book that
432
00:26:43 --> 00:26:49
does exactly that.
433
00:26:49 --> 00:26:51
So let's just pause a moment.
434
00:26:51 --> 00:26:54
This was the
eigenvector solution.
435
00:26:54 --> 00:27:04
And now I want to talk about
computational science.
436
00:27:04 --> 00:27:07
For which this problem is a
model, but the problem has
437
00:27:07 --> 00:27:13
probably got more complications
and things change
438
00:27:13 --> 00:27:19
in time, whatever.
439
00:27:19 --> 00:27:22
This was an exact solution.
440
00:27:22 --> 00:27:28
That's more than we hope
for in real computations.
441
00:27:28 --> 00:27:32
We hope for accuracy, but we
don't hope for 100% accuracy
442
00:27:32 --> 00:27:38
unless we have this
exact model problem.
443
00:27:38 --> 00:27:40
So we use finite differences.
444
00:27:40 --> 00:27:45
Ready for finite differences?
445
00:27:45 --> 00:27:50
So this is method one, which
allows you to understand
446
00:27:50 --> 00:27:52
the model problem.
447
00:27:52 --> 00:27:56
And now comes method two, which
allows you to compute much
448
00:27:56 --> 00:28:00
more, many more problems.
449
00:28:00 --> 00:28:03
And if I do it just for
this, let me come back
450
00:28:03 --> 00:28:06
to just this model.
451
00:28:06 --> 00:28:14
One spring and one mass.
452
00:28:14 --> 00:28:22
Shall we start the
mass at rest?
453
00:28:22 --> 00:28:26
Let me make that the answer.
454
00:28:26 --> 00:28:28
I want that to be the answer.
455
00:28:28 --> 00:28:34
So I'm going to start with u(0)
to be zero-- to be one, right?
456
00:28:34 --> 00:28:39
I'm going to start with u(0) to
be one and u'(0) to be zero.
457
00:28:39 --> 00:28:42
Because the cosine starts
at one and its derivative
458
00:28:42 --> 00:28:47
starts at zero.
459
00:28:47 --> 00:28:50
Apologies that this is
such a simple model.
460
00:28:50 --> 00:28:54
I'm just pulling this this
mass down and let go.
461
00:28:54 --> 00:28:58
I just pull it down an
amount one, I let go,
462
00:28:58 --> 00:29:05
it oscillates forever.
463
00:29:05 --> 00:29:14
How do I draw the motion
of a single mass?
464
00:29:14 --> 00:29:16
A picture is important
and we have to think
465
00:29:16 --> 00:29:18
what's in the picture.
466
00:29:18 --> 00:29:26
So I think the picture should
be, maybe u in this direction
467
00:29:26 --> 00:29:29
and maybe u' in this direction.
468
00:29:29 --> 00:29:33
So it's the u, u' plane.
469
00:29:33 --> 00:29:36
The position velocity plane.
470
00:29:36 --> 00:29:38
Sometimes called
the phase plane.
471
00:29:38 --> 00:29:42
I'll just write that
word down, phase plane.
472
00:29:42 --> 00:29:45
And what will be
our picture here?
473
00:29:45 --> 00:29:50
It starts there, right? u(0) is
one, no velocity, starts there.
474
00:29:50 --> 00:29:52
Oh, and then it just
follows cos(t).
475
00:29:52 --> 00:29:55
And of course, I
know what u' is.
476
00:29:55 --> 00:29:59
If u is cos(t), u' is -sin(t).
477
00:29:59 --> 00:30:03
478
00:30:03 --> 00:30:07
Let me put that somewhere where
my picture won't run over it,
479
00:30:07 --> 00:30:10
up here maybe. u' is -sin(t).
480
00:30:10 --> 00:30:16
481
00:30:16 --> 00:30:18
Help me out.
482
00:30:18 --> 00:30:24
If as time increases, I
follow this point, the
483
00:30:24 --> 00:30:25
x coordinate is cos(t).
484
00:30:26 --> 00:30:28
The y coordinate is -sin(t).
485
00:30:29 --> 00:30:31
It's a circle.
486
00:30:31 --> 00:30:33
It's a circle.
487
00:30:33 --> 00:30:36
So I just buzz along it,
because of that minus sign, the
488
00:30:36 --> 00:30:37
circle goes this way around.
489
00:30:37 --> 00:30:41
It's a unit circle.
490
00:30:41 --> 00:30:45
Goes on forever.
491
00:30:45 --> 00:30:49
And actually cos squared plus
sine squared is one, of
492
00:30:49 --> 00:30:52
course, that's why it's
a circle of radius one.
493
00:30:52 --> 00:30:58
And that represents the
total energy actually.
494
00:30:58 --> 00:31:00
I'm trying to draw a circle.
495
00:31:00 --> 00:31:07
Actually how would you get your
computer to draw a circle?
496
00:31:07 --> 00:31:09
I'm thinking not
just one circle.
497
00:31:09 --> 00:31:11
I want to refresh it.
498
00:31:11 --> 00:31:14
I mean, I want to
follow this in time.
499
00:31:14 --> 00:31:23
So I want to draw, I want to
follow a point around a circle.
500
00:31:23 --> 00:31:26
Well I suppose what you
would do would be to ask
501
00:31:26 --> 00:31:29
it to plot those points.
502
00:31:29 --> 00:31:31
That's not allowed.
503
00:31:31 --> 00:31:37
I want you to follow the
solution to the equation.
504
00:31:37 --> 00:31:41
I want to take that equation
whose exact solution is a
505
00:31:41 --> 00:31:45
circle and I want to use
finite differences in time.
506
00:31:45 --> 00:31:48
So I'm going to
take finite steps.
507
00:31:48 --> 00:31:50
Of course that's what
the computer has to do.
508
00:31:50 --> 00:31:54
The computer will
take finite steps.
509
00:31:54 --> 00:31:57
So I'd be very happy
if those finite steps
510
00:31:57 --> 00:32:00
keep me on the circle.
511
00:32:00 --> 00:32:03
And I'd be even happier if it
came around exactly right at
512
00:32:03 --> 00:32:12
2pi but it's not going to.
513
00:32:12 --> 00:32:15
I want finite differences.
514
00:32:15 --> 00:32:17
I want to replace that
differential equation
515
00:32:17 --> 00:32:19
by finite differences.
516
00:32:19 --> 00:32:24
And when we do this equation,
we're practically doing
517
00:32:24 --> 00:32:26
that one at the same time.
518
00:32:26 --> 00:32:33
So this model will be fine.
519
00:32:33 --> 00:32:35
So I'm going to introduce
finite differences
520
00:32:35 --> 00:32:37
for this equation.
521
00:32:37 --> 00:32:40
And everybody, we've been
talking about how to replace
522
00:32:40 --> 00:32:43
second derivatives by
second differences.
523
00:32:43 --> 00:32:49
So that looks good.
524
00:32:49 --> 00:32:53
Let me write down, maybe I can
write down three or four
525
00:32:53 --> 00:32:57
possible finite differences.
526
00:32:57 --> 00:33:01
For u'' I think I'll write
down u_(n+1)-2u_n+u_(n-1).
527
00:33:07 --> 00:33:10
That's the second difference
that replaces the second
528
00:33:10 --> 00:33:12
derivative divided by what?
529
00:33:12 --> 00:33:20
What do you think goes
in the denominator now?
530
00:33:20 --> 00:33:21
Square of what?
531
00:33:21 --> 00:33:23
Yeah, now what's
the right name?
532
00:33:23 --> 00:33:25
It's going to be
something squared.
533
00:33:25 --> 00:33:28
What do I put down there?
534
00:33:28 --> 00:33:31
The step size.
535
00:33:31 --> 00:33:33
I don't want to call
it delta x, right?
536
00:33:33 --> 00:33:35
What should I call that?
537
00:33:35 --> 00:33:39
Delta t would be
natural, delta t.
538
00:33:39 --> 00:33:44
Or h, I could again use h
just to have a shorthand.
539
00:33:44 --> 00:33:47
That's my second difference.
540
00:33:47 --> 00:33:50
It's centered, it's
the natural choice.
541
00:33:50 --> 00:33:52
And not the only
choice, by the way.
542
00:33:52 --> 00:33:55
Oh.
543
00:33:55 --> 00:33:59
If I'm an astronomer or like,
Professor Wisdom here.
544
00:33:59 --> 00:34:05
So he followed Pluto
around, right?
545
00:34:05 --> 00:34:08
And discovered that Pluto
wasn't a planet, right?
546
00:34:08 --> 00:34:12
He discovered that the motion
of Pluto was not like,
547
00:34:12 --> 00:34:14
regular, like the Earth.
548
00:34:14 --> 00:34:18
I think the Earth is regular.
549
00:34:18 --> 00:34:20
But Pluto has chaotic motion.
550
00:34:20 --> 00:34:22
So it does crazy stuff.
551
00:34:22 --> 00:34:24
Right.
552
00:34:24 --> 00:34:30
Well we're looking for very
periodic, circular motion here.
553
00:34:30 --> 00:34:35
So what I was going to say,
Professor Wisdom, he would
554
00:34:35 --> 00:34:37
sneer at this second
differences, right?
555
00:34:37 --> 00:34:43
I mean, as we'll see that's
got decent accuracy, but not
556
00:34:43 --> 00:34:46
accuracy that would allow you
to follow Pluto for
557
00:34:46 --> 00:34:50
100 million years.
558
00:34:50 --> 00:34:51
But it'll do for us.
559
00:34:51 --> 00:34:54
We're not going 100
million years here.
560
00:34:54 --> 00:34:58
So now comes -u.
561
00:34:58 --> 00:35:02
Now, question.
562
00:35:02 --> 00:35:05
I'm taking this, the u
term and I'm putting
563
00:35:05 --> 00:35:07
it on the other side.
564
00:35:07 --> 00:35:10
So mass is one.
565
00:35:10 --> 00:35:12
Acceleration is this.
566
00:35:12 --> 00:35:15
The force is -u.
567
00:35:15 --> 00:35:23
Now comes the big question.
568
00:35:23 --> 00:35:26
For -u, do I use
the newest value?
569
00:35:26 --> 00:35:28
Do I use the middle value?
570
00:35:28 --> 00:35:31
Or do I use the oldest value?
571
00:35:31 --> 00:35:34
Or some combination?
572
00:35:34 --> 00:35:38
If I want high accuracy, oh,
then I go way up in these
573
00:35:38 --> 00:35:44
differences and I find tricky
combinations that get me above
574
00:35:44 --> 00:35:46
first and second
order accuracy.
575
00:35:46 --> 00:35:48
But let's not go there.
576
00:35:48 --> 00:35:53
So I've gotta make one
of those choices.
577
00:35:53 --> 00:35:58
And that choice is going to
decide, so I mean, this is a
578
00:35:58 --> 00:36:01
choice you have to make when
you have such a problem.
579
00:36:01 --> 00:36:04
Where are you going
to evaluate this?
580
00:36:04 --> 00:36:08
I guess one choice
jumps out as natural.
581
00:36:08 --> 00:36:16
What would you think of doing?
582
00:36:16 --> 00:36:20
Well how many would
do that one?
583
00:36:20 --> 00:36:22
So those are people,
that's the conservative
584
00:36:22 --> 00:36:28
choice, the stable.
585
00:36:28 --> 00:36:30
And then, how many
would do this?
586
00:36:30 --> 00:36:33
Yeah, you would do that, right.
587
00:36:33 --> 00:36:35
I would call that the
leapfrog choice.
588
00:36:35 --> 00:36:39
Somehow, this second
difference is leaping
589
00:36:39 --> 00:36:41
over this middle point.
590
00:36:41 --> 00:36:43
So that would be the
leapfrog method.
591
00:36:43 --> 00:36:52
And then if I use a very
low, to use the first
592
00:36:52 --> 00:36:54
value would be.
593
00:36:54 --> 00:36:58
So my point is those are
three different choices.
594
00:36:58 --> 00:37:04
And each one has these
questions of accuracy.
595
00:37:04 --> 00:37:08
The middle one will be more
accurate because it's centered.
596
00:37:08 --> 00:37:11
Each one has a question
of stability.
597
00:37:11 --> 00:37:14
Ah, you don't see
stability yet.
598
00:37:14 --> 00:37:18
So that's the point of the rest
of the lecture is to see.
599
00:37:18 --> 00:37:21
What is this
stability question?
600
00:37:21 --> 00:37:30
The issue of speed, speed would
be, this is the slow one.
601
00:37:30 --> 00:37:35
Because it involves, I have
to bring u_(n+1) over there.
602
00:37:35 --> 00:37:39
And if I've got a system of
equations and then I, it
603
00:37:39 --> 00:37:41
would take me some time.
604
00:37:41 --> 00:37:48
This would be called
an implicit method.
605
00:37:48 --> 00:37:52
The new u_(n+1) is only given
implicitly because it's showing
606
00:37:52 --> 00:37:54
up on the right-hand side.
607
00:37:54 --> 00:37:56
I've got to move it
over to the left.
608
00:37:56 --> 00:38:00
If it's non-linear I've
got a system to solve.
609
00:38:00 --> 00:38:05
It's safer, but more expensive.
610
00:38:05 --> 00:38:11
But this would be
the natural choice.
611
00:38:11 --> 00:38:19
I want to get eigenvalues
into this picture.
612
00:38:19 --> 00:38:23
So I'm going to do the step
we often do when we see
613
00:38:23 --> 00:38:26
a second order equation.
614
00:38:26 --> 00:38:30
That is, reduce it to two
first order equations.
615
00:38:30 --> 00:38:33
Can I do that and
see the same thing?
616
00:38:33 --> 00:38:36
So what would be my two
first order equations.
617
00:38:36 --> 00:38:38
My unknowns will be u and u'.
618
00:38:38 --> 00:38:45
So it'll be first order, the
derivative of u and u'.
619
00:38:46 --> 00:38:49
Shall I call it v for velocity?
620
00:38:49 --> 00:38:53
Let me write down, I don't
need to write this fancy.
621
00:38:53 --> 00:38:56
The two equations
will be, u'=v.
622
00:38:56 --> 00:38:59
623
00:38:59 --> 00:39:05
And v' is what?
624
00:39:05 --> 00:39:10
So u'=v sort of told me what v
was. v is the velocity, the
625
00:39:10 --> 00:39:11
time derivative of du/dt.
626
00:39:12 --> 00:39:15
And now the derivative of the
velocity, now that should
627
00:39:15 --> 00:39:17
reflect my true equation.
628
00:39:17 --> 00:39:20
So this is all coming
from u''+u=0.
629
00:39:20 --> 00:39:24
630
00:39:24 --> 00:39:28
So v' is what? v' is u''.
631
00:39:29 --> 00:39:32
Do I want -u there?
632
00:39:32 --> 00:39:38
Yeah, that's my equation.
v', which is u'' is -u.
633
00:39:39 --> 00:39:39
Good.
634
00:39:39 --> 00:39:41
That's my system.
635
00:39:41 --> 00:39:44
So I have a matrix here.
636
00:39:44 --> 00:39:51
I have a two by two system with
a matrix [0, 1; -1, 0] I guess.
637
00:39:51 --> 00:39:57
So while thinking about
difference methods, so this,
638
00:39:57 --> 00:39:59
I was thinking about a
second order equation.
639
00:39:59 --> 00:40:01
I went right to it.
640
00:40:01 --> 00:40:04
Here I'm thinking about
a first order system.
641
00:40:04 --> 00:40:07
A lot of problems come
in first order systems.
642
00:40:07 --> 00:40:13
Gas dynamics comes as a big
first order system for those
643
00:40:13 --> 00:40:20
components of mass, momentum,
energy, typically five
644
00:40:20 --> 00:40:23
non-linear equations
in gas dynamics.
645
00:40:23 --> 00:40:25
I mean, that's a very,
very serious problem,
646
00:40:25 --> 00:40:28
how to solve those.
647
00:40:28 --> 00:40:32
Here we've got a much simpler
model problem, linear,
648
00:40:32 --> 00:40:34
constant coefficients.
649
00:40:34 --> 00:40:40
But what do we do?
650
00:40:40 --> 00:40:43
Can I propose three
possibilities here?
651
00:40:43 --> 00:40:48
And actually they are identical
to these three possibilities.
652
00:40:48 --> 00:40:50
If I just switch over to
a first order system.
653
00:40:50 --> 00:40:58
So possibility one.
654
00:40:58 --> 00:41:04
Possibility one is that u_(n+1)
and v_(n+1), n so the u_ (n+1)
655
00:41:04 --> 00:41:09
would be, how do I model that
equation? u_(n+1) is
656
00:41:09 --> 00:41:09
u_n+delta t*v_n.
657
00:41:09 --> 00:41:15
658
00:41:15 --> 00:41:19
That would be a natural, right?
659
00:41:19 --> 00:41:23
Let me just pause there,
because that's meant to be the
660
00:41:23 --> 00:41:27
natural forward difference.
661
00:41:27 --> 00:41:30
So I replaced u' by
u_(n+1)-u_n over delta t.
662
00:41:30 --> 00:41:33
663
00:41:33 --> 00:41:38
And I replaced v by
the value I know.
664
00:41:38 --> 00:41:45
And do you know whose name
is associated with that?
665
00:41:45 --> 00:41:49
I'm almost going
to say his name.
666
00:41:49 --> 00:41:52
It's just I replace a
differential equation by a
667
00:41:52 --> 00:41:56
difference equation where each
step I just could figure out
668
00:41:56 --> 00:42:02
what the slope is and I go
along a straight line.
669
00:42:02 --> 00:42:05
Delta t further along.
670
00:42:05 --> 00:42:07
Do you know his name?
671
00:42:07 --> 00:42:08
Euler, Euler, right.
672
00:42:08 --> 00:42:10
It's Euler's method.
673
00:42:10 --> 00:42:14
The very first method
you would think of.
674
00:42:14 --> 00:42:17
Right?
675
00:42:17 --> 00:42:22
So I'm doing the same for v.
676
00:42:22 --> 00:42:27
This is my position
I've reached.
677
00:42:27 --> 00:42:30
So Euler's method is replace
this by v_(n+1)-v_n over
678
00:42:30 --> 00:42:30
delta t equals -u_n.
679
00:42:30 --> 00:42:38
680
00:42:38 --> 00:42:39
So this is Euler's idea.
681
00:42:39 --> 00:42:45
Predict the new value by a
forward difference starting
682
00:42:45 --> 00:42:50
from where the slope of
the line is this old one.
683
00:42:50 --> 00:42:55
That's the first difference
method you would think of.
684
00:42:55 --> 00:42:58
So now if I multiply up
by delta t I have a
685
00:42:58 --> 00:42:59
minus delta t*u_n.
686
00:42:59 --> 00:43:09
687
00:43:09 --> 00:43:14
Good time to pay attention.
688
00:43:14 --> 00:43:18
This is forward Euler.
689
00:43:18 --> 00:43:27
Forward Euler.
690
00:43:27 --> 00:43:32
The time step, you follow the
derivative, what the derivative
691
00:43:32 --> 00:43:34
is at the start of the step.
692
00:43:34 --> 00:43:36
You follow it for
a little step.
693
00:43:36 --> 00:43:39
Ah, let me draw what
the thing would do.
694
00:43:39 --> 00:43:42
Shall I draw what
forward Euler will do?
695
00:43:42 --> 00:43:43
So we're here.
696
00:43:43 --> 00:43:46
And what's our first step?
697
00:43:46 --> 00:43:47
Our first step.
698
00:43:47 --> 00:43:51
So we're starting at v
is zero at the start.
699
00:43:51 --> 00:43:53
So u won't change in one step.
700
00:43:53 --> 00:43:58
But u is one, so
v will go down.
701
00:43:58 --> 00:44:04
I think the first step
will take us there.
702
00:44:04 --> 00:44:08
That point was (1, 0).
703
00:44:08 --> 00:44:13
And I think the second step,
u came from the zero, so
704
00:44:13 --> 00:44:15
it'll still be a one.
705
00:44:15 --> 00:44:19
And v came from a zero,
but minus delta.
706
00:44:19 --> 00:44:27
I think it'll just
be 1-delta t.
707
00:44:27 --> 00:44:33
It's left the
circle, of course.
708
00:44:33 --> 00:44:42
And what do you think happens
when I do 1,000 steps?
709
00:44:42 --> 00:44:45
This is two lines in MATLAB.
710
00:44:45 --> 00:44:52
Just write those equations for
the next u coming from the
711
00:44:52 --> 00:44:57
previous u starting at (1,
0) and see what happens.
712
00:44:57 --> 00:45:00
And what do you think?
713
00:45:00 --> 00:45:03
It's going to sort
of go around.
714
00:45:03 --> 00:45:05
I mean, it's a reasonable,
Euler wasn't dumb.
715
00:45:05 --> 00:45:10
Anything true here is the fact
that Euler was not dumb.
716
00:45:10 --> 00:45:13
He had a bunch of ideas.
717
00:45:13 --> 00:45:16
This wasn't his greatest.
718
00:45:16 --> 00:45:22
But it's the starting point of
any, if you have to code some
719
00:45:22 --> 00:45:26
complicated problem, start
with Euler, forward Euler.
720
00:45:26 --> 00:45:28
What'll happen?
721
00:45:28 --> 00:45:30
It'll spiral out.
722
00:45:30 --> 00:45:30
Exactly.
723
00:45:30 --> 00:45:32
It'll spiral out.
724
00:45:32 --> 00:45:38
And if I was an eigenvalue
person, I guess which I
725
00:45:38 --> 00:45:45
probably am, I would look at
the, there's some matrix
726
00:45:45 --> 00:45:48
here that's finding the new
u, v from the old u, v.
727
00:45:49 --> 00:45:52
Maybe we can even find
out, write above
728
00:45:52 --> 00:45:57
what that matrix is.
729
00:45:57 --> 00:46:01
The new u, v come from the
old u, v, well there is
730
00:46:01 --> 00:46:02
the identity matrix.
731
00:46:02 --> 00:46:04
Here is a delta t.
732
00:46:04 --> 00:46:07
And here is our minus delta t.
733
00:46:07 --> 00:46:11
That would be my matrix.
734
00:46:11 --> 00:46:13
That's the forward
Euler matrix.
735
00:46:13 --> 00:46:19
The growth matrix for
forward Euler is that.
736
00:46:19 --> 00:46:24
Now what do you think about
its eigenvalues before?
737
00:46:24 --> 00:46:29
Spiral out was the
right answer.
738
00:46:29 --> 00:46:35
If delta t is small it'll stay
close to the circle, but
739
00:46:35 --> 00:46:38
it'll spiral out from it.
740
00:46:38 --> 00:46:45
So if we were following a
planet, it would just be off.
741
00:46:45 --> 00:46:48
So I can't immediately, well
you could almost see the
742
00:46:48 --> 00:46:51
eigenvalues of this matrix.
743
00:46:51 --> 00:46:54
What's the determinant
of that matrix?
744
00:46:54 --> 00:46:55
What's the determinant of G?
745
00:46:55 --> 00:46:58
I'm just asking you
to think through.
746
00:46:58 --> 00:47:00
You see a matrix like this.
747
00:47:00 --> 00:47:02
You wonder about its
eigenvalues so you
748
00:47:02 --> 00:47:03
take the determinant.
749
00:47:03 --> 00:47:06
And what do you get?
750
00:47:06 --> 00:47:10
One, something bigger than
one or less than one?
751
00:47:10 --> 00:47:11
Bigger.
752
00:47:11 --> 00:47:12
So that tells me what?
753
00:47:12 --> 00:47:16
Since the determinant is the
product of lambda_1*lambda_2,
754
00:47:16 --> 00:47:20
whatever those guys are, at
least one of those eigenvalues
755
00:47:20 --> 00:47:23
is bigger than one.
756
00:47:23 --> 00:47:27
This has an eigenvalue,
actually, the actual
757
00:47:27 --> 00:47:32
eigenvalues of this happen to
be one from the identity matrix
758
00:47:32 --> 00:47:38
and then this guy is what was
our example of an imaginary
759
00:47:38 --> 00:47:41
eigenvalue, i*delta t.
760
00:47:41 --> 00:47:45
Bigger than one in
absolute value.
761
00:47:45 --> 00:47:49
It goes out.
762
00:47:49 --> 00:47:57
So that's forward Euler
corresponding to that choice.
763
00:47:57 --> 00:48:00
Can I change to a second,
to Euler's other choice?
764
00:48:00 --> 00:48:03
And you can guess
what its name is.
765
00:48:03 --> 00:48:09
Instead of forward Euler
it's going to be?
766
00:48:09 --> 00:48:14
Backward Euler.
767
00:48:14 --> 00:48:17
Euler said, if you don't like
it, I've got another one.
768
00:48:17 --> 00:48:19
And what will happen then?
769
00:48:19 --> 00:48:25
Backward Euler is going to use,
instead of the old value, it
770
00:48:25 --> 00:48:29
will use the new values here.
771
00:48:29 --> 00:48:37
So these differences, you
see that this is now a
772
00:48:37 --> 00:48:39
backward difference?
773
00:48:39 --> 00:48:40
So I'm at time n+1.
774
00:48:41 --> 00:48:42
I'm at the new time.
775
00:48:42 --> 00:48:45
And this goes back from that.
776
00:48:45 --> 00:48:49
So this is at the new time and
this is a difference backward.
777
00:48:49 --> 00:48:52
And what happens now?
778
00:48:52 --> 00:48:55
Well, we could follow
backward Euler.
779
00:48:55 --> 00:48:57
We could follow its first step.
780
00:48:57 --> 00:48:59
What would its first
step be? u_n+u_1.
781
00:48:59 --> 00:49:02
782
00:49:02 --> 00:49:07
No, I can't follow that
first step so easily.
783
00:49:07 --> 00:49:09
Help me out by just
making a guess.
784
00:49:09 --> 00:49:14
What do you think
backward Euler does?
785
00:49:14 --> 00:49:16
It spirals in.
786
00:49:16 --> 00:49:23
Right, the cover of the book
has forward Euler, I think, is
787
00:49:23 --> 00:49:24
on the cover of the textbook.
788
00:49:24 --> 00:49:26
Spiraling out in the front.
789
00:49:26 --> 00:49:31
But maybe the back cover also
has Euler spiraling in.
790
00:49:31 --> 00:49:32
So here is backward Euler.
791
00:49:32 --> 00:49:37
And after I take a bunch of
steps it comes back somewhere
792
00:49:37 --> 00:49:40
there, but it's spiraled in.
793
00:49:40 --> 00:49:41
No good.
794
00:49:41 --> 00:49:47
Have we got one minute
for the good method?
795
00:49:47 --> 00:49:52
You're guessing that
it's this one?
796
00:49:52 --> 00:49:53
How does that work?
797
00:49:53 --> 00:49:58
Actually, it's pretty neat.
798
00:49:58 --> 00:50:03
So my question is, where
do I evaluate these guys?
799
00:50:03 --> 00:50:07
So this is like
forward and backward.
800
00:50:07 --> 00:50:14
The u equation, so I know
u_n and v_n, right?
801
00:50:14 --> 00:50:14
Everybody with me?
802
00:50:14 --> 00:50:18
I've got to time in, I know u_n
and v_n, the position and
803
00:50:18 --> 00:50:22
velocity approximately
at point n.
804
00:50:22 --> 00:50:23
And I want to go n+1.
805
00:50:23 --> 00:50:26
806
00:50:26 --> 00:50:27
what I know is v_n.
807
00:50:28 --> 00:50:31
So that's great.
808
00:50:31 --> 00:50:33
But what's going
to change here?
809
00:50:33 --> 00:50:39
Now I know u_(n+1) from
the first equation.
810
00:50:39 --> 00:50:42
What shall I do, where
shall I evaluate u in
811
00:50:42 --> 00:50:45
the second equation?
812
00:50:45 --> 00:50:45
At n+1.
813
00:50:46 --> 00:50:47
Why not?
814
00:50:47 --> 00:50:52
I know it already from the
first step, from the forward
815
00:50:52 --> 00:50:54
step with this equation.
816
00:50:54 --> 00:50:56
I've taken u forward.
817
00:50:56 --> 00:51:03
So I'll use that forward value
of u in the v equation.
818
00:51:03 --> 00:51:08
So it's forward n backward
at the same time.
819
00:51:08 --> 00:51:17
And the net result if I go back
from my system of two first
820
00:51:17 --> 00:51:20
order equations back to one
second order equation, I'll
821
00:51:20 --> 00:51:24
find the centered guy.
822
00:51:24 --> 00:51:34
And what do you think happens
with that center difference?
823
00:51:34 --> 00:51:39
I think MATLAB would
have to show it.
824
00:51:39 --> 00:51:43
And look in Section 2.2 for
the picture that shows
825
00:51:43 --> 00:51:44
the leapfrog method.
826
00:51:44 --> 00:51:49
So what do you think?
827
00:51:49 --> 00:51:57
I'm having to depend
on my memory now.
828
00:51:57 --> 00:52:01
It stays almost on the circle.
829
00:52:01 --> 00:52:05
I think it maybe, it stays
very close to the circle.
830
00:52:05 --> 00:52:11
And comes back almost at the
right time, but not exactly.
831
00:52:11 --> 00:52:13
Because we have an error here.
832
00:52:13 --> 00:52:15
No, we haven't got
the real equation.
833
00:52:15 --> 00:52:19
We've got a
difference equation.
834
00:52:19 --> 00:52:23
So maybe that's
the point to say.
835
00:52:23 --> 00:52:25
When you start with a
differential equation you've
836
00:52:25 --> 00:52:29
got various choices.
837
00:52:29 --> 00:52:32
Backward is actually safer.
838
00:52:32 --> 00:52:37
Probably for a car crash.
839
00:52:37 --> 00:52:42
No, actually I think for a car
crash they would use forward
840
00:52:42 --> 00:52:45
differences because they
have to take so many.
841
00:52:45 --> 00:52:47
Let me check on that.
842
00:52:47 --> 00:52:54
But for our problem, this
one, leapfrog is the winner.
843
00:52:54 --> 00:52:57
Actually, if anybody does
computation in computational
844
00:52:57 --> 00:53:00
chemistry, I mean that's a
subject that uses giant
845
00:53:00 --> 00:53:05
supercomputers to follow
molecular dynamics, to follow
846
00:53:05 --> 00:53:09
what molecules are doing at
high speed, very high speed.
847
00:53:09 --> 00:53:13
So this, there's a giant
number multiplying that u.
848
00:53:13 --> 00:53:19
And they use leapfrog method.
849
00:53:19 --> 00:53:23
I'll say a little more
Friday if I can about this
850
00:53:23 --> 00:53:26
general subject, because
it's so important.