CoCalc Public FilesLab 4 / Lab4-turnin.ipynb
Author: Nushrat Esha
Views : 142
Compute Environment: Ubuntu 20.04 (Default)
In [ ]:
# Name: Nushrat Esha
# I worked on this code with: Xavier Herrera, Fujia Guo

# Please do all of your work for this week's lab in this worksheet. If
# you wish to create other worksheets for scratch work, you can, but
# this is the one that will be graded. You do not need to do anything
# to turn in your lab. It will be collected by your TA at the beginning
# of (or right before) next week’s lab.

# Be sure to clearly label which question you are answering as you go and to
# use enough comments that you and the grader can understand your code.


In [ ]:
#1
# X represents the plants because X's differential equations has one outflow which represents the predation rate.
# Y represents the herbivores because there is a natural death rate outflow and an outflow from Z, the carnivore, population.
# Z represents the carnivores because there is no predation outflow in this differential equation and it contains a natural death rate outflow.

In [18]:
#2 time series
var("X", "Y", "Z")
a1 = 5
b1 = 3
a2 = 0.1
b2 = 2
d1 = 0.4
d2 = 0.01

X0 = 1
Y0 = 0.5
Z0 = 8

Xprime = X*(1-X)-((a1*X)/(1+b1*X))*Y
Yprime = ((a1*X)/(1+b1*X))*Y - d1*Y - ((a2*Y)/(1+b2*Y))*Z
Zprime = ((a2*Y)/(1+b2*Y))*Z - d2*Z

t = srange(0,1000,0.1)
sol1 = desolve_odeint([Xprime, Yprime, Zprime], ics = [X0, Y0, Z0], times = t, dvars = [X, Y, Z])

timeseries1 = list_plot(list(zip(t, sol1[:, 0])), plotjoined = True, color = "magenta", legend_label = "plants") + list_plot(list(zip(t, sol1[:, 1])), plotjoined = True, color = "blue", legend_label = "herbivores") + list_plot(list(zip(t, sol1[:, 2])), plotjoined = True, color = "green", legend_label = "carnivores", axes_labels = ["Time", "State Variable"])
timeseries1

In [21]:
#2 trajectory
trajectory1 = list_plot(sol1, plotjoined = True, axes_labels = ["plants", "herbivores", "carnivores"])
trajectory1

In [19]:
#3 time series
var("X", "Y", "Z")
a1 = 5
b1 = 3
a2 = 0.1
b2 = 2
d1 = 0.4
d2 = 0.01

X0 = 2
Y0 = 2.5
Z0 = 10

Xprime = X*(1-X)-((a1*X)/(1+b1*X))*Y
Yprime = ((a1*X)/(1+b1*X))*Y - d1*Y - ((a2*Y)/(1+b2*Y))*Z
Zprime = ((a2*Y)/(1+b2*Y))*Z - d2*Z

t = srange(0,1000,0.1)
sol2 = desolve_odeint([Xprime, Yprime, Zprime], ics = [X0, Y0, Z0], times = t, dvars = [X, Y, Z])

timeseries2 = list_plot(list(zip(t, sol2[:, 0])), plotjoined = True, color = "red", legend_label = "plants2") + list_plot(list(zip(t, sol2[:, 1])), plotjoined = True, color = "purple", legend_label = "herbivores2") + list_plot(list(zip(t, sol2[:, 2])), plotjoined = True, color = "orange", legend_label = "carnivores2", axes_labels = ["Time", "State Variable"])
timeseries1 + timeseries2

In [22]:
#3 trajectory
trajectory2 = list_plot(sol2, plotjoined = True, axes_labels = ["plants", "herbivores", "carnivores"], color = "orange")
trajectory1 + trajectory2