Callahan
    8.1: 1(a)(b)(c), 2, 5(a)(b)(c)

1)  Predator-prey models:

    a)  The original Lotka-Volterra model:

    (1) Why is it a reasonable model and what is it saying about how the populations change?

    The equation:
    The rate of change of the rabbit population is proportional to the rabbit population, and to the negative product of the rabbit and fox populations.  We can think of the product term as an interaction term, where the frequency with which foxes encounters and eat rabbits is proportional to how big both the rabbit and fox populations.  This makes since because the more animals their are, the more encounters there will be.

    The equation:
    The rate of change of the fox population is proportional to the negative of the fox population, and to the product of the rabbit and fox populations.  The product term works the same as in the rabbit rate equation, except here it has a positive effect on the rate of change of the fox population.  The rate of change of the fox population being proportional to the negative of the fox population models the fact that in the absence of rabbits, the fox population can not survive and so must have a negative rate of change.

    (2) Plot (in red) the set of points where , and mark the regions where and .

    (3) Plot (in green) the set of points where , and mark the regions where and .

    (4) Mark the equilibrium points.  What color are they?

    (5) Plot the vector field and a few trajectories.    (p412)

    (6) On the basis of your plots make a conjecture about the stability of the equilibrium points.

    The equilibrium point is totally stable (center), all trajectories are stable solutions, they don't change from cycle to cycle.  We know this because as we increase the number of steps used in Euler's method, the curve gets closer and closer to returning to where it started after one cycle.  The equilibrium point is unstable in the rabbit axis, but stable in the fox axis.

    b)  The Leslie-Gower model:

    (1) Why is it a reasonable model and what is it saying about how the populations change?

    The equation:
    This is the same as for the Lotka-Volterra model.  The rate of change of the rabbit population is proportional to the rabbit population, and to the negative product of the rabbit and fox populations.  We can think of the product term as an interaction term, where the frequency with which foxes encounters and eat rabbits is proportional to how big both the rabbit and fox populations are in the given region.  This makes since because the more animals their are, the more encounters there will be.

    The equation:
    The rate of change of the fox population is proportional to the fox population scaled by a carrying capacity like term.  In the carrying capacity term, the rabbit population is proportional to the carrying capacity; the constants and act to scale the carrying capacity.  It makes sense to have the carrying capacity be proportional to the rabbit population because the fox population is dependent of the rabbits as a food source.

    (2) Plot (in red) the set of points where , and mark the regions where and .

    (3) Plot (in green) the set of points where , and mark the regions where and .

    (4) Mark the equilibrium points.  What color are they?

    The equilibrium points occur anywhere the and lines cross:

    To find the coordinate of an equilibrium point, we set the equations of a line equal to a line:

    So:

    And we already had the coordinate of the equilibrium point:

    So our equilibrium points are:

    (5) Plot the vector field and a few trajectories.

    (6) On the basis of your plots make a conjecture about the stability of the equilibrium points.

    The equilibrium point is totally stable (spiral sink).  The equilibrium point is unstable in the rabbit axis, but stable in the fox axis.

    c)  The Leslie-Gower model with carrying capacity for rabbits:

    (1) Why is it a reasonable model and what is it saying about how the populations change?

    The equation:
    The rate of change of the rabbit population is proportional to the rabbit population scaled by a carrying capacity term, and to the negative product of the rabbit and fox populations.  We can think of the product term as an interaction term, where the frequency with which foxes encounters and eat rabbits is proportional to how big both the rabbit and fox populations are in the given region.  This makes since because the more animals their are, the more encounters there will be.  The carrying capacity term simply equals zero when , is greater than zero when , and is less than zero when .

    The equation:
    The rate of change of the fox population is proportional to the fox population scaled by a carrying capacity like term.  In the carrying capacity term, the rabbit population is proportional to the carrying capacity; the constants and act to scale the carrying capacity.  It makes sense to have the carrying capacity be proportional to the rabbit population because the fox population is dependent of the rabbits as a food source.

    (2) Plot (in red) the set of points where , and mark the regions where and .

    (3) Plot (in green) the set of points where , and mark the regions where and .

    (4) Mark the equilibrium points.  What color are they?

    The equilibrium points occur anywhere the and lines cross:

    To find the coordinate of an equilibrium point, we set the equations of a line equal to a line:

    So:

    And it is easy to find the coordinate:

    Also we have a equilibrium point on the axes, where the line intersects it:

    So:

    So our equilibrium points are:

    (5) Plot the vector field and a few trajectories.

    (6) On the basis of your plots make a conjecture about the stability of the equilibrium points.

    The equilibrium point is totally stable (spiral sink).  The equilibrium point is unstable in the rabbit axis, but stable in the fox axis.  The equilibrium point is stable in the rabbit axis, but unstable in the fox axis.

from visual.graph import *
from visual import *
field = display(background=color.white)
graph1 = gdisplay()

a = 0.1
b = 0.005
c = 0.00004
d = 0.04
e = 0.001
f = 0.05
g = 0.04
h = 0.004
K = 10000.0


# Create the vector field (auto scaled):    
arrowsize = 3.5
dt = arrowsize
xmax = 2000.0
xmin = 100.0
rangeX = xmax - xmin
ymax = 23.0
ymin = 15.0
rangeY = ymax - ymin
z = 0

xstep = rangeX/10.0
ystep = rangeY/10.0

if rangeX > rangeY:
    normalizeX = rangeX/rangeY
else:
    normalizeX = 1

if rangeY > rangeX:
    normalizeY = rangeY/rangeX
else:
    normalizeY = 1

field.center=((xmax-(rangeX/2.0))/normalizeX,(ymax-(rangeY/2.0))/normalizeY,0)
field=[]
for x in arange(xmin,xmax,xstep):
    for y in arange(ymin,ymax,ystep):
        Xprime = a*x*(1-(x/K))-b*x*y
        Yprime = (e-f*(y/x))*y
        
        deltaX = Xprime * dt
        deltaY = Yprime * dt
        
        vec=arrow(pos=(x/normalizeX,y/normalizeY,z),axis=(deltaX/normalizeX,deltaY/normalizeY,0),shaftwidth=arrowsize/5.0,color=(0,0,0))
        field.append(vec)

# Plot the function:    
tinitial = 0
tfinal = 2000

t = tinitial
x = 1000
y = 22

numberofsteps = 10000
deltat = (tfinal - tinitial) / float(numberofsteps)

Xcurve = gcurve(gdisplay=graph1)
Ycurve = gcurve(gdisplay=graph1)
XYcurve = curve(color=(0,0,0))

for k in range(numberofsteps):
    Xprime = a*x*(1-(x/K))-b*x*y
    Yprime = (e-f*(y/x))*y
    
    deltaX = Xprime * deltat
    deltaY = Yprime * deltat
    
    t = t + deltat
    x = x + deltaX
    y = y + deltaY
    
    Xcurve.plot(pos=(t,x))
    Ycurve.plot(pos=(t,y))
    XYcurve.append(pos=(x/normalizeX,y/normalizeY,0))

# Plot the X'=0 lines(RED):
numberofsteps = 2

Rzero = curve(color=(1,0,0))
yinitial = 0
yfinal = ymax
y = yinitial
deltay = (yfinal - yinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    x = 0
    Rzero.append(pos=(x/normalizeX,y/normalizeY,0))
    y = y + deltay

xinitial = 0
xfinal = xmax
x = xinitial
deltax = (xfinal - xinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    y = (a/b)-(a/(b*K))*x
    Rzero.append(pos=(x/normalizeX,y/normalizeY,0))
    x = x + deltax

# Plot the Y'=0 lines(GREEN):
Fzero = curve(color=(0,1,0))
xinitial = 0
xfinal = xmax
x = xinitial
deltax = (xfinal - xinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    y = 0
    Fzero.append(pos=(x/normalizeX,y/normalizeY,0))
    x = x + deltax

yinitial = 0
yfinal = ymax
y = yinitial
deltay = (yfinal - yinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    x = (f/e)*y
    Fzero.append(pos=(x/normalizeX,y/normalizeY,0))
    y = y + deltay

5(a)(b)(c):  More on Lotka-Volterra model

    a) Sketch the vector field, together with some typical trajectories, in the rest of the R-F plane, include negative values.  What happens to any trajectory starting at a state with a negative R or F value?

    Any trajectory that starts with a negative R or F value diverges to infinity, these solutions are not phisically meaningful

    b) Using Euler's method draw the first 500 segments of the trajectory, use the values of a, b, c, and d given above, start at the point , and use a time step of .

    Does the trajectory look like a closed loop?

    No.

    How small does have to be before the trajectory looks like it closes?  Can you explain the phenomenon?

    In order for it to look like a closed loop we must have of about 0.5.  The problem is due to the inherent error involved in the process of Euler's method.

    c) Using Euler's method draw the first 1000 segments of the trajectory, use the values of a, b, c, and d given above, start at the point , and use a time step of .

    What happens?

    The computer gives the error "Result too large" because the solution goes to infinity very quickly.

2)  Symbiosis and mutualism:  In this model is the number of bees per acre, measured in hundreds of bees, while C is the weight of clover per acre, in thousands of pounds.  Time is measured in months.

    a)  Do these equations describe symbiosis?

    Yes.  The interaction terms are both positive.

    b)  Each equation has a negative term in it.  What aspect of reality is this term capturing?

    Rewriting each equation as above, we can see that the negative term is part of the "carrying capacity".  If the population is above the carrying capacity, then some will die off because their is insufficient food to sustain that great of a population.

    c)  Sketch the vector field for this system in the plane, mark the equilibrium points.

    The set of points where :

    The set of points where :

from visual.graph import *
from visual import *
field = display(background=color.white)
graph1 = gdisplay()

a = 0.1
b = 0.0005
c = 0.03
d = 0.0012
Kb = 100
Kc = 10


# Create the vector field (auto scaled):    
arrowsize = 1.5
dt = arrowsize
xmax = 200.0
xmin = 0
rangeX = xmax - xmin
ymax = 120.0
ymin = 0
rangeY = ymax - ymin
z = 0

xstep = rangeX/10.0
ystep = rangeY/10.0

if rangeX > rangeY:
    normalizeX = rangeX/rangeY
else:
    normalizeX = 1

if rangeY > rangeX:
    normalizeY = rangeY/rangeX
else:
    normalizeY = 1

field.center=((xmax-(rangeX/2.0))/normalizeX,(ymax-(rangeY/2.0))/normalizeY,0)
field=[]
for x in arange(xmin,xmax,xstep):
    for y in arange(ymin,ymax,ystep):
        Xprime = a*x*(1-(x/Kb))+b*x*y
        Yprime = c*y*(1-(y/Kc))+d*x*y
        
        deltaX = Xprime * dt
        deltaY = Yprime * dt
        
        vec=arrow(pos=(x/normalizeX,y/normalizeY,z),axis=(deltaX/normalizeX,deltaY/normalizeY,0),shaftwidth=arrowsize/5.0,color=(0,0,0))
        field.append(vec)

# Plot the function:    
tinitial = 0
tfinal = 100

start = [5,10,50,120,90,5,200,5,200,50,150,120]

for i in range(len(start)/2):
    t = tinitial
    x = start[0+(2*i)]
    y = start[1+(2*i)]
    print x, y

    numberofsteps = 10000
    deltat = (tfinal - tinitial) / float(numberofsteps)

    Xcurve = gcurve(gdisplay=graph1)
    Ycurve = gcurve(gdisplay=graph1)
    XYcurve = curve(color=(0,0,0))

    for k in range(numberofsteps):
        Xprime = a*x*(1-(x/Kb))+b*x*y
        Yprime = c*y*(1-(y/Kc))+d*x*y
        
        deltaX = Xprime * deltat
        deltaY = Yprime * deltat
        
        t = t + deltat
        x = x + deltaX
        y = y + deltaY
        
        Xcurve.plot(pos=(t,x))
        Ycurve.plot(pos=(t,y))
        XYcurve.append(pos=(x/normalizeX,y/normalizeY,0))

# Plot the X'=0 lines(RED):
numberofsteps = 10

Rzero = curve(color=(1,0,0))
yinitial = 0
yfinal = ymax
y = yinitial
deltay = (yfinal - yinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    x = 0
    Rzero.append(pos=(x/normalizeX,y/normalizeY,0))
    y = y + deltay

xinitial = 0
xfinal = xmax
x = xinitial
deltax = (xfinal - xinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    y = (a/(b*Kb))*x-(a/b)
    Rzero.append(pos=(x/normalizeX,y/normalizeY,0))
    x = x + deltax

# Plot the Y'=0 lines(GREEN):
Fzero = curve(color=(0,1,0))
xinitial = 0
xfinal = xmax
x = xinitial
deltax = (xfinal - xinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    y = 0
    Fzero.append(pos=(x/normalizeX,y/normalizeY,0))
    x = x + deltax

yinitial = 0
yfinal = ymax
y = yinitial
deltay = (yfinal - yinitial)/float(numberofsteps)
for k in range(0,numberofsteps+1) :
    x = (c/(d*Kc))*y-(c/d)
    Fzero.append(pos=(x/normalizeX,y/normalizeY,0))
    y = y + deltay

    d)  Determine the stability of the equilibrium points:

    The equilibrium points occur anywhere the and lines cross:

    To find the coordinate of an equilibrium point, we set the equations of a line equal to a line:

    So:

    And it is easy to find the coordinate:

    Also we have a equilibrium point on the axes, where the line intersects it:

    So:

    And we have a equilibrium point on the axes, where the line intersects it:

    So:

    So our equilibrium points are:

    The equilibrium point is totally stable (sink).  All other equilibrium points tend toward this equilibrium point.

    e)  Suppose an acre of land has pounds of clover on it, and a hive of bees is introduced.  What happens?

    f)

    g)