%   =============
%   Simple EGT
%   2013 August 5
%   =============

%   Play with the following parameters
T = 0.5;  R = 0.25;    P = -0.25;   S = -0.5;
p = [T R P S]';

%   Initial conditions
C0 = 1000;
D0 = 100;

%   Time-interval of integration
MaxTime = 21;

%   ==================================================
%   Numerically integrate differential equations
%   ==================================================
[T,Y] = ode45(@SimpleEGTDEs,[0:3:MaxTime],[C0 D0],[],p);

figure;
plot(T,Y(:,1),'-s','MarkerSize',16,'LineWidth',5,'Color',[222 125 0]/255,'MarkerFaceColor',[255 177 100]/255);
hold on;
plot(T,Y(:,2),'--o','MarkerSize',16,'LineWidth',5,'Color',[11 132 199]/255,'MarkerFaceColor',[173 235 255]/255);
hold off;
set(gca,'FontSize',20,'FontName','Times New Roman','LineWidth',3);
xlabel('Time {\itt} (days)','FontSize',20);
ylabel('Numbers of cells','FontSize',20);
hl = legend('{\itx}({\itt})','{\ity}({\itt})');
ymax = 1.1*max(max(Y));
ylim([0 ymax]);
set(hl,'LineWidth',3,'FontSize',20);
grid on;

figure;
subplot(2,1,1);
plot(T,Y(:,1),'-s','MarkerSize',16,'LineWidth',5,'Color',[222 125 0]/255,'MarkerFaceColor',[255 177 100]/255);
set(gca,'FontSize',20,'FontName','Times New Roman','LineWidth',3);
xlabel('Time {\itt} (days)','FontSize',20);
ylabel('{\itx}({\itt}) (cells)','FontSize',20);
ymax1 = 1.1*max(Y(:,1));
ylim([0 ymax1]);
grid on;

subplot(2,1,2);
plot(T,Y(:,2),'--o','MarkerSize',16,'LineWidth',5,'Color',[11 132 199]/255,'MarkerFaceColor',[173 235 255]/255);
set(gca,'FontSize',20,'FontName','Times New Roman','LineWidth',3);
xlabel('Time {\itt} (days)','FontSize',20);
ylabel('{\ity}({\itt}) (cells)','FontSize',20);
ymax2 = 1.1*max(Y(:,2));
ylim([0 ymax2]);
grid on;

r = sqrt(2)*1000/0.75;
theta = [0 atan(0.5)/2 atan(0.5) pi/4 atan(2) atan(2)+(atan(0.5)/2) pi/2];
Cnor = r*cos(theta);
Dnor = r*sin(theta);
C = [0.25*Cnor 0.5*Cnor 0.75*Cnor Cnor 1.25*Cnor 1.5*Cnor 1.75*Cnor 2*Cnor];
D = [0.25*Dnor 0.5*Dnor 0.75*Dnor Dnor 1.25*Dnor 1.5*Dnor 1.75*Dnor 2*Dnor];

pC = C./(C + D);
pD = D./(C + D);

T = p(1);

DC = (R*pC + S*pD).*C;
DD = (T*pC + P*pD).*D;

figure;
q = quiver(C,D,DC,DD,'MaxHeadSize',0.1,'AutoScaleFactor',1,'AutoScale','off','Color',[0.4 0.4 0.4]);
set(q,'LineWidth',3);
set(gca,'FontSize',20,'FontName','Times New Roman','LineWidth',3);
hold on;
plot(Y(:,1),Y(:,2),'-o','MarkerSize',16,'LineWidth',5,'Color',[88 42 4]/255,'MarkerFaceColor',[152 72 7]/255);
hold on;
grid on;
xlabel('{\itx} (cells)','FontSize',20);
ylabel('{\ity} (cells)','FontSize',20);
axis equal;

figure;
q = quiver(C,D,DC,DD,'MaxHeadSize',0.1,'AutoScaleFactor',0.4,'AutoScale','off','Color',[0.4 0.4 0.4]);
set(q,'LineWidth',3);
set(gca,'FontSize',20,'FontName','Times New Roman','LineWidth',3);
grid on;
xlabel('{\itx} (cells)','FontSize',20);
ylabel('{\ity} (cells)','FontSize',20);
axis equal;