tspan = [0 20];                      % time span  
m0 = [.3  .3  .3 1  1  1    1.7 1.7 1.7 .7;...
      .3 1.3 1.9 .3 .3 1.7   1   1.7  .3 1.7]; % initial condition
N = size(m0);
for i=1:N(1,2)
[t,mi] = ode23('ricci3d',tspan,m0(:,i));   % run ode23
x = mi(:,1); y = mi(:,2)               % x = col 1 of z, y = col 2
plot(x,y)
hold on
axis([0 2 -2 2])
end
xlabel('m2 = a2/a1'), ylabel('m3 = a3/a1')
title ('Phase Plane Plot')