SharedLab 10.sagewsOpen in CoCalc
#6
f(x,y)=x^2+y^2
f(y)=1+y^2
f(x)=x^2+1

dfdx = diff(f(x),x)
k1 = dfdx.subs(x=1)

dfdy = diff(f(y),y)
k2 = dfdy.subs(y=-1)

g(x,y) = k1*x+(k2)*y
z=k1*(x-1)+(k2)*(y+4)+2

plot3d(f(x,y), (x,-4,4), (y,-4,4))+plot3d(z, (x,-4,4), (y,-4,4), color="red")

3D rendering not yet implemented
#8
@interact
def Zoom(a=(0.5,2,0.5),b=(0.5,2,0.5)):
var("x,y")
t=srange(0,100,0.1)
r = 0.2
x_prime = r*x*(1-(x/7))-(x/(1+x))*y
y_prime = (1/10)*y*(1-(y/x))
sol=desolve_odeint([x_prime,y_prime],ics=[0.5,0.5],times=t,dvars=[x,y])
p=plot_vector_field([x_prime,y_prime], (x,0,a), (y,0,b))+point([0.239,0.24],size=50)+list_plot(sol,color="red",plotjoined=True)+point([.5,.5],color="red",size=30)
show(p)

#The equilibrium point in this model is a stable spiral

#9
f1 = 0.2*x*(1-(x/7))-(x/(1+x))*y
f2 = (1/10)*y*(1-(y/x))

Jacobian_matrix = jacobian([f1,f2], [x,y])
jacobian([f1,f2], [x,y])

[-0.0571428571428571*x - y/(x + 1) + x*y/(x + 1)^2 + 0.200000000000000 -x/(x + 1)] [ 1/10*y^2/x^2 -1/5*y/x + 1/10]
#10
Jacobian_matrix.subs(x=0.239, y=0.24)

#Jacobian = matrix(RDF, [[0.0300033894396597,-0.192897497982244],[0.100838570753313,-0.100836820083682]])
Jacobian_matrix.subs(x=0.239, y=0.24).eigenvectors_right()#Eigenvalues of -1/27240646557740*I*sqrt(11258216119645872304850111) - 2336015071/65957981980 and 1/27240646557740*I*sqrt(11258216119645872304850111) - 2336015071/65957981980

[0.0300033894396597 -0.192897497982244] [ 0.100838570753313 -0.100836820083682] [(-1/27240646557740*I*sqrt(11258216119645872304850111) - 2336015071/65957981980, [], 1), (1/27240646557740*I*sqrt(11258216119645872304850111) - 2336015071/65957981980, [], 1)]
#This is an unstable node or a saddle point since one of the eigenvalues is positive

#11
var("x,y")
x_prime = 0.0300033894396597*x-0.192897497982244*y
y_prime = 0.100838570753313*x-0.100836820083682*y
plot_vector_field([x_prime,y_prime], (x,-10,10),(y,-10,10))

(x, y)
#The graph of this system is uniform throughout and only shows off one equilibrium point. The graph of the nonlinear system wasn't uniform and contained more equilibrium points. This graph also shows a stable spiral equilibrium.

@interact
def interactive(v=(0.48,0.6,0.005)):
var("S,P")
t=srange(0,100,0.1)
S_prime = v-(0.23*S*P^2)
P_prime =(0.23*S*P^2)-0.4*P
sol=desolve_odeint([S_prime,P_prime],ics=[1.5,1.5],times=t,dvars=[S,P])
show(list_plot(sol,plotjoined=True))