ORDINARY DIFFERENTIAL EQUATIONS POWERTOOLLesson 7 -- Application: The Spruce BudwormProf. 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 7</Text-field>7.A Biological Background7.B Bifurcation Analysis7.C References<Text-field layout="Heading 1" style="Heading 1">Initialization</Text-field>restart;with( DEtools ):with( plots ):with( PDEtools ):<Text-field bookmark="7.A" layout="Heading 1" style="Heading 1">7.A Biological Background</Text-field>A classical example of bifurcation in nature is the interaction between the spruce budworm and balsam fir forests in North America. The basic model for the budworm population is a logistic model with a predation term:f[logistic] := r*B*(1-B/K);f[predation] := beta*B^2/(alpha^2+B^2);where, in the absence of predation, NiMlInJH is the intrinsic growth rate and NiMlIktH is the carrying capacity. All four parameters, NiMlInJH, NiMlIktH, NiMlJmFscGhhRw==, and NiMlJWJldGFH, are positive. The corresponding differential equation isode := Diff( B, t ) = f[logistic] - f[predation];Via dimensional analysis, the problem can be reduced to one involving only two parameters. Introduce new, dimensionless, dependent and independent variables, NiMvJSJ5Ry1GJDYjJSR0YXVH, defined bynew_var_eq := B(t)=alpha*y(tau), t=alpha/beta*tau;and new (dimensionless) parameters NiMlIlJH and NiMlIlFH defined so thatnew_par_eq := r=beta/alpha*R, K=alpha*Q;The corresponding differential equation in terms of the new variables and parameters isode2 := dchange( {new_var_eq, new_par_eq}, subs(B=B(t), ode/beta), [y,tau,Q,R] );which can be simplified toode2 := map( collect, ode2, {R,alpha} );The advantage of this ODE is that it involves only two parameters.<Text-field bookmark="7.B" layout="Heading 1" style="Heading 1">7.B Bifurcation Analysis</Text-field>To begin the bifurcation analysis, it would be nice to be able to identify all equilibrium solutions for the dimensionless ODEode2;While it is easy to see that NiMvJSJ5RyIiIQ== is an equilibrium solution, no other equilibrium solution is immediately obvious.equil_eq := eval( rhs(ode2), y(tau)=y ) = 0:factor( equil_eq );It is apparent from the numerator that, in addition to the trivial equilibrium, there can be up to three additional equilibria. While Maple is capable of finding explicit formulas for the remaining equilibria, these expressions are so cumbersome that they are not likely to be very useful.solve( equil_eq, {y} );However, they do serve to show that the equation equil_eq;defines NiMvJSJ5Ry1GJDYkJSJRRyUiUkc= either implicitly, or explicitly. All equilibrium values NiMlInlH are a function of the two parameters NiMlIlFH and NiMlIlJH, and therefore lie on a surface above the NiMlI1FSRw==-plane. This surface is seen in Figure 7.1.implicitplot3d(equil_eq, Q=0..100, R=0..1, y=0..10, axes=box, grid=[20,20,20], title="Figure 7.1");Above a point in the NiMlI1FSRw==-plane where this surface folds over itself, there are multiple equilibria. Careful inspection of this surface shows there are regions in the NiMlI1FSRw==-plane where only one equlibrium solution exists, and regions where there are clearly three equilibrium solutions. Exactly along the leading edge of a fold in the surface, there will be just two equilibrium solutions. It is along such folds that the bifurcation points are found.The bifurcation analysis presented in Lesson 6, Section B can be used to identify bifurcation points in terms of the parameters NiMlIlFH and NiMlIlJH. The two necessary conditions arebif_eq1 := equil_eq;bif_eq2 := diff( eval(rhs(ode2),y(tau)=y), y ) = 0;The parametric solutions to these two equations arebif_sol := solve( {bif_eq1,bif_eq2}, {R,Q} );Since the original unknowns and parameters are positive, NiMlInlH, NiMlIlFH and NiMlIlJH should also be positive. Thus, there are physically realistic equilibria only when NiMqKCUieUciIiIlIj5HRiVGJUYl.Figure 7.2 shows the regions I and II in the first quadrant of the NiMlI1FSRw==-plane where there will be either one or three nontrivial equilibrium solutions, respectively.p1 := plot( eval([Q,R,y=1.001..100],bif_sol), labels=['Q','R'], view=[0..200,0..1], numpoints=200 ):p2 := textplot([[15,0.1,`Region I`],[50,0.7,`Region I`], [100,0.25,`Region II`]]):display([p1,p2], title="Figure 7.2");In Region I there is only one non-trivial equilibria while in Region II there are three non-trivial equilibria. At every point of the boundary between Regions I and II, i.e., the red curve, there are two non-trivial equilibria. To understand this conclusion, note that the equilibria must satisfyNiMvKiYlIlJHIiIiLCZGJkYmKiYlInlHRiYlIlFHISIiRitGJiomRilGJiwmRiZGJiokRikiIiNGJkYrThe function on the left-hand side is linear while the one on the right-hand side does not depend on the parameters. From a graph of these functions it is possible to see how the number of equilibria change with the parameter values.P := proc(R,Q,T) if nargs=2 then T := `` end if; plot( [R*(1-y/Q),y/(1+y^2)], y=0..10, title=T )end proc:display( array([P(0.5,5,`Region I`),P(0.5,7.3,`Boundary`),P(0.5,10,`Region II`)]) );To see the transition from three, to two, to one equilibrium, the animation in Figure 7.3 is helpful:display( seq( P(0.5, 10-q/10, sprintf("P= 0.5, Q=%4.1f",10.-q/10)), q=0..50), view=0..0.5, insequence=true, title="Figure 7.3" );For additional information about this problem -- both biological and mathematical -- please consult the references.<Text-field bookmark="7.C" layout="Heading 1" style="Heading 1">7.C References</Text-field>1. Ludwig, Jones, and Holling, ``Qualitative Analysis of Insect Breakout Systems: The Spruce Budworm and Forest,'' J. Animal Ecology (1978), 47, 315-332.2. Murray, J., Mathematical Biology, Springer-Verlag, 1993.3. Strogatz, S., Nonlinear Dynamics and Chaos, Addison-Wesley, 1994.4. McKelvey, S., Spruce Budworm Model, URL: http://www.stolaf.edu/people/mckelvey/envision.dir/spruce.html.[Back to ODE Powertool Table of Contents ]