# Plotting Limit Cycle

I have Lotka Volterra Equations where the prey equation is modified. The Equations are

x'[t]=a*x[t]*(1-x[t]/k)-b*x[t]*y[t]
y'[t]=-c*y[t]+d*x[t]*y[t]


To plot the limit cycles, i used NDSolve to solve the coupled nonlinear differential equations for x[t] and y[t]. I then used ParametricPlot to plot for x[t] and y[t] with respect to t, but I'm not getting appropriate limit cycles.

I did a steady state analysis of the equations and got certain conditions on the value of a, b, c, d and k. The equation has a non-trivial steady state, at $(x,y)=(c/d,(a/b)-(ac/bdk))$. The commands are:

Equations = {x'[t] == a*x[t]*(1 - x[t]/k) - b*x[t]*y[t], x == 2, y == 1, y'[t] == -c*y[t] + d*x[t]*y[t]};
s = NDSolve[Equations, {x, y}, {t, 0, time}];
ParametricPlot[Evaluate[{x[t], y[t]} /. s], {t, 0, 30}]


For the value of a = 2; b = 1; c = 2.5; d = 1.2; k = 2.8, I'm supposed to get a unstable spiral. But instead, the trajectory approaches the steady state of $(x,y) = (2.08333 , 0.511905)$.

Any help will be greatly appreciated.

• Please give the values of the constants you have used, and also the NDSolve and Plotting commands you used. Also, you say you don't get appropriate limit cycles, so what would you have expected instead of what you get? – Marius Ladegård Meyer Jun 11 '16 at 20:51
• @MariusLadegårdMeyer - I did a steady state analysis of the equations and got certain conditions on the value of a,b,c,d and k. Now, the equation has a non trivial steady state, at (x,y)=(c/d,(a/b)-(ac/bdk)).S s The commands are - Equations = {x'[t] == ax[t]*(1 - x[t]/k) - bx[t]*y[t], x == 2, y == 1, y'[t] == -cy[t] + dx[t]*y[t]}; s = NDSolve[{Equations}, {x, y}, {t, 0, time}]; ParametricPlot[Evaluate[{x[t], y[t]} /. s], {t, 0, 30}] – Sum-Meister Jun 11 '16 at 20:53
• Thanks. Just to have a reasonable starting point, what values did you use for a, b, c, d, k? NDSolve will not do anything without specifying these. – Marius Ladegård Meyer Jun 11 '16 at 20:57
• For the value of a = 2; b = 1; c = 2.5; d = 1.2; k = 2.8, I'm supposed to get a unstable spiral. But instead, the Parametric Plot goes to the steady state of (x,y)=(2.08333 , 0.511905) – Sum-Meister Jun 11 '16 at 20:58

I'm afraid that you have calculated the stability incorrectly. Here is the Jacobian of your system:

a = 2; b = 1; c = 2.5; d = 1.2; k = 2.8;
jac[x_, y_] := {{D[a x*(1 - x/k) - b x*y, x], D[a x*(1 - x/k) - b x*y, y]},
{D[-c y + d x*y, x], D[-c y + d x*y, y]}};


At the equilibrium, this is:

jacEq = jac[x, y] //. {x -> 2.08333, y -> 0.511905}


The eigenvalues of this are:

Re[Eigenvalues[jacEq]]
{-0.744047, -0.744047}


Hopf bifurcation analysis

xp = a*x*(1 - x/k) - b*x*y;
yp = (-c)*y + d*x*y;


The differential system

sys = {xp, yp}

{a*x*(1 - x/k) - b*x*y, (-c)*y + d*x*y}


The equilibrium points

eqp = Solve[sys == 0, {x, y}]

{{x-> c/d, y-> -((a*(c - d*k))/(b*d*k))}, {x -> 0, y -> 0}, {x -> k, y -> 0}}


The Jacobian in eqp[] (coexistence equilibrium)

J1 = FullSimplify[D[sys, {{x, y}}] /.eqp[]]
MatrixForm[J1]

{{-((a*c)/(d*k)), -((b*c)/d)}, {(a*(-c + d*k))/(b*k), 0}}


The Hopf bifurcation conditions for limit cycles require that the Jacobian matrix have a null trace for some parameter, but the trace of J1 is null if and only if a=0 or c=0, but all the parameters of the system are strictly positive. For thus, the Hopf bifurcation not take place.