How to Solve and Plot Lotka-Volterra Differential Equations in Matlab

Question

I was wondering if someone might be able to help me solve the Lotka-Volterra equations using MatLab. My code doesn't seem to be working. I do the following:

Step 1 -

I created a file entitled pred_prey_odes.m containing the following code:

% the purpose of this program is to model a predator prey relationship 
% I will be using the Lotka-Volterra equations 

% Program consists of the following differential equations: 
% dY1/dt = a * Y1 - c * Y1 * Y2
% dY2/dt = b * Y2 - d * Y1 * Y2 

function dy = pred_prey_odes(t, y) 
% function that is to be integrated 

%select constants 
a = 1;  
b = 2;  
c = 3; 
d = 4; 

%set up differential equations 
dy = zeros(2,1); 
dy(1) = a * y(1) - c * y(1) * y(2); 
dy(2) = b * y(2) - d * y(1) * y(2); 

I saved the file and made sure it was in the current directory before typing the following code into the command window:

clc
tspan = [0, 20]; 
y0 = [10; 10]; 
ode = @(t, y) pred_prey_odes(t, y); 
[t, y] = ode45(ode, tspan, y0); 
plot (t,y)

However, no plot pops up. In fact, nothing happens in matlab and I can't even clear the command window. If I type clc nothing happens...

Any help would be appreciated!

Thanks!

-Sneha Inguva

Answer 1

Your code is fine. With this choice of the parameters, the predators simply die out and the prey experience exponential growth. Since the value of the predator variable is tiny, ode45 thinks it needs to take really tiny steps to resolve it correctly. For example I ran your code with tspan = [0,10] and the number of time steps was 227485.

You might have better luck either changing the parameters, or if these are really the parameters you care about, you could use ode23 which should be able to take larger time steps due to it's stability properties. However that might be slow too. Fundamentally your problem is that MATLAB is trying to control the relative error based on the predator variable which is tiny and basically meaningless. If you rolled your own RK4 method that did not have adaptive time stepping it would probably work fine.

EDIT There is a sign error in the second equation. Fixing that should give oscillations instead of the predators going extinct.