CoCalc Public FilesLab 5 Group.sagews
Views : 172
##Problem one Romeo and Juliet

s=var("R,J")         # R and J are the variables
t=srange(0,20,0.05)   # we will use t for desolve_odeint.  It is the range of time values
# for which we will solve the differential equations.
ic1=[.5,.6]           # our 4 initial conditions.  Play around with different values for these.
ic2=[.5,-.2]
ic3=[-.9,.6]
ic4=[-.5,-.7]

@interact
def Romeo_Juliet(a=slider(-2,2,0.5,default=-1),
b=slider(-2,2,0.5,default=0),
c=slider(-2,2,0.5,default=0),
d=slider(-2,2,0.5,default=-1)):
# use 4 sliders for the coefficients, set the default values
# for a,d to -1, b and c to 0.
titleeq= ("$R'= "+'%2.1f'%(a)+" R "+"%+2.1f"%(b)+" J$\n"
+"$J' = "+'%2.1f'%(c)+" R "+"%+2.1f"%(d)+" J$")
# some fancy formatting stuff to create a title with our differential equations
vector_field=[a*R+b*J,c*R+d*J]
# our differential equations
pv=plot_vector_field(vector_field,(R,-1,1),(J,-1,1),plot_points=11,axes_labels=["R","J"])
# plot the vector field,R and J both going from -1 to 1, using 11x11 arrows
sol1=desolve_odeint(vector_field,times=t,ics=ic1,dvars=[R,J])
# solve the differential equation starting at ic1, for times in t
ps1=list_plot(sol1,plotjoined=True,color='blue')+disk(ic1,.03,(0,2*pi),color='blue')
# plot a trajectory starting at ic2, with a starting dot
sol2=desolve_odeint(vector_field,times=t,ics=ic2,dvars=[R,J])
# solve the differential equation starting at ic2, for times in t
ps2=list_plot(sol2,plotjoined=True,color='red')+disk(ic2,.03,(0,2*pi),color='red')
# plot a trajectory starting at ic2, with a starting dot
sol3=desolve_odeint(vector_field,times=t,ics=ic3,dvars=[R,J])
# solve the differential equation starting at ic3, for times in t
ps3=list_plot(sol3,plotjoined=True,color='green')+disk(ic3,.03,(0,2*pi),color='green')
#plot a trajectory starting at ic3, with a starting dot
sol4=desolve_odeint(vector_field,times=t,ics=ic4,dvars=[R,J])
# solve the differential equation starting at ic4, for times in t
ps4=list_plot(sol4,plotjoined=True,color='purple')+disk(ic4,.03,(0,2*pi),color='purple')
#a trajectory starting at ic2, with a starting dot
show(pv+ps1+ps2+ps3+ps4,xmin=-1,xmax=1,ymin=-1,ymax=1,aspect_ratio=1,
title=titleeq,figsize=8,fontsize=20)#"R'=aR+bJ\n J'=cR+dJ")
# combine and show the plots but limit the range of x and of y, add the title etc.












vs=var("D,M")
t=srange(0,50,0.05)
Dmax=3.1
Mmax=2.5
@interact
def Deer_Moose(r_D=slider(0,5,0.5,default=3,label="$r_D$"),
r_M=slider(0,5,0.5,default=2,label="$r_M$"),
alpha=slider(0,2,0.5,default=0.5,label="$\\alpha$"),
D0=slider(0,Dmax,0.1,default=1,label="$D_0$"),
M0=slider(0,Mmax,0.1,default=1,label="$M_0$")):
titleeq=("$D'= " + '%2.1f'%(r_D) + " D - M D - D^2$\n " +
"$M'= " + '%2.1f'%(r_M) + " M - " + '%2.1f'%(alpha) + " M D - M^2$")
vector_field=[r_D*D-M*D-D^2,r_M*M-alpha*M*D-M^2]
initial_value=[D0,M0]
pv=plot_vector_field(vector_field,(D,0,Dmax),(M,0,Mmax),plot_points=11,axes_labels=["D","M"])
sol1=desolve_odeint(vector_field,times=t,ics=initial_value,dvars=[D,M])
ps1=list_plot(sol1,plotjoined=True,color='green',thickness=2)
d1=disk(initial_value,.035,(0,2*pi),color='green')
nullclines=( plot(r_D-D,(D,0,r_D),color="red",thickness=3)
+plot(r_M-alpha*D,(D,0,r_M/alpha),color="blue",thickness=3))
show(pv+ps1+d1+nullclines,xmin=0,xmax=Dmax,ymin=0,ymax=Mmax,
aspect_ratio=1,title=titleeq,figsize=8,fontsize=20,svg=False)

def Scale_vf(xy):
# a function that takes  a vector input and returns a scaled
# version of the vector to keep arrows from getting too long
#
scale=sqrt(5+xy[0]^2+xy[1]^2)
return [xy[0]/scale,xy[1]/scale]


#red line is null line for D
#where they cross D' and M' are zero and that is an equilibrium point


86f69a6a-92be-4793-8680-37e8ab5a35︠











#Deer Moose norm





##I have made the comments on this so we can put it into the answer PDF when ready. Please do not delete! --Nicole

def Scale_vf(xy):
# a function that takes  a vector input and returns a scaled
# version of the vector to keep arrows from getting too long
#
scale=sqrt(5+xy[0]^2+xy[1]^2)
return [xy[0]/scale,xy[1]/scale]

vs=var("D,M")#tells the program what the names of the variables are
t=srange(0,50,0.05)#sets up "t" as a list of numbers from 0-50 in steps of 0.05, this will be used later int he code
Dmax=4.0 #sets up Dmax as a variable to be used later
Mmax=3.5 #sets up Mmax as a variable to be used later
@interact #starts the code as an interaction
def Deer_Moose_norm(r_D=slider(0,5,0.5,default=3,label="$r_D$"), #set up the 4 sliders in the interaction so that we can manipulate the values of the coefficients in the equations
r_M=slider(0,5,0.5,default=2,label="$r_M$"),
alpha=slider(0,2,0.5,default=0.5,label="$\\alpha$"),
D0=slider(0,Dmax,0.1,default=1,label="$D_0$"),
M0=slider(0,Mmax,0.1,default=1,label="$M_0$")):
titleeq=("$D'= " + '%3.1f'%(r_D) + " D - M D - D^2$\n " +
"$M'= " + '%3.1f'%(r_M) + " M - " + '%3.1f'%(alpha) + " M D - M^2$") #these two lines list the equations of the system as a title to the graph/vector field
vector_field=[r_D*D-M*D-D^2,r_M*M-alpha*M*D-M^2] #tells it to plot a vector field using the system equation. the coefficients in the original system of equations have been replaced with the terms used in the slider function
initial_value=[D0,M0] #tells the program what the initial values are for the starting point
pv=plot_vector_field(Scale_vf(vector_field),(D,0,Dmax),(M,0,Mmax), #puts the vector field into a named variable to be used later
plot_points=11,axes_labels=["D","M"])
# plot the scaled vector field
sol1=desolve_odeint(vector_field,times=t,ics=initial_value,dvars=[D,M]) #this next section of code is to set up the trajectory line
ps1=list_plot(sol1,plotjoined=True,color='green',thickness=2) #stores the trajectory line in to the variable ps1
d1=disk(initial_value,.035,(0,2*pi),color='green') #this line is just to put a dot ("disk") at the beginning of the trajectory line to make it easier to see
nullclines=(plot(r_D-D,(D,0,r_D),color="red",thickness=3) #these lines are to plot the nullclines of the system as red and blue lines on the vector graph
+plot(r_M-alpha*D,(D,0,r_M/alpha),color="blue",thickness=3))
show(pv+ps1+d1+nullclines,xmin=0,xmax=Dmax,ymin=0,ymax=Mmax,
aspect_ratio=1,title=titleeq,figsize=8,fontsize=20,svg=False)#tells it to show all the plots on one graph with various guidelines #square, no frame, size eight







##Shark tuna







def Scale_vf(xy):
# a function that takes  a vector input and returns a scaled
# version of the vector to keep arrows from getting too long
#
scale=sqrt(5+xy[0]^2+xy[1]^2)
return [xy[0]/scale,xy[1]/scale]

vs=var('T,S')
t=srange(0,150,0.05)
Tmax=100
Smax=40
@interact
def Shark_Tuna_norm(r_T=slider(0,0.4,0.01,default=0.1,label="$r_T$"),
r_S=slider(-0.5,0,0.1,default=-0.2,label="$r_S$"),
beta=slider(0,0.1,0.005,default=0.01,label="$\\beta$"),
m=slider(0,1,0.01,default=0.5,label="$m$"),
T0=slider(0,Tmax,0.1,default=20,label="$T_0$"),
S0=slider(0,Smax,0.1,default=10,label="$S_0$")):
titleeq=("$T'= "+'%4.3f'%(r_T)+" T - "'%4.3f'%(beta)+" S T$\n "
+"$S' = "+'%4.3f'%(r_S)+" S + "+'%4.3f'%(beta*m)+" S T$")
vector_field=[r_T*T-beta*T*S,r_S*S+m*beta*T*S]
initial_value=[T0,S0]
pv=plot_vector_field(Scale_vf(vector_field),(T,0,Tmax),
(S,0,Smax),plot_points=11,axes_labels=["T","S"])
sol1=desolve_odeint(vector_field,times=t,ics=initial_value,dvars=[T,S])
ps1=list_plot(sol1,plotjoined=True,color='green',thickness=2)
d1=disk(initial_value,1,(0,2*pi),color='green')
show(pv+ps1+d1,xmin=0,xmax=Tmax,ymin=0,ymax=Smax,
aspect_ratio=1,title=titleeq,figsize=8,fontsize=20,svg=False)




def Scale_vf(xy):
# a function that takes  a vector input and returns a scaled
# version of the vector to keep arrows from getting too long
#
scale=sqrt(5+xy[0]^2+xy[1]^2)
return [xy[0]/scale,xy[1]/scale]

vs=var('T,S') #tells the program we will be using S and T as the variables
t=srange(0,150,0.05) #sets up "t" as a list from 0 to 150 in steps of 0.05. This will be used later to make a trajectory line
Tmax=100 #setting up variables Tmax and Smax. These will be used to define the window size of the graph
Smax=40
@interact #starts the interaction program
def Shark_Tuna_norm(r_T=slider(0,0.4,0.01,default=0.1,label="$r_T$"), #these lines set up each of the sliders which will allow us to change the values of the coefficeints in the system of equations
r_S=slider(-0.5,0,0.1,default=-0.2,label="$r_S$"),
beta=slider(0,0.1,0.005,default=0.01,label="$\\beta$"),
m=slider(0,1,0.01,default=0.5,label="$m$"),
T0=slider(0,Tmax,0.1,default=20,label="$T_0$"),
S0=slider(0,Smax,0.1,default=10,label="$S_0$")):
titleeq=("$T'= "+'%4.3f'%(r_T)+" T - "'%4.3f'%(beta)+" S T$\n "
+"$S' = "+'%4.3f'%(r_S)+" S + "+'%4.3f'%(beta*m)+" S T$") #these two lines set up the title of the graph to display the system of equations with the current coefficient values displayed (based on what value the slider is at)
vector_field=[r_T*T-beta*T*S-0.001*T*T,r_S*S+m*beta*T*S] #this line sets up the vector field with the T prime and S prime equations. Notice the new added T squared term for tuna competition
initial_value=[T0,S0] #setting initial values for T and S
pv=plot_vector_field(Scale_vf(vector_field),(T,0,Tmax),
(S,0,Smax),plot_points=11,axes_labels=["T","S"]) #these lines set up the rest of the vector field and store it as a variable "pv" and define figure size, labels etc
sol1=desolve_odeint(vector_field,times=t,ics=initial_value,dvars=[T,S]) #this line sets up for plotting the trajectory line
ps1=list_plot(sol1,plotjoined=True,color='green',thickness=2) #stores the trajectroy line as a variable "ps1"
d1=disk(initial_value,1,(0,2*pi),color='green') #sets up a dot ("disk") to mark the start of the trajecctory line
show(pv+ps1+d1,xmin=0,xmax=Tmax,ymin=0,ymax=Smax,
aspect_ratio=1,title=titleeq,figsize=8,fontsize=20,svg=False)#these lines set up the system to display all the plots on one graph and make the figure square, add a nice title and show the size