1
00:00:03,000 --> 00:00:06,000
In 5.4, youll actually get your hands dirty in writing an ODE solver
2
00:00:06,000 --> 00:00:08,000
Youll write forward Euler and backward Euler
3
00:00:08,000 --> 00:00:13,000
Forward Euler is also very commonly called explicit Euler because its an explicit solver
4
00:00:13,000 --> 00:00:19,000
You will then use the explicit Euler solver you wrote to generate trajectories of the simple harmonic oscillator for different parameter sets
5
00:00:19,000 --> 00:00:22,000
Heres my Matlab implementation of forward Euler
6
00:00:22,000 --> 00:00:30,000
It takes in an initial condition x0, a step size h, a final time, and a showPlot parameter which I just need later on as a flag whether or not to plot the trajectories
7
00:00:30,000 --> 00:00:32,000
You can pretty much just ignore this one if you need to
8
00:00:32,000 --> 00:00:39,000
The first step is just to make xCurrent, which is a helper variable, equal to the initial condition, and make the x0 the first element in the trajectory
9
00:00:39,000 --> 00:00:44,000
Notice that, since Matlab starts array indexing at one, x0 will actually be the first element in the trajectory
10
00:00:44,000 --> 00:00:48,000
That is, x0 will be trajectory[1], so x1 will be trajectory[2]
11
00:00:48,000 --> 00:00:50,000
These off-by-one errors can often be very confusing
12
00:00:50,000 --> 00:00:54,000
We then just have this main loop
13
00:00:54,000 --> 00:01:02,000
All it does is take xCurrent, which is that helper variable I used above, starts from that position, and takes one step in the direction of the differential equation
14
00:01:02,000 --> 00:01:09,000
So the next element in the trajectory is simply the current position plus a tiny step in the direction of the ODE
15
00:01:09,000 --> 00:01:14,000
We then update xCurrent and loop for as many times as we need to get to the final time
16
00:01:14,000 --> 00:01:18,000
We then may or may not plot the trajectory, depending on that showPlot parameter
17
00:01:18,000 --> 00:01:24,000
The simple harmonic oscillator is given by these two differential equations, and has three parameters: k, m, and beta
18
00:01:24,000 --> 00:01:31,000
We then calculate x based on these differential equations evaluated at xCurrent and return x
19
00:01:31,000 --> 00:01:41,000
Notice that k = 2, m = 0.5, and beta = 0 are precisely the parameters that we used in both Quiz 5.3 and will use in Homework 5.4 for the first problem
20
00:01:41,000 --> 00:01:45,000
Its always a good idea when you write new software to desk check the software
21
00:01:45,000 --> 00:01:50,000
That means to do the first iteration by hand and make sure that what you got by hand is what you get with the software
22
00:01:50,000 --> 00:02:00,000
In Quiz 5.3, you used forward Euler to solve the simple harmonic oscillator for one time step using a time step of 0.1, starting at [-1, -2]
23
00:02:00,000 --> 00:02:08,000
And if we do that, we see that we get -1.2 and -1.6, which is precisely what we got doing the desk check
24
00:02:08,000 --> 00:02:09,000
So I would say that our software is good to go
25
00:02:09,000 --> 00:02:21,000
The first problem on Homework 5.4 after weve developed this software of course is to generate a trajectory for 0.5 seconds starting at k = 2, m = 0.5, and beta = 0 from the initial condition [-1, -2]
26
00:02:21,000 --> 00:02:22,000
Lets see what we get
27
00:02:22,000 --> 00:02:43,000
So, starting from initial condition [-1, -2], using a step size of 0.1, and going to time, or final time, 0.5, and in this case we want to not plot the final trajectory, we get that the last point in the trajectory, or the trajectory at time 0.5 using this solver, is position -1.5283 and velocity 0.6246
28
00:02:43,000 --> 00:02:46,000
That is the answer to question 1(a)
29
00:02:46,000 --> 00:02:50,000
For part b, we want to keep all the parameters the same, and the same initial condition
30
00:02:50,000 --> 00:02:58,000
We want to generate a 200-point trajectory instead of a 5-point trajectory, and well generate this trajectory using a time step of 0.1 as well as a time step of 0.11
31
00:02:58,000 --> 00:03:00,000
If we do this, we get the following plot
32
00:03:00,000 --> 00:03:06,000
In this plot, the blue curve is generated using a time step of 0.1, and the red curve is using a time step of 0.11
33
00:03:06,000 --> 00:03:10,000
As you can see, and as youd expect from the lecture, the red curve spirals out faster
34
00:03:10,000 --> 00:03:13,000
That is, if you take a bigger time step, youre going to step farther off the trajectory
35
00:03:13,000 --> 00:03:16,000
And in this case, youll spiral out faster
36
00:03:16,000 --> 00:03:19,000
See the lecture for a very nice visualization of this
37
00:03:19,000 --> 00:03:24,000
From this plot, we can see that using the delta t of 0.11 spiraled out faster
38
00:03:24,000 --> 00:03:27,000
Question 2 asks us to write backward Euler
39
00:03:27,000 --> 00:03:30,000
This is also called implicit Euler, as it is an implicit solver
40
00:03:30,000 --> 00:03:39,000
My implementation is a Matlab implementation that takes the exact same inputs as forward Euler did: namely, an initial condition, a step size, a final time, and this plot parameter
41
00:03:39,000 --> 00:03:43,000
My backward Euler code is identical to my forward Euler code except for this one key step
42
00:03:43,000 --> 00:03:48,000
Remember, with an implicit solver, you approximate the next step and then you use that as your current guess
43
00:03:48,000 --> 00:04:02,000
So in the case of backward Euler, this line is identical, except for, instead of evaluating the simple harmonic oscillator at xCurrent, you use xCurrent to do a forward Euler step, save that, and then you update the trajectory by using the slope at the forward Euler step
44
00:04:02,000 --> 00:04:10,000
This can be a little bit confusing, but youre approximating a step into the future, pulling the direction at that new step back in time, and then using the direction of the new step to go forward in time
45
00:04:10,000 --> 00:04:17,000
Liz provided a really nice visualization of whats happening here in the lecture, and if youre confused about this I encourage you to go back and see that
46
00:04:17,000 --> 00:04:25,000
While implicit solvers seem like they may be a little bit confusing, and you dont seem like youre gaining much, theyll become a lot more important when we do industrial-strength integrators in Unit 6
47
00:04:25,000 --> 00:04:35,000
Question 2, part a asks us to find what x at time 0.5 would be using a time step of 0.1, starting at [-1, -2] using our backwards Euler solver
48
00:04:35,000 --> 00:04:44,000
If we do this, we see that x at time 0.5 will be -1.2451 and the velocity will be 0.6136
49
00:04:44,000 --> 00:04:46,000
This is precisely the answer to question 2(a)
50
00:04:46,000 --> 00:04:51,000
Question 2(b) asks us to repeat problem 1(b), but instead using the backward Euler solver we just created
51
00:04:51,000 --> 00:04:58,000
Again, in this plot, the blue trajectory uses a step size of 0.1, and the red trajectory uses a step size of 0.11
52
00:04:58,000 --> 00:05:02,000
As we would expect, the red trajectory using a larger step size spirals in faster
53
00:05:02,000 --> 00:05:09,000
So this is the answer to question 2(b): the delta t of 0.11 spirals in faster than using a time step of 0.1
54
00:05:09,000 --> 00:05:19,000
Problem 2(c) asks us to use both a forward and backward Euler solver, using a time step of 0.1, to generate two 50-point trajectories, and plot them both together to see what the difference is
55
00:05:19,000 --> 00:05:24,000
Heres that plot: the blue uses backward Euler, and the red uses forward Euler
56
00:05:24,000 --> 00:05:27,000
Clearly the difference here is that one is spiraling in and one is spiraling out
57
00:05:27,000 --> 00:05:28,000
This should not be surprising from what we just did
58
00:05:28,000 --> 00:05:36,000
Whats interesting here is that this is an undamped simple harmonic oscillator, yet one of these is growing in amplitude, and one of these is shrinking in amplitude
59
00:05:36,000 --> 00:05:39,000
So the numerical solver is introducing dynamics that dont actually exist
60
00:05:39,000 --> 00:05:44,000
But all you need to know for this problem is that the difference between these two is that one is spiraling in and one is spiraling out
61
00:05:44,000 --> 00:05:46,000
This is the answer to question 2(c)