This laboratory will investigate two topics:
We have seen initial-value problems where we have that the first-derivative of one function of one variable is described by a function of time and the value of the function:
How do we deal with a differential equation with a higher order, or a system of coupled differential equations.
So as to not alarm you, the beauty is that the code you have already written will be reused with only three small changes to three line of code.
Suppose you have rabbits: left to themselves, they will grow exponentially (assuming adequate food), but if you introduce foxes into the equation, the more foxes, the two will play one off another: an increase in rabbits will result in an increase in foxes, but at some point, the foxes will overpopulate and decimate the rabbits, at which point, the foxes die out, once again allowing the rabbits to begin increasing exponentially.
Suppose that Consider the following system of differential equations:
Here represents the density of rabbit families per acre while represents the density of fox families per acre. If you begin with the density of rabbits of 0.5 families per acre and a density of foxes of 0.3 families per acre, the relatively low number of foxes allows the rabbit population to grow exponentially and the equally relative low number of rabbits causes the foxes to begin to die off. However, as the population of rabbits increases, the population of foxes rebounds until they begin to decimate the rabbits. This cycle is shown in Figure 1.
Figure 1. The population densities of foxes and rabbits over time.
To understand this, we note that at time , , and therefore the population density of rabbits is increasing. The fox population density, however, is decreasing: . Using Euler's method, we could therefore predict the population at time : and . We could carry on and now predict and .
The differential equation described above falls into a category of differential equations called predator-prey equations (see the wikipedia article on the Lotka-Volterra equation). There is no analytic solution for these equations and they must be solved numerically.
Another way of looking at the interaction of these two species is to plot versus over time. Such a plot is shown in Figure 2.
Figure 2. An animation of the population densities of foxes and rabbits over time.
Again, when both population densities are relatively small, the number of rabbits increases, but this increase is followed by an increase in foxes which in turn decreases the rabbit population density. Without a food supply, the foxes themselves begin to die off allowing the rabbits a chance to recover.
Another interesting differential equation are the Lorenz equations. You will have heard of the butterfly effect and this description is, as we will see, a direct response to the solution to this system of three differential equations:
with, for example, the initial values , , and .
The solutions to these coupled initial-value problems must be solved numerically. The solutions the these three are shown in Figure 3.
Again, like before, we can make a plot of the three points in space and watch how the three change over time. This is shown in Figure 4 which gives us the familiar butterfly.
This is reasonably straight-forward to implement: Modify dp45 so that it accepts a column vector as an initial condition. This requires three changes:
Now, create the function
function [dy] = f10a( t, y ) dy = [y(1)*(1 - 2*y(2)), -y(2)*(1 - y(1))]'; end
With a little bit of thought, you should notice that this is the predator-prey model described above. If you wanted, you could write it as:
function [dy] = f10a( t, y ) r = y(1); % rabbits f = y(2); % foxes dy = [r*(1 - 2*f), -f*(1 - r)]'; end
Next, execute the following to give plots similar to Figure 1 and Figure 2:
hold off [ts, pts] = dp45( @f10a, [0, 6], [0.5 0.3]', 0.1, 1e-8 ); plot( pts(2,:), pts(1,:) ) % rabbits on the vertical axis, foxes on the horizontal plot( ts, pts(1,:) ) % plot population density of rabbits over time hold on plot( ts, pts(2,:) ) % plot population density of foxes over time
Similarly, write a function [dv] = f10c( t, v ) which, in a similar manner, implements the three coupled Lorenz equations. In this case, v(1), v(2) and v(3) will represent , and , respectively. When you run the following code, you should get an output similar to Figure 5.
[ts, pts] = dp45( @f10c, [0, 100], [-5 0 28]', 0.1, 1e-5 ); size( pts ) ans = 3 34145 plot3( pts(1,:), pts(2,:), pts(3,:) )
Figure 5. The output in Matlab of the given commands.
This is your first use of the plot3 command in Matlab for plotting points in three dimensions. You will notice a few entries in the toolbar, highlighted in Figure 6.
Figure 6. Five tools for 3-dimensional plots: Zoom In, Zoom Out, Pan, Rotate 3D, Data Cursor, and Brush/Select Data.
These five tools do, in order:
Assignment: use the pan tool to move the image around, then rotate the image to change the angle, and finally select the Zoom In tool and click twice on one region of the plot to zoom in on the plot. Save this image and print it. Each student must do their own plot.
Given a real-valued function of a real variable , suppose we want to find a root of the function, that is, a point such that .
We will look at three techniques for approximating
You may have already seen Newton's method.
First, if is a root, then :
Next, assuming that the derivative at the point is not zero, that is, , we can divide the equation by that derivative:
The next step is to isolate the error term for on one side and group the other:
The expression in the parentheses is the formula for Newton's method, and thus, this equals :
That is, the error of the approximation is a constant multiplied by the square of the error of the approximation. Thus, the claim is, if the error is at one step, the error will be on the order of with the next step, and on the order of with the next.
Now, assuming that we are close to the root and the second derivative is continuous, then it is a reasonably good approximation that
and we will use this to demonstrate how well Newton's method works.
Consider the function and let us start with .
|Absolute Error||Estimated Absolute Error|