ORDINARY DIFFERENTIAL EQUATIONS POWERTOOLLesson 13 -- Direction Fields and Phase PortraitsProf. Douglas B. MeadeIndustrial Mathematics InstituteDepartment of MathematicsUniversity of South CarolinaColumbia, SC 29208
URL: http://www.math.sc.edu/~meade/E-mail: meade@math.sc.eduCopyright \251 2001 by Douglas B. MeadeAll rights reserved-------------------------------------------------------------------<Text-field layout="Heading 1" style="Heading 1">Outline of Lesson 13</Text-field>13.A Direction Fields and Phase Portraits13.B Equilibrium Points13.C Phase Portraits for 2d Linear Autonomous Systems13.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 Portrait13.C-6 Representative Examples13.D Phase Portraits for 2d Nonlinear Autonomous Systems13.D-1 Example 513.D-2 Example 6 (The Nonlinear Pendulum)13.D-3 Example 713.E Nonautonomous Systems13.E-1 Example 813.F Stability of an Equilibrium Point13.F-1 Example 913.F-2 Example 1013.F-3 Example 1113.F-4 Summary13.G Linearization of Nonlinear Systems13.G-1 Tangent-Line Approximation13.G-2 Tangent-Plane Approximation13.G-3 Linearizing Autonomous Systems13.G-3A Example 1213.G-4 Interpretation Rules13.G-5 Example 13<Text-field layout="Heading 1" style="Heading 1">Initialization</Text-field>restart;with( DEtools ):
with( plots ):
with( plottools ):
with( LinearAlgebra ):
with( PDEtools ):
with( Student[Calculus1] ):
with( VectorCalculus ):<Text-field bookmark="13.A" layout="Heading 1" style="Heading 1">13.A Direction Fields and Phase Portraits</Text-field>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 representationNiMvKiYlI2R4RyIiIiUjZHRHISIiLSUiZkc2JCUieEclInlHNiMvKiYlI2R5RyIiIiUjZHRHISIiLSUiZ0c2JCUieEclInlHThe direction field, a figure in the phase plane (which is the NiMlI3h5Rw==-plane), contains arrows whose direction atNiMtJSFHNiQlInhHJSJ5Rw== is determined by the vector NiMtJSdtYXRyaXhHNiM3JDcjLSUiZkc2JCUieEclInlHNyMtJSJnR0Yq. The slope of the shaft of the arrow is given byNiMvKiYlI2R5RyIiIiUjZHhHISIiKiYtJSJnRzYkJSJ4RyUieUdGJi0lImZHRixGKA==at any point for which NiMwLSUiZkc2JCUieEclInlHIiIh. If NiMvLSUiZkc2JCUieEclInlHIiIh, 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.)<Text-field bookmark="13.B" layout="Heading 1" style="Heading 1">13.B Equilibrium Points</Text-field>An equilibrium point in the phase plane is a point at which NiMvLSUjeCdHNiMlInRHLSUjeSdHRiY= = 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, NiMtJSFHNiQiIiFGJg==, 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 byCENTER: Orbits are closed paths around the equilibrium pointNODE: All paths either enter or leave the equilibrium pointSPIRAL: Each path spirals around the equilibrium point.============================================= Table 13.1Section 13.C contains examples of the phase portraits for each type of equilibrium point listed in Table 13.1<Text-field bookmark="13.C" layout="Heading 1" style="Heading 1">13.C Phase Portraits for 2d Linear Autonomous Systems</Text-field>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.<Text-field bookmark="13.C-1" layout="Heading 2" style="Heading 2">13.C-1 Example 1 (Saddle)</Text-field>Consider the system of differential equationsq1 := diff(x(t),t) = x(t) + 2*y(t);
q2 := diff(y(t),t) = 4*x(t) - y(t);for which solutions through the eight initial pointsP1 := [1,1];
P2 := [-1,-1];
P3 := [5,-10];
P4 := [-5,10];
P5 := [6,-10];
P6 := [4,-10];
P7 := [-6,10];
P8 := [-4,10];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 NiQtJSJ4RzYjJSJ0Ry0lInlHRiU=, packaged in a list of the formNiM3JS0lInhHNiMlInRHLSUieUdGJi9GJzsiIiEiIiY=ready for parametric plotting.for k from 1 to 8 do
init := x(0)=P||k[1], y(0)=P||k[2];
qq := dsolve({q1,q2, init},{x(t),y(t)});
f||k := eval([x(t),y(t),t=0..5], qq);
print(f||k);
od:Next, plot all eight trajectories on a single set of axes, forming Figure 13.1.pp1 := plot([f||(1..8)], x = -10..10, y = -10..10, color=black, tickmarks=[[-10,-5,5,10],[-10,-5,5,10]], labels=[x,`y `], labelfont=[TIMES,ITALIC,12]):
display(pp1, title="Figure 13.1");Adding some arrows and labels by executing the following Maple commands, gives Figure 13.2.pp2 := textplot({[3.5,.6,`A`], [.42,5.9,`B`], [-3.7,-.8,`C`], [ .6,-4.7,`D`]}):
for k from 1 to 8 do
pp||(k+2) := disk(P||k,.25,color=black):
od:
pp11 := textplot({[1.5,.7,`1`],[-.8,-2,`2`], [4,-11,`6`], [5,-11,`3`], [6,-11,`5`], [-6,10.6,`7`], [-5,10.6,`4`], [-4,10.6,`8`]}):
pp12 := arrow([5,5],vector([.4,.4]),.2,.4,1,color=black):
pp13 := arrow([-5,-5],vector([-.4,-.4]),.2,.4,1,color=black):
pp14 := arrow([2.5,-5],vector([-.2,.4]),.2,.4,1,color=black):
pp15 := arrow([-2.5,5],vector([.2,-.4]),.2,.4,1,color=black):
pp16 := arrow([-9.408, -8.265],vector([-.4,-.4]),.2,.4,1,color=black):
pp17 := arrow([-9.4, -10.4],vector([-.4,-.4]),.2,.4,1,color=black):
pp18 := arrow([9.408, 8.265],vector([.4,.4]),.2,.4,1,color=black):
pp19 := arrow([9.4, 10.4],vector([.4,.4]),.2,.4,1,color=black):
display([pp||(1..19)],view=[-10..10,-11..11], scaling=constrained, labels=[x,y],labelfont=[TIMES,ITALIC,12], title="Figure 13.2");The directions on the trajectories correspond to the motion of the system point NiMtJSFHNiQtJSJ4RzYjJSJ0Ry0lInlHRig= as the parameter NiMlInRH 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 NiYlIkFHJSJCRyUiQ0clIkRH in Figure 13.2.NiMtJSdNQVRSSVhHNiM3KDcnJSZQb2ludEclJih4LHkpRyUjeCdHJSN5J0clK0NvbmNsdXNpb25HNyclJi0tLS0tRyUqLS0tLS0tLS0tRyUwLS0tLS0tLS0tLS0tLS0tR0YwJT4tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLUc3JyUiQUclKChwb3MsMClHJSt4J349fnh+Pn4wRyUseSd+PX40eH4+fjBHJTZib3Rofnh+YW5kfnl+aW5jcmVhc2VHNyclIkJHJSgoMCxwb3MpRyUseCd+PX4yeX4+fjBHJSx5J349fi15fjx+MEclOXh+aW5jcmVhc2VzLH55fmRlY3JlYXNlc0c3JyUiQ0clKChuZWcsMClHJSp4J349fnh+PDBHJSx5J349fjR4fjx+MEclNmJvdGh+eH5hbmR+eX5kZWNyZWFzZUc3JyUiREclKCgwLG5lZylHJSt4J349fjJ5fjwwRyUseSd+PX4teX4+fjBHJTl4fmRlY3JlYXNlcyx+eX5pbmNyZWFzZXNH Table 13.2The 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.<Text-field bookmark="13.C-2" layout="Heading 2" style="Heading 2">13.C-2 Example 2 (Center)</Text-field>The equilibrium point at the origin is a center for the differential equationsq5 := diff(x(t),t) = x(t) - 2*y(t);
q6 := diff(y(t),t) = x(t) - y(t);A set of initial points isinits := [seq([0,k,0],k=1..3)];The initial points have three values, the first, representing the initial time. Here, that is NiMvJSJ0RyIiIQ==, the time at which the system point is at that initial point whose NiMlInhH and NiMlInlH coordinates are the second and third value in the initial condition.The resulting phase portrait, Figure 13.3, isDEplot([q5,q6],[x(t),y(t)], t=0..2*Pi, x=-4..4, y=-4..4, inits, stepsize=.1, linecolor=black, arrows=medium, title="Figure 13.3");A center is surrounded by closed trajectories.<Text-field bookmark="13.C-3" layout="Heading 2" style="Heading 2">13.C-3 Example 3 (Node)</Text-field>The following linear system of differential equations has an equilibrium point at the origin, and it is a node.q3 := diff(x(t),t) = 7*x(t) - 2*y(t);
q4 := diff(y(t),t) = 15*x(t) - 4*y(t);With the following initial points,inits := [seq([0,5-k,0],k=3..7)];the behavior of solutions near the node can be determined.The phase portrait, Figure 13.4, is drawn viaDEplot([q3,q4],[x(t),y(t)], t=-10..1/10, x=-4..4, y=-8..8, inits, stepsize=.1, arrows=medium, linecolor=black, title="Figure 13.4");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.<Text-field bookmark="13.C-4" layout="Heading 2" style="Heading 2">13.C-4 Example 4 (Spiral)</Text-field>The equilibrium point at the origin is a spiral for the differential equationsq7 := diff(x(t),t) = -x(t) - 2*y(t);
q8 := diff(y(t),t) = x(t) - 3*y(t);A set of initial points isinits := [seq([0,k,0],k=7..9), seq([0,k,0],k=-9..-7)];and the phase portrait, Figure 13.5, isDEplot([q7,q8],[x(t),y(t)], t=-2..2, x=-12..12, y=-30..30, inits, stepsize=.1, linecolor=black, arrows=medium, title="Figure 13.5");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.<Text-field bookmark="13.C-5" layout="Heading 2" style="Heading 2">13.C-5 Eigenvalues and the Phase Portrait</Text-field>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 signsCENTER: Pure imaginaryNODE: Same signSPIRAL: Complex conjugates=================================== Table 13.3<Text-field bookmark="13.C-6" layout="Heading 2" style="Heading 2">13.C-6 Representative Examples</Text-field>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 => NiQtJSRleHBHNiMsJCUidEchIiItRiQ2IywkKiYiIiMiIiJGJ0YuRig=Improper Node In: -1,-1 => NiQtJSRleHBHNiMsJCUidEchIiIqJkYnIiIiRiNGKg==Proper Node Out: 1,2 => NiQtJSRleHBHNiMlInRHLUYkNiMqJiIiIyIiIkYmRis=Improper Node Out: 1,1 => NiQtJSRleHBHNiMlInRHKiZGJiIiIkYjRig=Saddle: -1,2 => NiQtJSRleHBHNiMsJCUidEchIiItRiQ2IyomIiIjIiIiRidGLQ==Center: NiQqJiIiJCIiIiUiaUdGJSwkRiMhIiI= => NiQtJSRzaW5HNiMqJiIiJCIiIiUidEdGKC0lJGNvc0dGJQ==Spiral In: NiQsJiIiIiEiIiomIiIjRiQlImlHRiRGJCwmRiRGJUYmRiU= => NiQqJi0lJGV4cEc2IywkJSJ0RyEiIiIiIi0lJGNvc0c2IyomIiIjRipGKEYqRioqJkYkRiotJSRzaW5HRi1GKg==Spiral Out: NiQsJiIiIkYkKiYiIiNGJCUiaUdGJEYkLCZGJEYkRiUhIiI= => NiQqJi0lJGV4cEc2IyUidEciIiItJSRjb3NHNiMqJiIiI0YoRidGKEYoKiZGJEYoLSUkc2luR0YrRig====================================================================== Table 13.4<Text-field bookmark="13.D" layout="Heading 1" style="Heading 1">13.D Phase Portraits for 2d Nonlinear Autonomous Systems</Text-field>The two-dimensional nonlinear autonomous systemNiMvKiYlI2R4RyIiIiUjZHRHISIiLSUiZkc2JCUieEclInlHNiMvKiYlI2R5RyIiIiUjZHRHISIiLSUiZ0c2JCUieEclInlHcan have more than one equilibrium point. Each real solution of the simultaneous equations NiMvLSUiZkc2JCUieEclInlHIiIh = NiMtJSJnRzYkJSJ4RyUieUc= 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.<Text-field bookmark="13.D-1" layout="Heading 2" style="Heading 2">13.D-1 Example 5</Text-field>A direction field for the nonlinear autonomous system of differential equationssys1 := { diff( x(t), t ) = x(t) * ( 6 - 3*x(t) - 2*y(t) ), diff( y(t), t ) = y(t) * ( 5 - x(t) - y(t) ) };appears in Figure 13.6.dir_field1 := DEplot( sys1, {x(t),y(t)}, t=0..1, x=-5..5, y=-5..10, scene=[x,y] ):display(dir_field1, title="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,equil_soln1 := solve( map(rhs,sys1), {x(t),y(t)} );which are easier to read as the list of pointsequil_pts1 := map( subs, [equil_soln1], [x(t),y(t)] );Figure 13.7 shows these points as blue circles on the direction field of Figure 13.6.equil_plot1 := pointplot( equil_pts1, symbol=CIRCLE, symbolsize=18, color=BLUE ):display( [dir_field1,equil_plot1], axes=BOXED, title="Figure 13.7" );Figure 13.7 suggests that the equilibrium solution at the origin is a source, the one at NiMtJSFHNiQiIiEiIiY= is a sink, and the ones at NiMtJSFHNiQiIiMiIiE= and NiMtJSFHNiQsJCIiJSEiIiIiKg== 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 listpt_list := [ op(equil_pts1), [5,0], [1,10], [-2,10], [0,.2], [0,.5], [0,-.5] ];Perturb these points a small amount in each component to obtain the final listh := 0.1:q1 := [seq( seq( X+dX, dX=[ [h,h], [h,-h], [-h,h], [-h,-h] ] ), X=pt_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.IC1 := map( p->[x(0)=p[1],y(0)=p[2]], q1 ):phase_port1 := DEplot( sys1, {x(t),y(t)}, t=0..2, x=-5..5, y=-5..10, IC1, scene=[x,y], arrows=NONE, linecolor=green ):display( [dir_field1, equil_plot1, phase_port1], axes=BOXED, title="Figure 13.8" );<Text-field bookmark="13.D-2" layout="Heading 2" style="Heading 2">13.D-2 Example 6 (The Nonlinear Pendulum)</Text-field>The equationssys2 := { diff( x(t), t ) = y(t), diff( y(t), t ) = -sin(x(t)) };formulate Newton's Second Law of Motion for a simple plane pendulum in which NiMtJSJ4RzYjJSJ0Rw== is the angle between the stable hanging position and the pendulum's position at time NiMlInRH; NiMvJSJ5RyomJSNkeEciIiIlI2R0RyEiIg== is the angular velocity of the pendulum bob. Thus, for every integer NiMlImtH, the point NiMvLSUhRzYkJSJ4RyUieUctRiU2JCooIiIjIiIiJSJrR0YtJSNQaUdGLSIiIQ== corresponds to the stable hanging position and the point NiMtJSFHNiQqJiwmKiYiIiMiIiIlImtHRipGKkYqISIiRiolI1BpR0YqIiIh 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.dir_field2 := DEplot( sys2, {x(t),y(t)}, t=0..1, x=-10..10, y=-5..5, scene=[x,y], dirgrid=[25,25], scaling=constrained ):display(dir_field2, title="Figure 13.9");Although solve finds only one equilibrium solution for this systemsolve( map(rhs,sys2), {x(t),y(t)} );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_EnvAllSolutions := true:temp_soln := solve( map(rhs,sys2), {x(t),y(t)} ):_EnvAllSolutions := false:
unassign('k');
Q := remove(has,indets(temp_soln),{t})[1]:
equil_soln2 := eval(temp_soln, Q=k);(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, NiMlJC4uLkc=. To avoid such complications, the operation of solving for the equilibrium solutions has been coded to return the symbol NiMlImtH for the arbitrary integer.)To plot the equilibrium solutions in the direction field, convert the general expression to a pointq1 := op(map( subs, [equil_soln2], [x(t),y(t)] ));and substitute appropriate integers to get the list of pointsequil_pts2 := [ q1 $ k=-3..3 ];The desired plot then appears in Figure 13.10.equil_plot2 := pointplot( evalf(equil_pts2), symbol=CIRCLE, symbolsize=18, color=BLUE ):display( [dir_field2, equil_plot2], title="Figure 13.10" );The direction field suggests different behaviors for the equilibria at even and odd integer multiples of NiMlI1BpRw==. Solution curves will help to further elaborate on these differences. With initial conditions selected along the NiMlInhH-axis, and from the "inflow" portions of the window for the direction field, i.e., from NiQvJSJ4RyIjNTIsJCIiJiEiIiUieUc=NiMyJSFHIiIh and NiMvJSJ4RywkIiM1ISIi, NiMyIiIhJSJ5Rw==NiMyJSFHIiIm, the totality of initial points isIC2 := [ [x(0)=10,y(0)=i] $ i=-4..-1, [x(0)=-10,y(0)=i] $ i=1..4, [x(0)=2*i,y(0)=0] $ i=-4..4 ];
The combined direction field and phase portrait so generated is seen in Figure 13.11.phase_port2 := DEplot( sys2, {x(t),y(t)}, t=0..15, x=-10..10, y=-5..5, IC2, scene=[x,y], arrows=NONE, linecolor=CYAN ):display( [dir_field2, equil_plot2, phase_port2], title="Figure 13.11" );
Thus, for any integer NiMlImtH, the point NiMtJSFHNiQqKCIiIyIiIiUia0dGKCUjUGlHRigiIiE= is a center and NiMtJSFHNiQqJiwmKiYiIiMiIiIlImtHRipGKkYqISIiRiolI1BpR0YqIiIh is a saddle point.<Text-field bookmark="13.D-3" layout="Heading 2" style="Heading 2">13.D-3 Example 7</Text-field>If damping is added to the simple plane pendulum of Example 6, the system equations becomesys3 := { diff( x(t), t ) = y(t), diff( y(t), t ) = -sin(x(t)) - y(t) };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.dir_field3 := DEplot( sys3, {x(t),y(t)}, t=0..1, x=-10..10, y=-5..5, scene=[x,y], dirgrid=[30,30], scaling=constrained ):display(dir_field3, title="Figure 13.12");
The differences in the direction fields are obvious. Nonetheless, the system has the same equilibrium solutions, namely,_EnvAllSolutions := true:temp_soln := solve( map(rhs,sys3), {x(t),y(t)} ):_EnvAllSolutions := false:
unassign('k');
Q := remove(has,indets(temp_soln),{t})[1]:
equil_soln3 := eval(temp_soln, Q=k);where the integer NiMlImtH has the same explanation as appeared in Example 6.Again writing this solution as the pointq3 := op(map( subs, [equil_soln3], [x(t),y(t)] ));and generating the following list of equilibrium points,equil_pts3 := [ q3 $ k=-3..3 ];Figure 13.13 showing both the equilibrium points and the direction field, is generated.equil_plot3 := pointplot( evalf(equil_pts3), symbol=CIRCLE, symbolsize=18, color=BLUE ):display( [dir_field3,equil_plot3], title="Figure 13.13" );For the phase portrait, select initial conditions from the top and bottom edges of the viewing window. With the initial pointsIC3 := [ [x(0)=i,y(0)=5] $ i=-10..10, [x(0)=i,y(0)=-5] $ i=-10..10 ];the phase portrait in Figure 13.14 is generated.sol_traj3 := DEplot( sys3, {x(t),y(t)}, t=0..5, x=-10..10, y=-5..5, IC3, scene=[x,y], arrows=NONE, linecolor=green ):display( [dir_field3,equil_plot3,sol_traj3], title="Figure 13.14" );From the direction field and equilibrium solution, it appears as though the equilibrium solutions at NiMtJSFHNiQqKCIiIyIiIiUia0dGKCUjUGlHRigiIiE= that were centers in the undamped system are inward nodes (attractors, or sinks) in the damped system while the equilibria at NiMtJSFHNiQqJiwmKiYiIiMiIiIlImtHRipGKkYqISIiRiolI1BpR0YqIiIh are saddle points for both systems.<Text-field bookmark="13.E" layout="Heading 1" style="Heading 1">13.E Nonautonomous Systems</Text-field>Any nonautonomous system can be reformulated as an autonomous system via the introduction of the independent variable, say NiMlInRH, and augmenting the system with the (trivial) differential equationNiMvLSUlRGlmZkc2JCUidEdGJyIiIg==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 NiMvJSJzRyUidEc=. Then all occurrences of NiMlInRH in the original system are replaced with NiMlInNH and the auxiliary ODE becomesNiMvKiYlI2R0RyIiIiUjZHNHISIiRiY=<Text-field bookmark="13.E-1" layout="Heading 2" style="Heading 2">13.E-1 Example 8</Text-field>The nonautonomous systemsys4 := { diff( x(t), t ) = y(t), diff( y(t), t ) = -sin(x(t)) - sin(t) };can be converted to an autonomous system with the help of the dchange command from the PDEtools package. This results in the three-dimensional systemsys4a := dchange( t=s, sys4, s ) union { diff( t(s), s ) = 1 };Because this is a three-dimensional system, it is not realistic to draw a direction field. Instead, a sample of solution curves with initial conditionsIC4a := [ [x(0)=1,y(0)=1,t(0)=0], [x(0)=1,y(0)=-1,t(0)=0], [x(0)=-1,y(0)=1,t(0)=0], [x(0)=-1,y(0)=-1,t(0)=0] ]:obtained with the DEplot3d command, appears in Figure 13.15.sol_traj4a := DEplot3d( sys4a, {x(s),y(s),t(s)}, s=0..10, IC4a, scene=[t,x,y], stepsize=0.2 ):display(sol_traj4a, title="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 withdisplay( sol_traj4a, orientation=[0,90] );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 NiMlInRH-NiMlInhH view is obtained withdisplay( sol_traj4a, orientation=[-90,0] );and the t-y view is obtained withdisplay( sol_traj4a, orientation=[-90,90] );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 useIC4 := [ [x(0)=1,y(0)=1], [x(0)=1,y(0)=-1], [x(0)=-1,y(0)=1], [x(0)=-1,y(0)=-1] ]:DEplot3d( sys4, {x(t),y(t)}, t=0..10, IC4, scene=[t,x,y], stepsize=0.2 );<Text-field bookmark="13.F" layout="Heading 1" style="Heading 1">13.F Stability of an Equilibrium Point</Text-field>An isolated equilibrium point for an autonomous system of differential equations is said to beOrbitally (Neutrally) Stable if all solutions sufficiently close to the equilibrium point remain close to that equilibrium point as NiMsJiUidEciIiIqJiUiPkdGJSUpaW5maW5pdHlHRiUhIiI=,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 NiMsJiUidEciIiIqJiUiPkdGJSUpaW5maW5pdHlHRiUhIiI=,Unstable if at least one solution starting close to the equilibrium point does not remain close to that equilibrium point as NiMsJiUidEciIiIqJiUiPkdGJSUpaW5maW5pdHlHRiUhIiI=.Special noteSome 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)NiMvLSUjeCdHNiMlInRHLSUieUdGJg==NiMvJSJBRy0lJ01BVFJJWEc2IzckNyQiIiEiIiI3JEYrRio= + NiMlImlH Figure 13.16NiMvLSUjeSdHNiMlInRHLCQtJSJ4R0YmISIi___________________________________________________________________________(13.10)NiMvLSUjeCdHNiMlInRHLCYqJiIiKSIiIi0lInhHRiZGKyEiIiomIiIkRistJSJ5R0YmRitGKw==NiMvJSJBRy0lJ01BVFJJWEc2IzckNyQsJCIiKSEiIiIiJDckLCQiIzVGLEYtNiQsJCIiJCEiIiwkIiIjRiU= Figure 13.17NiMvLSUjeSdHNiMlInRHLCYqJiIjNSIiIi0lInhHRiZGKyEiIiomIiIkRistJSJ5R0YmRitGKw==___________________________________________________________________________(13.11)NiMvLSUjeCdHNiMlInRHLSUieEdGJg==NiMvJSJBRy0lJ01BVFJJWEc2IzckNyQiIiIiIiE3JCIiIywkRiohIiI= + 1 Figure 13.18NiMvLSUjeSdHNiMlInRHLCYqJiIiIyIiIi0lInhHRiZGK0YrLSUieUdGJiEiIg==____________________________________________________________________________ Table 13.5<Text-field bookmark="13.F-1" layout="Heading 2" style="Heading 2">13.F-1 Example 9</Text-field>The system in Example 9, Table 13.5, with equationsq1 := diff(x(t),t) = y(t);
q2 := diff(y(t),t) = -x(t);has an orbitally stable equilibrium point at the origin. To see what this means, form the coefficient matrix, A,(Aneg,b) := GenerateMatrix(eval([q1,q2],[x(t)=_x,y(t)=_y]),[_x,_y]):A := -Aneg;and compute its eigenvalues.Eigenvalues(A);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 pointsinits := [seq([0,k,0],k=1..3)];and use Maple's DEplot command to obtain a phase portrait.DEplot([q1,q2], [x(t),y(t)], t=0..2*Pi, inits, linecolor=black, stepsize=.1, scaling=constrained, arrows=medium, title="Figure 13.16");<Text-field bookmark="13.F-2" layout="Heading 2" style="Heading 2">13.F-2 Example 10</Text-field>The system in Example 10, Table 13.5, whose equations areq3 := diff(x(t),t) = -8*x(t) + 3*y(t);
q4 := diff(y(t),t) = -10*x(t) + 3*y(t);has an asymptotically stable equilibrium point at the origin. To see what this means, form the coefficient matrix, A(Aneg,b) := GenerateMatrix(eval([q3,q4],[x(t)=_x,y(t)=_y]),[_x,_y]):A := -Aneg;and compute the eigenvalues of A.Eigenvalues(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 NiMsJiUidEciIiIqJiUiPkdGJSUpaW5maW5pdHlHRiUhIiI=. For example, consider the solution starting at (NiQlImFHJSJjRw==), obtained viaq5 := dsolve({q3,q4, x(0)=a,y(0)=c},{x(t),y(t)});Extract the component solutionsxt := eval(x(t), q5);
yt := eval(y(t), q5);and compute the limitsLimit(x,t=infinity) = limit(xt,t=infinity);
Limit(y,t=infinity) = limit(yt,t=infinity);As promised, an arbitrary solution starting near the origin, not only remains near the origin, but tends to the origin as NiMsJiUidEciIiIqJiUiPkdGJSUpaW5maW5pdHlHRiUhIiI=. Hence, all solutions, not only the ones near the origin, tend towards the origin with increasing NiMlInRH, as seen in the phase portrait, Figure 13.17. To obtain this figure, pick the initial pointsinits := [seq([0,-3,k],k=-3..3), seq([0,3,k],k=-3..3)];and use Maple's DEplot command to obtainDEplot([q3,q4], [x(t),y(t)], t=0..1, inits, linecolor=black, stepsize=.1, scaling=constrained, dirgrid=[10,20], arrows=medium, title="Figure 13.17");<Text-field bookmark="13.F-3" layout="Heading 2" style="Heading 2">13.F-3 Example 11</Text-field>The system in Example 11, Table 13.5, whose equations areq6 := diff(x(t),t) = x(t);
q7 := diff(y(t),t) = 2*x(t) - y(t);has an unstable equilibrium point at the origin. To see what this means, form the coefficient matrix, A,(Aneg,b) := GenerateMatrix(eval([q6,q7],[x(t)=_x,y(t)=_y]),[_x,_y]):A := -Aneg;and compute the eigenvalues of A.Eigenvalues(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 NiMlInRH 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 preciselyq8 := dsolve({q6,q7, x(0)=0,y(0)=a},{x(t),y(t)}):
x1 := eval(x(t), q8);
y1 := eval(y(t), q8);for whichLimit(x,t=infinity) = limit(x1,t=infinity);
Limit(y,t=infinity) = limit(y1,t=infinity);Any other solution with NiMqKCUiYUciIiIlIj5HRiUiIiFGJQ==, given byq9 := dsolve({q6,q7, x(0)=a,y(0)=c},{x(t),y(t)}):
x2 := eval(x(t), q9);
y2 := eval(y(t), q9);will not remain near the origin, as shown by the following limits. Indeed, Maple computes the limitsLimit(x,t=infinity) = limit(x2,t=infinity) assuming a>0;
Limit(y,t=infinity) = limit(y2,t=infinity) assuming a>0;for NiMqKCUiYUciIiIlIj5HRiUiIiFGJQ==, whereas if NiMyJSJhRyIiIQ==, Maple computesLimit(x,t=infinity) = limit(x2,t=infinity) assuming a<0;
Limit(y,t=infinity) = limit(y2,t=infinity) assuming a<0;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 pointsinits := [seq([0,k,2],k=-2..2), seq([0,k,-2],k=-2..2)];and use the DEplot command.DEplot([q6,q7], [x(t),y(t)], t=-2..2, x=-5..5, y=-5..5, inits, linecolor=black, stepsize=.1, arrows=medium, scaling=constrained, title="Figure 13.18");Solutions starting in the half-plane NiMqKCUieEciIiIlIj5HRiUiIiFGJQ== tend to NiMlKWluZmluaXR5Rw==, whereas solutions starting in the half-plane NiMyJSJ4RyIiIQ== tend to NiMsJCUpaW5maW5pdHlHISIi.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.<Text-field bookmark="13.F-4" layout="Heading 2" style="Heading 2">13.F-4 Summary</Text-field>Table 13.6 summarizes the relationships between phase portraits, eigenvalues, and stability.Equilibrium Point Eigenvalues Conditions Stability=========== ======== ============ ===============Node (in) real and negative Asymptotically stableSpiral (in) NiMlImFH + NiMqJiUiYkciIiIlImlHRiU=NiMyJSJhRyIiIQ== Asymptotically stableCenter + NiMqJiUiYkciIiIlImlHRiU= pure imaginary Orbitally stableNode (out) real and positive UnstableSpiral (out) NiMlImFH + NiMqJiUiYkciIiIlImlHRiU=NiMqKCUiYUciIiIlIj5HRiUiIiFGJQ== UnstableSaddle real, with opposite signs Unstable============================================================================== Table 13.6<Text-field bookmark="13.G" layout="Heading 1" style="Heading 1">13.G Linearization of Nonlinear Systems </Text-field><Text-field bookmark="13.G-1" layout="Heading 2" style="Heading 2">13.G-1 Tangent-Line Approximation</Text-field>Near the point of contact, the curve NiMvJSJ5Ry0lImZHNiMlInhH is approximated by its tangent line. For example, the functionf := (x-1)^2+2;and its tangent at the point NiMtJSFHNiQiIiMiIiQ= appear in Figure 13.19.Student[Calculus1][Tangent]( f, x=2, output=plot, title="Figure 13.19", view=[0..3,0..5], labels=[x,y], scaling=constrained );Near the point NiMtJSFHNiQiIiMiIiQ= the function NiMtJSJmRzYjJSJ4Rw== very much resembles its tangent line NiMvJSJ5RywmKiYiIiMiIiIlInhHRihGKEYoISIi. Analytically, the tangent line is the Taylor polynomial of degree 1, NiMsJiomIiIjIiIiJSJ4R0YmRiZGJiEiIg==.The tangent line can be constructed by the following calculations. First, the slope at NiMvJSJ4RyIiIw==,m := eval(diff(f,x), x=2);then, the equation of the tangent line via the point-slope formulaNiMvJSJ5RywmKiYlIm1HIiIiLCYlInhHRigmRio2I0YoISIiRihGKCZGJEYsRig=Y := m*(x - 2) + 3;Analytically, the tangent line is a Taylor polynomial of degree 1, obtained in Maple viay = convert(taylor(f, x = 2, 2), polynom);Replacing a curve NiMtJSJmRzYjJSJ4Rw== with its tangent-line approximation is known as linearizing in one dimension.<Text-field bookmark="13.G-2" layout="Heading 2" style="Heading 2">13.G-2 Tangent-Plane Approximation</Text-field>The tangent plane approximates a function of two variables. This is the two-dimensional version of linearization.For example, at the point where NiMvLSUhRzYkJSJ4RyUieUctRiU2JCIiI0Yr, the tangent plane approximating the functionf := 3*x^2 + 2*y^2 - 5*x -7*y -9;is NiMvJSJ6RywoKiYiIigiIiIlInhHRihGKCUieUdGKCIjSCEiIg==. This first-degree polynomial approximation of NiMtJSJmRzYkJSJ4RyUieUc= is just the multivariable Taylor polynomial of degree one.To construct this tangent plane, first find the corresponding z-coordinate at NiMtJSFHNiQiIiNGJg==, the function value NiMtJSJmRzYkIiIjRiY=. Thus,z1 := eval(f, {x = 2, y = 2});To construct a plane tangent at NiMtJSFHNiUiIiNGJiwkIiM4ISIi obtain a vector normal to the surface NiMvJSJ6Ry0lImZHNiQlInhHJSJ5Rw==. The general normal vector isNiMvJSJORy0lJ01BVFJJWEc2IzclNyMsJCYlImZHNiMlInhHISIiNyMsJCZGLDYjJSJ5R0YvNyMiIiI=so computeN := < -diff(f,x), -diff(f,y), 1 >;and evaluate it at NiMtJSFHNiQiIiNGJg== to obtainN1 := eval(N, {x = 2, y = 2});Once the normal for the tangent plane is known, write the equation of the plane asq := N1 . <x,y,z> = d;To determine the value of d, force the plane to go through the point NiMtJSFHNiUiIiNGJiwkIiM4ISIi. Hence,eval(q, {x = 2, y = 2, z = -13});from which the implicit form of the tangent plane is found to beq1 := eval(q, d = -29);Explicitly solving this equation for NiMvJSJ6Ry1GJDYkJSJ4RyUieUc= yieldsZ := solve(q1,z);That is one way to find the equation of the tangent plane. Alternatively, at the point where NiMvLSUhRzYkJSJ4RyUieUctRiU2JCIiI0Yr, generate the multivariable Taylor polynomial of degree one. In Maple, use the mtaylor command to create a linear approximation to the function NiMtJSJmRzYkJSJ4RyUieUc=, thereby getting the explicit form of the equation of the tangent plane, namely, NiMvJSJ6Ry1GJDYkJSJ4RyUieUc=. This results inz = mtaylor(f,[x=2, y=2], 2);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.<Text-field bookmark="13.G-3" layout="Heading 2" style="Heading 2">13.G-3 Linearizing Autonomous Systems</Text-field>Near an equilibrium point of the autonomous nonlinear systemNiMvKiYlI2R4RyIiIiUjZHRHISIiLSUiRkc2JC0lInhHNiMlInRHLSUieUdGLg==NiMvKiYlI2R5RyIiIiUjZHRHISIiLSUiR0c2JC0lInhHNiMlInRHLSUieUdGLg==the behavior of solutions can be studied by linearizing the functions NiMtJSJGRzYkJSJ4RyUieUc= and NiMtJSJHRzYkJSJ4RyUieUc= at each equilibrium point.Let (NiQlImFHJSJiRw==) be an equilibrium point determined by the equations NiMvLSUiRkc2JCUiYUclImJHLSUiR0dGJg== = 0. Then, lumping the higher-order remainder terms into the functions NiMtJSJmRzYkJSJ4RyUieUc= and NiMtJSJnRzYkJSJ4RyUieUc=, the first-degree Taylor expansionsNiMvLSUiRkc2JCUieEclInlHLCotRiU2JCUiYUclImJHIiIiKiYtJkYlNiNGJ0YrRi4sJkYnRi5GLCEiIkYuRi4qJi0mRiU2I0YoRitGLiwmRihGLkYtRjRGLkYuLSUiZkdGJkYuNiMvLSUiR0c2JCUieEclInlHLCotRiU2JCUiYUclImJHIiIiKiYtJkYlNiNGJ0YrRi4sJkYnRi5GLCEiIkYuRi4qJi0mRiU2I0YoRitGLiwmRihGLkYtRjRGLkYuLSUiZ0dGJkYubecome NiMvLSUiRkc2JCUieEclInlHLCgqJi0mRiU2I0YnNiQlImFHJSJiRyIiIiwmRidGMUYvISIiRjFGMSomLSZGJTYjRihGLkYxLCZGKEYxRjBGM0YxRjEtJSJmR0YmRjE=NiMvLSUiR0c2JCUieEclInlHLCgqJi0mRiU2I0YnNiQlImFHJSJiRyIiIiwmRidGMUYvISIiRjFGMSomLSZGJTYjRihGLkYxLCZGKEYxRjBGM0YxRjEtJSJnR0YmRjE=because NiMtJSFHNiQlImFHJSJiRw== is an equilibrium point at which NiMvLSUiRkc2JCUiYUclImJHLSUiR0dGJg== = 0. The nonlinear autonomous system now readsDefine the change of variablesNiMvLSUidUc2IyUidEcsJi0lInhHRiYiIiIlImFHISIiNiMvLSUidkc2IyUidEcsJi0lInlHRiYiIiIlImJHISIiso thatNiMvKiYlI2R1RyIiIiUjZHRHISIiKiYlI2R4R0YmRidGKA==NiMvKiYlI2R2RyIiIiUjZHRHISIiKiYlI2R5R0YmRidGKA==and the system now readsNiMpLSUnTUFUUklYRzYjNyQ3IyUidUc3IyUidkclIidH = NiMtJSdNQVRSSVhHNiM3JDckLSYlIkZHNiMlInhHNiQlImFHJSJiRy0mRio2IyUieUdGLTckLSYlIkdHRitGLS0mRjdGMkYtNiMtJSdNQVRSSVhHNiM3JDcjJSJ1RzcjJSJ2Rw== + NiMtJSdNQVRSSVhHNiM3JDcjKiYpJSJmRyUiKkciIiI2JCUidUclInZHRiw3IyomKSUiZ0dGK0YsRi1GLA==or u' = Au + Pwhere u = NiMtJSdNQVRSSVhHNiM3JDcjJSJ1RzcjJSJ2Rw==, the matrix A is the Jacobian MatrixNiMvJSJBRy0lJ01BVFJJWEc2IzckNyQtJiUiRkc2IyUieEc2JCUiYUclImJHLSZGLDYjJSJ5R0YvNyQtJiUiR0dGLUYvLSZGOUY0Ri8=and the functions NiMqJiklImZHJSIqRyIiIjYkJSJ1RyUidkdGJw== and NiMqJiklImdHJSIqRyIiIjYkJSJ1RyUidkdGJw== are just the functions NiMtJSJmRzYkJSJ4RyUieUc= and NiMtJSJnRzYkJSJ4RyUieUc= under the change of variables NiQvJSJ4RywmJSJ1RyIiIiUiYUdGJy8lInlHLCYlInZHRiclImJHRic=.The equilibrium point (NiQlImFHJSJiRw==) in the xy-plane has been translated to the equilibrium point NiMtJSFHNiQiIiFGJg== 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 NiMtJSFHNiQiIiFGJg==.
The behavior of the nonlinear xy-system at (NiQlImFHJSJiRw==) can often be deduced by analyzing the behavior of the linear uv-system at NiMtJSFHNiQiIiFGJg==, as will be shown shortly.<Text-field bookmark="13.G-3A" layout="Heading 3" style="Heading 3">13.G-3A Example 12</Text-field>The differential equations NiMvJSN4J0cqJiUieEciIiIsKCIjZ0YnKiYiIiRGJ0YmRichIiIqJiIiJUYnJSJ5R0YnRixGJw==NiMvJSN5J0cqJiUieUciIiIsKCIjVUYnKiYiIiRGJ0YmRichIiIqJiIiI0YnJSJ4R0YnRixGJw==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 NiMtJSJGRzYkJSJ4RyUieUc= and NiMtJSJHRzYkJSJ4RyUieUc= areF := x*(60 - 3*x - 4*y);
G := y*(42 - 3*y - 2*x);The equilibrium points are obtained by solving the equations NiMvLSUiRkc2JCUieEclInlHLSUiR0dGJg== = 0, done in Maple viaq := solve({F=0, G=0}, {x,y});Extract and label these four equilibrium points asfor k from 1 to 4 do
P||k := eval([x,y], q[k]);
od;The Jacobian matrix isNiMvJSJBRy0lJ01BVFJJWEc2IzckNyQsKCIjZyIiIiomIiInRiwlInhHRiwhIiIqJiIiJUYsJSJ5R0YsRjAsJComRjJGLEYvRixGMDckLCQqJiIiI0YsRjNGLEYwLCgiI1VGLComRi5GLEYzRixGMComRjlGLEYvRixGMA==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 commandA := Jacobian([F,G],[x,y]);At each of the equilibrium points, the Jacobian matrix evaluates tofor k from 1 to 4 do
A||k := eval(A, q[k]);
od;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 NiYmJSJBRzYjIiIiJkYkNiMiIiMmRiQ2IyIiJCZGJDYjIiIl.for k from 1 to 4 do
Eigenvalues(A||k);
od;To determine the signs of the eigenvalues from the matrix NiMmJSJBRzYjIiIl, obtain floating point equivalents of its exact eigenvalues.evalf(Eigenvalues(A4));For the linearized systems we haveP1 = (0, 0) - node outP2 = (20, 0) - saddleP3 = (0, 14) - saddleP4 = (12, 6) - node inRules 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.<Text-field bookmark="13.G-4" layout="Heading 2" style="Heading 2">13.G-4 Interpretation Rules</Text-field>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 (NiMlJ2xhbWJkYUc='s real and distinct) in the linearized system remain the same in the nonlinear system.4. Improper nodes (NiMlJ2xhbWJkYUc='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 systemThis 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.<Text-field bookmark="13.G-5" layout="Heading 2" style="Heading 2">13.G-5 Example 13</Text-field>To analyze the stability of the nonlinear autonomous system NiMvJSN4J0cqKCIiJiIiIiUieEdGJywmRidGJyUieUchIiJGJw==NiMvJSN5J0cqJiUieUciIiIsJiomIiIkRiclInhHRidGJyIjNyEiIkYnobtain the equilibrium points NiMtJSFHNiQiIiFGJg== and NiMtJSFHNiQiIiUiIiI=, and the Jacobian matrixNiMvJSJBRy0lJ01BVFJJWEc2IzckNyQsJiIiJiIiIiomRitGLCUieUdGLCEiIiwkKiZGK0YsJSJ4R0YsRi83JComIiIkRixGLkYsLCYqJkY1RixGMkYsRiwiIzdGLw==Evaluate A at each equilibrium point, and compute the corresponding eigenvalues, obtaining the matricesNiMtJSdNQVRSSVhHNiM3JDckIiImIiIhNyRGKSwkIiM3ISIi and NiMtJSdNQVRSSVhHNiM3JDckIiIhLCQiIz8hIiI3JCIiJEYoand the respective eigenvalues NiQiIiYsJCIjNyEiIg==, and + NiMqKCIiIyIiIiUiaUdGJS0lJXNxcnRHNiMiIzpGJQ==. The point NiMtJSFHNiQiIiFGJg== is a saddle, and unstable, for the linearized system; therefore, it is a saddle, and unstable, for the nonlinear system. The point NiMtJSFHNiQiIiUiIiI= 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 NiMtJSFHNiQiIiUiIiI= may well be a center in the nonlinear case.To implement these calculations in Maple, enter the funtions NiMtJSJGRzYkJSJ4RyUieUc= and NiMtJSJHRzYkJSJ4RyUieUc=.F := (x,y) -> 5*x*(1 - y):
G := (x,y) -> y*(3*x - 12):
'F'(x,y) = F(x,y);
'G'(x,y) = G(x,y);Then, solve the equations x' = y' = 0 for equilibrium points, obtainingq := solve({F(x,y) = 0, G(x,y) = 0}, {x,y});There are two equilibrium points that can now be carefully identified and labeled.for k from 1 to 2 do
P||k := eval([x,y], q[k]);
od;Compute the generic Jacobian Matrix of first partial derivatives.A := Jacobian([F(x,y),G(x,y)],[x,y]);Evaluate the Jacobian Matrix at each of the two equilibrium points.for k from 1 to 2 do
A||k := eval(A, q[k]);
od;Obtain the eigenvalues of each of the matrices NiQmJSJBRzYjIiIiJkYkNiMiIiM=.for k from 1 to 2 do
Eigenvalues(A||k);
od;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.q1 := diff(x(t),t) = F(x(t),y(t));
q2 := diff(y(t),t) = G(x(t),y(t));and select the initial conditionsinits := [seq([0,10,2*k],k=-3..3), seq([0,-3,2*k],k=-3..3)];Then, use the DEplot command to construct a phase portrait.display(DEplot([q1,q2],[x(t),y(t)],t=-1..1, x = -3..12, y = -6..8, inits, stepsize=.01, linecolor = black, arrows=medium), labels=[` x`,y], labelfont=[TIMES,ITALIC,12], title="Figure 13.20");The phase portrait suggests that (4, 1) may indeed be a center in the nonlinear case.[Back to ODE Powertool Table of Contents]