Mathematical Principles for Scientific Computing and Visualization
http : // www.farinhansford.com/books/scv
by Gerald Farin and Dianne Hansford

This notebook demonstrates the predator-prey model where two ODEs describe two populations. This particular example is called the fox and rabbit model.

Figures 9.8 and  9.9 were created with this notebook.
(Errata: The figure was generated, as you will find in this notebook,  with equation (9.12) omitting the 0.5 factor.)

See Chapter 9 : Computing Dynamic Processes for details.

Updated July 2008

(* Load the vectorfield package *)

Needs["VectorFieldPlots`"];

(* Specify the output folder *)

SetDirectory[ToFileName[Extract["FileName"/.NotebookInformation[EvaluationNotebook[]],{1},FrontEnd`FileName]]];

(* Create a vectorfield *)
(* The first equation describes rabbit growth, so the x-axis will represent the number of rabbits. The y-axis represents the number of foxes *)
vecfield=VectorFieldPlot[{2 r-r f,-f+0.9` r f},{r,0,6},{f,0,7},PlotPoints->20,Frame->True,BaseStyle->{"FontFamily"->"Helvetica","FontSize"->11}];

(* If you want the vectorfield plotted on its own .. *)
(* This is Figure 9.8 *)
(* Export["rabfoxfield.eps",vecfield]; *)

(* Solve the dynamical system *)
(* Given: two governing ODEs, starting values, and a domain for t *)
(* We want the solution over time = 0 to time=10 *)
(* Solution: comes in definitions for r(t) and f(t) *)
rabfox1=NDSolve[{r'[t]==2 r[t]-r[t] f[t],
                                f'[t]==-f[t]+0.9` r[t] f[t],
                                r[0]==5,f[0]==2},
                                 {r[t],f[t]},{t,0,10}]

(* This code loads the InterpolationFunction objects into  rab1 and fox1. A shorter way to do this is below *)
(*
rab1=rabfox1[[1,1,2]]
fox1=rabfox1[[1,2,2]]
rabfoxplot1=ParametricPlot[{rab1,fox1},{t,0,10},PlotStyle->{Thickness[0.02`],GrayLevel[0.4`]},BaseStyle->{"FontFamily"->"Helvetica","FontSize"->11}];
*)

(* r(t) and f(t) are defined in rabfox1. Use the ReplaceAll command to get them for plotting *)
rabfoxplot1=ParametricPlot[Evaluate[{r[t],f[t]}/. rabfox1],{t,0,10},PlotStyle->{Thickness[0.015], RGBColor[0.3, 0.3, 0.8]},BaseStyle->{"FontFamily"->"Helvetica","FontSize"->11}];

(* Create a plot of the parametric plot *)
points=Graphics[rabfoxplot1];

(* Create the start point -- five rabbits and 2 foxes *)
startpoint=ListPlot[{{5,2}},Axes->False,PlotStyle->PointSize[0.05`],BaseStyle->{"FontFamily"->"Helvetica","FontSize"->11}];

(* Plot it all together *)
pic=Show[vecfield,points,startpoint]

Export["foxrab.eps",pic]