CoCalc -- Collaborative Calculation in the Cloud
SharedLab 10 / Lab 10-turnin.sagewsOpen in CoCalc
# Lab 9

# Name: Angelica Jensen
# I worked on this code with:

# 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.
#1
#Xprime represents plants
#Yprime represents herbivores
#Zpime represents carnivores

#2
var("X,Y,Z")
a1= 5
b1= 3
a2= 0.1
b2= 2
d1=0.4
d2=0.01
f1= X*(1-X) -(a1*X)/ (1+b1*X)*Y
f2= (a1*X)/(1+b1*X)*Y -d1*Y - (a2*Y)/(1+b2*Y)*Z
f3 = (a2*Y)/(1+b2*Y)*Z- d2*Z
t = srange(0,1000)
sol = desolve_odeint([f1,f2,f3], [1,0.5,8], t, [X,Y,Z])
plot = list_plot(sol , plotjoined = True, color="hotpink")
show(plot)
(X, Y, Z)
3D rendering not yet implemented

p=list_plot(zip(t,sol[:,0]),plotjoined=True, color="violet",axes_labels=["X","f(X)"])+list_plot(zip(t,sol[:,1]),plotjoined=True, color="gold"),list_plot(zip(t,sol[:,2]),plotjoined=True, color="pink",axes_labels=["X","f(X)"])
show(p)
()
#3
var("X,Y,Z")
a1= 5
b1= 3
a2= 0.1
b2= 2
d1=0.4
d2=0.01
f1= X*(1-X) -(a1*X)/ (1+b1*X)*Y
f2= (a1*X)/(1+b1*X)*Y -d1*Y - (a2*Y)/(1+b2*Y)*Z
f3 = (a2*Y)/(1+b2*Y)*Z- d2*Z
t = srange(0,1000)
sol1 = desolve_odeint([f1,f2,f3], [2,1,12], t, [X,Y,Z])
plot = list_plot(sol1 , plotjoined = True, color="hotpink") + list_plot(sol , plotjoined = True, color="turquoise")
show(plot)
(X, Y, Z)
3D rendering not yet implemented
p=list_plot(zip(t,sol[:,0]),plotjoined=True, color="violet",axes_labels=["X","f(X)"])+list_plot(zip(t,sol[:,1]),plotjoined=True, color="gold"),list_plot(zip(t,sol[:,2]),plotjoined=True, color="pink")+list_plot(zip(t,sol1[:,0]),plotjoined=True, color="purple")+list_plot(zip(t,sol1[:,1]),plotjoined=True, color="red") +list_plot(zip(t,sol1[:,2]),plotjoined=True, color="deeppink",axes_labels=["X","f(X)"])
show(p)
()
#4
#The first (top) time series looks the same however the bottom time series changes as the inital conditions change

#5
@interact
def hpg(b1=(2,3,0.1)):
    var("X,Y,Z")
    a1= 5
    a2= 0.1
    b2= 2
    d1=0.4
    d2=0.01
    f1= X*(1-X) -(a1*X)/ (1+b1*X)*Y
    f2= (a1*X)/(1+b1*X)*Y -d1*Y - (a2*Y)/(1+b2*Y)*Z
    f3 = (a2*Y)/(1+b2*Y)*Z- d2*Z
    t = srange(0,1000)
    sol1 = desolve_odeint([f1,f2,f3], [2,1,12], t, [X,Y,Z])
    plot = list_plot(sol1 , plotjoined = True, color="hotpink") + list_plot(sol , plotjoined = True, color="turquoise")
    show(plot)
Interact: please open in CoCalc
#6
#changing the b1 value from 2 to 3 shifts the trajectory up and down

#7
f(X)=2.5*X*(1-X)
vals = [.5]
for i in srange(0,20):
    a= vals[i]
    b= f(a)
    vals.append(b)
vals
[0.500000000000000, 0.625000000000000, 0.585937500000000, 0.606536865234375, 0.596624740865082, 0.601659148631890, 0.599163543748598, 0.600416478978050, 0.599791326874127, 0.600104227701753, 0.599947858990589, 0.600026063707993, 0.599986966447711, 0.600006516351461, 0.599996741718113, 0.600001629114403, 0.599999185436164, 0.600000407280259, 0.599999796359456, 0.600000101820168, 0.599999949089890]
list_plot(vals, plotjoined=True,axes_labels=["X","f(X)"])
f(X)=3.1*X*(1-X)
vals = [0.5]
for i in srange(0,20):
    a= vals[i]
    b= f(a)
    vals.append(b)
vals
[0.500000000000000, 0.775000000000000, 0.540562500000000, 0.769899519140625, 0.549178173659744, 0.767502672430025, 0.553171192752663, 0.766235755209903, 0.555267420208218, 0.765531088016937, 0.556429048019278, 0.765128863872878, 0.557090725178579, 0.764896012205610, 0.557473318424452, 0.764760134774773, 0.557696420205537, 0.764680481595945, 0.557827152252630, 0.764633663433289, 0.557903974951418]

list_plot(vals, plotjoined=True,axes_labels=["X","f(X)"])
f(X)=3.5*X*(1-X)
vals = [0.5]
for i in srange(0,20):
    a= vals[i]
    b= f(a)
    vals.append(b)
vals
[0.500000000000000, 0.875000000000000, 0.382812500000000, 0.826934814453125, 0.500897694844753, 0.874997179503880, 0.382819903774472, 0.826940887670016, 0.500883795893397, 0.874997266166866, 0.382819676285819, 0.826940701069839, 0.500884222943868, 0.874997263524249, 0.382819683222636, 0.826940706759848, 0.500884209921798, 0.874997263604850, 0.382819683011062, 0.826940706586302, 0.500884210318974]
list_plot(vals, plotjoined=True,axes_labels=["X","f(X)"])
f(X)=3.9*X*(1-X)
vals = [0.5]
for i in srange(0,20):
    a= vals[i]
    b= f(a)
    vals.append(b)
vals
[0.500000000000000, 0.975000000000000, 0.0950625000000001, 0.335499922265625, 0.869464925259000, 0.442633109113109, 0.962165255336889, 0.141972779361614, 0.475084386199614, 0.972578927536905, 0.104009713267468, 0.363447601972599, 0.902278426112569, 0.343871064749139, 0.879932646751986, 0.412039617334918, 0.944825587217508, 0.203307768130735, 0.631697506238894, 0.907357490716780, 0.327833511552022]
list_plot(vals, plotjoined=True,axes_labels=["X","f(X)"])
#This r condition at 3.9 results in a chaotic behavior

#8
f(X)=3.9*X*(1-X)
valsa = [0.4]
for i in srange(0,20):
    a= valsa[i]
    b= f(a)
    valsa.append(b)
valsa
[0.400000000000000, 0.936000000000000, 0.233625600000000, 0.698274248196096, 0.821680557758864, 0.571434313163790, 0.955098841720988, 0.167251672630438, 0.543186347467759, 0.967726263630336, 0.121805355010580, 0.417178360955174, 0.948248246813121, 0.191386685992957, 0.603555507428601, 0.933177401836699, 0.243193620298226, 0.717796885043427, 0.790001615774983, 0.647006345106807, 0.890717624543520]

f(X)=3.9*X*(1-X)
valsb = [0.401]
for i in srange(0,20):
    a= valsb[i]
    b= f(a)
    valsb.append(b)
valsb
[0.401000000000000, 0.936776100000000, 0.230983890028281, 0.692758297045167, 0.830092531788964, 0.550051789782892, 0.965229791524164, 0.130888840205693, 0.443652111689677, 0.962617170383573, 0.140342879298556, 0.470522346562113, 0.971611164996404, 0.107573345297866, 0.374405150649239, 0.913481141884582, 0.308230046692706, 0.831574711533221, 0.546227021619479, 0.966665943641550, 0.125669298476054]
list_plot(vals, plotjoined=True) +list_plot(vals1, plotjoined=True, color="red",axes_labels=["X","f(X)"])
#9
f(X)=3.7*X*(1-X)
vals1 = [0.5]
for i in srange(0,20):
    a= vals1[i]
    b= f(a)
    vals1.append(b)
vals1
[0.500000000000000, 0.925000000000000, 0.256687500000000, 0.705956401171875, 0.768053255020420, 0.659145574149943, 0.831288939045395, 0.518916263804854, 0.923676047365561, 0.260844845488171, 0.713377804660565, 0.756538676169480, 0.681495258228080, 0.803120043590673, 0.585037484942278, 0.898243916772360, 0.338186596189095, 0.828120762684378, 0.526646030853062, 0.922372959447177, 0.264924007572982]
f(X)=3.8*X*(1-X)
vals2 = [0.5]
for i in srange(0,20):
    a= vals2[i]
    b= f(a)
    vals2.append(b)
vals2
[0.500000000000000, 0.950000000000000, 0.180500000000000, 0.562095050000000, 0.935347978108890, 0.229794124234705, 0.672557381867257, 0.836851009859847, 0.518819309194325, 0.948654167685504, 0.185095863710026, 0.573174462800368, 0.929652892376735, 0.248513889874763, 0.709667998373493, 0.782949455740601, 0.645770500885172, 0.869253652072407, 0.431876613638521, 0.932364976076450, 0.239630004357161]
f(X)=3.9*X*(1-X)
vals3 = [0.5]
for i in srange(0,20):
    a= vals3[i]
    b= f(a)
    vals3.append(b)
vals3
[0.500000000000000, 0.975000000000000, 0.0950625000000001, 0.335499922265625, 0.869464925259000, 0.442633109113109, 0.962165255336889, 0.141972779361614, 0.475084386199614, 0.972578927536905, 0.104009713267468, 0.363447601972599, 0.902278426112569, 0.343871064749139, 0.879932646751986, 0.412039617334918, 0.944825587217508, 0.203307768130735, 0.631697506238894, 0.907357490716780, 0.327833511552022]
list_plot(vals1, plotjoined=True)+list_plot(vals2, plotjoined=True, color="red")+list_plot(vals3, plotjoined=True, color="green",axes_labels=["X","f(X)"])
#There is no relationship, changing the initial condition will completely change the behavior

#10
emptylist=[]
f(X)=3.9*X*(1-X)
for i in srange(0,8):
    a= abs(valsa[i]-valsb[i])
    emptylist.append(a)
emptylist
[0.00100000000000000, 0.000776100000000168, 0.00264170997171956, 0.00551595115092951, 0.00841197403010030, 0.0213825233808989, 0.0101309498031757, 0.0363628324247447]
list_plot(emptylist, plotjoined=True, axes_labels=["X","f(X)"])
#I see a plot that starts off increasing slightly then showing chaotic behavior when X is 4