ORDINARY DIFFERENTIAL EQUATIONS POWERTOOL
Lesson 13 -- Direction Fields and Phase Portraits
Prof. Douglas B. Meade
Industrial Mathematics Institute
Department of Mathematics
University of South Carolina
Columbia, SC 29208
URL: http://www.math.sc.edu/~meade/
E-mail: meade@math.sc.edu
Copyright 2001 by Douglas B. Meade
All rights reserved
-------------------------------------------------------------------
Outline of Lesson 13
13.A Direction Fields and Phase Portraits
13.B Equilibrium Points
13.C Phase Portraits for 2d Linear Autonomous Systems
13.C-1 Example 1 (Saddle)
13.C-2 Example 2 (Center)
13.C-3 Example 3 (Node)
13.C-4 Example 4 (Spiral)
13.C-5 Eigenvalues and the Phase Portrait
13.C-6 Representative Examples
13.D Phase Portraits for 2d Nonlinear Autonomous Systems
13.D-1 Example 5
13.D-2 Example 6 (The Nonlinear Pendulum)
13.D-3 Example 7
13.E Nonautonomous Systems
13.E-1 Example 8
13.F Stability of an Equilibrium Point
13.F-1 Example 9
13.F-2 Example 10
13.F-3 Example 11
13.F-4 Summary
13.G Linearization of Nonlinear Systems
13.G-1 Tangent-Line Approximation
13.G-2 Tangent-Plane Approximation
13.G-3 Linearizing Autonomous Systems
13.G-3A Example 12
13.G-4 Interpretation Rules
13.G-5 Example 13
Initialization
Warning, the name changecoords has been redefined
Warning, the assigned names arrow and translate now have a global binding
Warning, these names have been rebound: `&x`, ArcLength, CrossProduct, DotProduct, Tangent
Warning, the assigned names `<,>` and `<|>` now have a global binding
Warning, these protected names have been redefined and unprotected: `*`, `+`, `.`, D, Vector, diff, int, limit, series
A direction field for a two-dimensional system of first-order ODEs, drawn in the phase plane for the system, is similar to the direction field for a single first-order ODE (see Lesson 1, Lesson 3, or Lesson 11).
For the autonomous system, the independent variable does not explicitly appear, and therefore has the representation
The direction field, a figure in the phase plane (which is the -plane), contains arrows whose direction at is determined by the vector . The slope of the shaft of the arrow is given by
at any point for which . If , the direction field will have an arrow pointing either straight up or straight down. Since the direction field shows only the direction in which a trajectory moves, the length of the arrows in unimportant.
A graph of solution curves in the space of the dependent variables is a phase portrait for the system. (See Lesson 1 and Lesson 11.)
An equilibrium point in the phase plane is a point at which = 0. Such a point represents an equilibrium solution, a constant solution where the system point sits and does not travel along an orbit.
For the autonomous linear system, we find the only equilibrium point to be the origin, , provided the system matrix is not singular. In the two-dimensional case, there are four kinds of equilibrium points, and they have the names saddle, center, node, and spiral. Table 13.1 summarizes these names and briefly characterizes each type of equilibrium point.
=============================================
SADDLE: Two paths enter, two leave, all others pass by
CENTER: Orbits are closed paths around the equilibrium point
NODE: All paths either enter or leave the equilibrium point
SPIRAL: Each path spirals around the equilibrium point.
Table 13.1
Section 13.C contains examples of the phase portraits for each type of equilibrium point listed in Table 13.1
Corresponding to the four types of equilibrium points listed in Table 13.1, there are characteristic phase portraits. The geometric properties of these distinctive phase portraits are related to the analytic behaviors of the solutions of the corresponding systems.
The following four examples illustrate each of these characteristic phase portraits.
Consider the system of differential equations
for which solutions through the eight initial points
will be obtained. Something as repetitive as generating solutions for these eight points demands a for-loop structure. Each IVP will be solved, and the solution , packaged in a list of the form
ready for parametric plotting.
Next, plot all eight trajectories on a single set of axes, forming Figure 13.1.
Adding some arrows and labels by executing the following Maple commands, gives Figure 13.2.
The directions on the trajectories correspond to the motion of the system point as the parameter increases along the curve. These directions are obtained from the differential equations. Table 13.2 indicates how directions are obtained for the paths through the four points labeled in Figure 13.2.
Table 13.2
The origin is a saddle for the differential equations in this example. For the saddle, two opposing trajectories enter the equilibrium point, and two opposing trajectories leave it. One such trajectory is called a separatrix (plural, separatrices). Trajectories starting on one side of a separatrix remain on that same side.
Aside from the two separatrices that enter the equilibrium point, all other solutions move away from it. Two trajectories appear to exit the equilibrium point, and the remaining trajectories first move towards it, then ultimately move away.
The equilibrium point at the origin is a center for the differential equations
A set of initial points is
The initial points have three values, the first, representing the initial time. Here, that is , the time at which the system point is at that initial point whose and coordinates are the second and third value in the initial condition.
The resulting phase portrait, Figure 13.3, is
A center is surrounded by closed trajectories.
The following linear system of differential equations has an equilibrium point at the origin, and it is a node.
With the following initial points,
the behavior of solutions near the node can be determined.
The phase portrait, Figure 13.4, is drawn via
The equilibrium point at the origin is a node "out," sometimes called a repeller or source. Were the trajectories in this node to point inwards, towards the origin, the node would be a node "in," sometimes called an attractor or sink.
The equilibrium point at the origin is a spiral for the differential equations
and the phase portrait, Figure 13.5, is
The equilibrium point at the origin is an inward spiral. Like the inward node, the inward spiral is sometimes called an attractor also. Clearly, if the trajectories spiral outward, the phase portrait would be that of an outward spiral, and it would likewise be called a repeller.
For the linear system x' = Ax, the eigenvalues of the matrix A characterize the nature of the phase portrait at the origin. These relationships are summarized in the Table 13.3.
Equilibrium Point Eigenvalues
============ =============
SADDLE: Opposite signs
CENTER: Pure imaginary
NODE: Same sign
SPIRAL: Complex conjugates
===================================
Table 13.3
Table 13.4 below gives representative examples of the four types of equilibrium points, and their eigenvalues.
We distinguish between a proper node and an improper node. For the moment, the improper node is characterized by repeated eigenvalues. In the next lesson we will see that the improper node poses a minor difficulty when studying nonlinear systems.
Nodes and saddles have exponential solutions that tend either to zero or infinity. Hence, trajectories for these systems will either tend towards or away from the origin, depending on the sign of the eigenvalues. The center has trigonometric solutions that are the parametric representations of closed curves. The spirals have solutions with the trigonometric terms of the center, but have exponential factors that force the solutions either towards or away from the origin.
Equilibrium point Representative Eigenvalues Corresponding Solutions
===================================================================
Proper Node In: -1,-2 =>
Improper Node In: -1,-1 =>
Proper Node Out: 1,2 =>
Improper Node Out: 1,1 =>
Saddle: -1,2 =>
Center: =>
Spiral In: =>
Spiral Out: =>
=====================================================================
Table 13.4
The two-dimensional nonlinear autonomous system
can have more than one equilibrium point. Each real solution of the simultaneous equations = is an equilibrium point, so the phase portrait for such a system is much more complicated than that of the autonomous linear system. This additional complexity is illustrated in the following three examples.
A direction field for the nonlinear autonomous system of differential equations
appears in Figure 13.6.
The DEplot command was used to create this figure. The key differences in usage for this command are the specification of the scene and the dimensions of the viewing window; the "time" interval is required, but relevant only when solution curves are included in the plot:
This system appears to have several equilibrium solutions, namely,
which are easier to read as the list of points
Figure 13.7 shows these points as blue circles on the direction field of Figure 13.6.
Figure 13.7 suggests that the equilibrium solution at the origin is a source, the one at is a sink, and the ones at and are saddle points.
To create a phase portrait for this system it is necessary to specify an initial condition for each trajectory. To create a reasonable set of initial conditions quickly, start with the equilibrium points and a few additional points in the phase plane, thereby developing the list
Perturb these points a small amount in each component to obtain the final list
After converting each initial condition into a pair of equations, create the phase portrait, and display it, along with the direction field and equilibrium solutions, as Figure 13.8.
The equations
formulate Newton's Second Law of Motion for a simple plane pendulum in which is the angle between the stable hanging position and the pendulum's position at time ; is the angular velocity of the pendulum bob. Thus, for every integer , the point corresponds to the stable hanging position and the point corresponds to the unstable position with the pendulum pointing straight up from the pivot.
For this system, the direction field in Figure 13.9 appears better with the dirgrid option used to increase the resolution of the direction arrows, and with the constrained option used to impose equal scaling on both axes.
Although solve finds only one equilibrium solution for this system
it is readily apparent that there is an infinite collection of equilibrium solutions. Maple can find the entire family of equilibria when the _EnvAllSolutions environment variable is set to true. The result is then
(Note on syntax: Instead of the parameter k, Maple returns a symbol such as _Z1 to denote an arbitrary integer. Unfortunately, on subsequent executions of the solve command, this symbol will change sequentially to _Z2, _Z3, . To avoid such complications, the operation of solving for the equilibrium solutions has been coded to return the symbol for the arbitrary integer.)
To plot the equilibrium solutions in the direction field, convert the general expression to a point
and substitute appropriate integers to get the list of points
The desired plot then appears in Figure 13.10.
The direction field suggests different behaviors for the equilibria at even and odd integer multiples of . Solution curves will help to further elaborate on these differences. With initial conditions selected along the -axis, and from the "inflow" portions of the window for the direction field, i.e., from and , , the totality of initial points is
The combined direction field and phase portrait so generated is seen in Figure 13.11.
Thus, for any integer , the point is a center and is a saddle point.
If damping is added to the simple plane pendulum of Example 6, the system equations become
To facilitate comparisons with Example 6, create the direction field with the same window, but with a slightly higher resolution. The result appears in Figure 13.12.
The differences in the direction fields are obvious. Nonetheless, the system has the same equilibrium solutions, namely,
where the integer has the same explanation as appeared in Example 6.
Again writing this solution as the point
and generating the following list of equilibrium points,
Figure 13.13 showing both the equilibrium points and the direction field, is generated.
For the phase portrait, select initial conditions from the top and bottom edges of the viewing window. With the initial points
the phase portrait in Figure 13.14 is generated.
From the direction field and equilibrium solution, it appears as though the equilibrium solutions at that were centers in the undamped system are inward nodes (attractors, or sinks) in the damped system while the equilibria at are saddle points for both systems.
Any nonautonomous system can be reformulated as an autonomous system via the introduction of the independent variable, say , and augmenting the system with the (trivial) differential equation
The strange appearance of this differential equation that results from using the same name for the independent variable and one of the unknown functions can be avoided by introducing a new name for the independent variable, say . Then all occurrences of in the original system are replaced with and the auxiliary ODE becomes
The nonautonomous system
can be converted to an autonomous system with the help of the dchange command from the PDEtools package. This results in the three-dimensional system
Because this is a three-dimensional system, it is not realistic to draw a direction field. Instead, a sample of solution curves with initial conditions
obtained with the DEplot3d command, appears in Figure 13.15.
This plot can be rotated in real time by using the mouse to "grab" and move a corner of the box.
Specific views of the solution can be obtained by explicitly listing the viewing coordinates. For example, the x-y view is obtained with
Note that the crossing of solutions in this view is not reason for concern -- this is not an autonomous system in x and y -- and the solution curves do not intersect in the x-y-t space (the phase space for the augmented system).
The - view is obtained with
and the t-y view is obtained with
As a concluding note, it should be mentioned that the DEplot3d command is capable of handling nonautonomous systems. A simpler means of plotting the solution trajectories in this example would be to use
An isolated equilibrium point for an autonomous system of differential equations is said to be
Orbitally (Neutrally) Stable if all solutions sufficiently close to the equilibrium point remain close to that equilibrium point as ,
Asymptotically Stable if it is orbitally stable, and in addition, all solutions sufficiently close to the equilibrium point actually tend towards that equilibrium point as ,
Unstable if at least one solution starting close to the equilibrium point does not remain close to that equilibrium point as .
Special note
Some texts drop the modifier in conjunction with the orbitally (neutrally) stable equilibrium point, choosing to call such a point a stable equilibrium point. However, without insisting on some modifier in this case, it can be difficult to determine if the word "stable" is being used in a precise sense. Students who intend the term "asymptotic stability" but who carelessly drop the modifier then end up using the term "stable" ambiguously. They could have intended to use "stable" in the sense of orbital stability, or in the sense of asymptotic stability, improperly expressed. Hence, to avoid this possible ambiguity, we will always modify the term "stable" with either asymptotic or with orbital.
Table 13.5 listing Examples 9 - 11, which will be used to illustrate the terms given above.
Example System Matrix A Eigenvalues Phase Portrait
================= ============ ========= ==========
(13.9)
+ Figure 13.16
___________________________________________________________________________
(13.10)
Figure 13.17
(13.11)
+ 1 Figure 13.18
____________________________________________________________________________
Table 13.5
The system in Example 9, Table 13.5, with equations
has an orbitally stable equilibrium point at the origin. To see what this means, form the coefficient matrix, A,
and compute its eigenvalues.
The eigenvalues are a pure imaginary conjugate pair. Hence, the phase portrait is that of the center. The origin is the center of closed orbits that neither spread apart nor approach each other, the essential meaning of orbital stability. Any solution starting near the origin remains near the origin. It does not necessarily get arbitratily close to the origin, but it does not get far from the origin either. The following phase portrait, Figure 13.16, is useful in this regard since is shows that solutions starting near the origin remain near the origin, tending neither towards, nor away from the origin.
To draw the figure, pick the initial points
and use Maple's DEplot command to obtain a phase portrait.
The system in Example 10, Table 13.5, whose equations are
has an asymptotically stable equilibrium point at the origin. To see what this means, form the coefficient matrix, A
and compute the eigenvalues of A.
The eigenvalues of the system matrix are of the same sign, and negative. Hence, the phase portrait is that of the node in. That makes the origin an attractor, and hence, all solutions will be drawn into the origin, making the system asymptotically stable. Any solution that starts near the origin will tend towards the origin as . For example, consider the solution starting at ( ), obtained via
Extract the component solutions
and compute the limits
As promised, an arbitrary solution starting near the origin, not only remains near the origin, but tends to the origin as . Hence, all solutions, not only the ones near the origin, tend towards the origin with increasing , as seen in the phase portrait, Figure 13.17. To obtain this figure, pick the initial points
and use Maple's DEplot command to obtain
The system in Example 11, Table 13.5, whose equations are
has an unstable equilibrium point at the origin. To see what this means, form the coefficient matrix, A,
The eigenvalues are of opposite signs. Hence, the phase portrait is that of the saddle. There are two trajectories that will enter the origin, but all other solutions, no matter how close to the origin they begin, will tend towards infinity as increases. Since at least one solution starting near the origin does not remain near the origin, the origin is an unstable equilibrium point.
Solutions starting on the y-axis will tend towards the origin. In fact, such solutions are precisely
for which
Any other solution with , given by
will not remain near the origin, as shown by the following limits. Indeed, Maple computes the limits
for , whereas if , Maple computes
As promised, all other solutions, no matter how near the origin they may have started, tend towards infinity.
As seen in the following figure, Figure 13.18, the phase portrait, that of a saddle, makes the instability of this equilibrium clear.
To obtain this figure, pick the following initial points
and use the DEplot command.
Solutions starting in the half-plane tend to , whereas solutions starting in the half-plane tend to .
An equilibrium point that is not stable, is called unstable. Thus, if at least one solution starting near the equilibrium point does not remain near it, the equilibrium point is unstable.
Table 13.6 summarizes the relationships between phase portraits, eigenvalues, and stability.
Equilibrium Point Eigenvalues Conditions Stability
=========== ======== ============ ===============
Node (in) real and negative Asymptotically stable
Spiral (in) + Asymptotically stable
Center + pure imaginary Orbitally stable
Node (out) real and positive Unstable
Spiral (out) + Unstable
Saddle real, with opposite signs Unstable
==============================================================================
Table 13.6
Near the point of contact, the curve is approximated by its tangent line. For example, the function
and its tangent at the point appear in Figure 13.19.
Near the point the function very much resembles its tangent line . Analytically, the tangent line is the Taylor polynomial of degree 1, .
The tangent line can be constructed by the following calculations. First, the slope at ,
then, the equation of the tangent line via the point-slope formula
Analytically, the tangent line is a Taylor polynomial of degree 1, obtained in Maple via
Replacing a curve with its tangent-line approximation is known as linearizing in one dimension.
The tangent plane approximates a function of two variables. This is the two-dimensional version of linearization.
For example, at the point where , the tangent plane approximating the function
is . This first-degree polynomial approximation of is just the multivariable Taylor polynomial of degree one.
To construct this tangent plane, first find the corresponding z-coordinate at , the function value . Thus,
To construct a plane tangent at obtain a vector normal to the surface . The general normal vector is
so compute
and evaluate it at to obtain
Once the normal for the tangent plane is known, write the equation of the plane as
To determine the value of d, force the plane to go through the point . Hence,
from which the implicit form of the tangent plane is found to be
Explicitly solving this equation for yields
That is one way to find the equation of the tangent plane. Alternatively, at the point where , generate the multivariable Taylor polynomial of degree one. In Maple, use the mtaylor command to create a linear approximation to the function , thereby getting the explicit form of the equation of the tangent plane, namely, . This results in
The two computations of the equation of the tangent plane yield the same result. Hence, the first-degree Taylor polynomial gives the linear approximation known geometrically as the tangent plane.
Near an equilibrium point of the autonomous nonlinear system
the behavior of solutions can be studied by linearizing the functions and at each equilibrium point.
Let ( ) be an equilibrium point determined by the equations = 0. Then, lumping the higher-order remainder terms into the functions and , the first-degree Taylor expansions
become
because is an equilibrium point at which = 0. The nonlinear autonomous system now reads
Define the change of variables
so that
and the system now reads
= +
or
u' = Au + P
where u = , the matrix A is the Jacobian Matrix
and the functions and are just the functions and under the change of variables .
The equilibrium point ( ) in the xy-plane has been translated to the equilibrium point in the uv-plane. If the perturbation term P is dropped, the resulting system u' = Au is a linear system with an equilibrium point at . The behavior of the nonlinear xy-system at ( ) can often be deduced by analyzing the behavior of the linear uv-system at , as will be shown shortly.
The differential equations
model the dynamics of two species competing for the same resources. Such models will be studied at length in Lesson 27. For this model, the functions and are
The equilibrium points are obtained by solving the equations = 0, done in Maple via
Extract and label these four equilibrium points as
The Jacobian matrix is
which is calculate using Maple's built-in jacobian command.
Unfortunately, the Jacobian is defined mathematically as the determinant of the Jacobian matrix, but Maple returns just the matrix itself with the command
At each of the equilibrium points, the Jacobian matrix evaluates to
For each of the four matrices obtained by evaluating A at an equilibrium point, obtain the eigenvalues. Thus, compute the eigenvalues of each of the four matrices .
To determine the signs of the eigenvalues from the matrix , obtain floating point equivalents of its exact eigenvalues.
For the linearized systems we have
P1 = (0, 0) - node out
P2 = (20, 0) - saddle
P3 = (0, 14) - saddle
P4 = (12, 6) - node in
Rules are needed for deducing the nature of the equilibrium points of the nonlinear system from the information generated for the linear system. These rules appear in the next section.
The following five rules summarize how to convert information about the linearized system into information about the nonlinear system. They form a hierarchy, and are applied from top to bottom.
1. Equilibrium points that are unstable or asymptotically stable in the linearized system will be the same in the nonlinear system.
2. Saddles and spirals in the linearized system remain the same in the nonlinear system.
3. Proper nodes ( 's real and distinct) in the linearized system remain the same in the nonlinear system.
4. Improper nodes ( 's real and equal) in the linearized system can be nodes or spirals in the nonlinear system.
Such nodes are called improper because of the slight ambiguity in passing from the linearized to the nonlinear system. However, instability and asymptotic stability will persist by Rule 1. Thus, an improper outward node will be either an outward node or an outward spiral. An improper inward node will be either an inward node or an inward spiral. So even though there is a slight ambiguity, it is still possible to determine whether the nonlinear system is unstable or asymptotically stable.
5. Centers in the linearized system can be centers or spirals in the nonlinear system
This is the difficult case. Since centers are not asymptotically stable to begin with, the center that becomes a spiral can become an inward or an outward spiral. Centers are neither unstable nor asymptotically stable, so Rule 1 does not apply. Thus, the very nature of the behavior of the nonlinear system is in complete doubt if the linearized system has a center. No decision can be made on the basis of the linearization and far more powerful tools are then necessary to determine the behavior of the nonlinear system at such an equilibrium point.
It is often exceedingly difficult to distinguish between a spiral and a stable center in a nonlinear system. It is rarely clear if the spiral behavior is an artifact of numerical (roundoff) error in an approximating numerical scheme. Analyzing the stability of such a point requires a great deal more mathematical skill than will be developed here.
To analyze the stability of the nonlinear autonomous system
obtain the equilibrium points and , and the Jacobian matrix
Evaluate A at each equilibrium point, and compute the corresponding eigenvalues, obtaining the matrices
and
and the respective eigenvalues , and + .
The point is a saddle, and unstable, for the linearized system; therefore, it is a saddle, and unstable, for the nonlinear system.
The point is a center for the linearized system; therefore, the behavior for the nonlinear system cannot be determined from the linearization. The best that can be done with the present suite of tools is try to infer the behavior from the phase portrait in Figure 13.20, obtained below. The phase portrait suggests that may well be a center in the nonlinear case.
To implement these calculations in Maple, enter the funtions and .
Then, solve the equations x' = y' = 0 for equilibrium points, obtaining
There are two equilibrium points that can now be carefully identified and labeled.
Compute the generic Jacobian Matrix of first partial derivatives.
Evaluate the Jacobian Matrix at each of the two equilibrium points.
Obtain the eigenvalues of each of the matrices .
Point P1 = (0, 0) is a saddle for the linearized system; therefore it is a saddle for the nonlinear system.
Point P2 = (4, 1) is a center for the linearized system; therefore the behavior for the nonlinear system cannot be determined from the linearization. The best that can be done with the present suite of tools is to draw a phase portrait, hoping to infer the behavior from that plot.
For this, write the two differential equations.
and select the initial conditions
Then, use the DEplot command to construct a phase portrait.
The phase portrait suggests that (4, 1) may indeed be a center in the nonlinear case.
[Back to ODE Powertool Table of Contents]
