function f = cellular_dynamics(Tf,Off)

%Tf - total time of simulation
%Off - time after which the input is turned off

Vm = 1.0;
Vc = 10.0;
k12 = 1;
k2 = 1;

A =[-k12/Vm k12/Vm; k12/Vc -(k12+k2)/Vc];

disp('system matrix:')
disp(num2str(A))

disp(['eigenvalues: ' num2str(eig(A)')])

close all

rc = 1;
rm = 0.5;

T = linspace(0,Tf,2*Tf);
y0 = [0 0];
[t, y] = ode45(@(t,y) cell_dyn(t,y,Vm,Vc,k12,k2,Off),T,y0');


figure(2)
set(2,'Color','b')
x1fine = linspace(-rc,rc,100);
x2fine = linspace(-rm,rm,100);
y1fine = sqrt(rc^2-x1fine.^2);
y2fine = sqrt(rm^2-x2fine.^2);

plot(x1fine,y1fine,'r','LineWidth',3);
hold on
plot(x1fine,-y1fine,'r','LineWidth',3);

plot(x2fine,y2fine,'g','LineWidth',3);
hold on
plot(x2fine,-y2fine,'g','LineWidth',3);
set(gca,'DataAspectRatio',[1 1 1])

p1 = plot(0,0,'yo','MarkerFaceColor','y','MarkerSize',4);
p2 = plot(0,0,'yo','MarkerFaceColor','y','MarkerSize',4);
set([p1,p2],'Visible','off');
th = title('Production On','FontWeight','b','FontAngle','i');
set(th,'Color','y','FontSize',16);
drawnow
axis off

disp('hit any key to continue')
pause

Npt = 1000;
for i = 1:length(T)
    
    nc = round(y(i,2)*Npt);
    nm = round(y(i,1)*Npt);
    
    pc = rand(1,nc) * [rc - rm] + rm;
    phic = rand(1,nc) * 2*pi;
    pm = rand(1,nm) * rm;
    phim = rand(1,nm) * 2*pi;
    
    set(p1,'XData',pc.*cos(phic),'YData',pc.*sin(phic),'Visible','on');
    set(p2,'XData',pm.*cos(phim),'YData',pm.*sin(phim),'Visible','on');
    pause(0.05)
    if T(i) > Off
        set(th,'Color','r','String','Production Off');
    end
end



function dy = cell_dyn(t,y,Vm,Vc,k12,k2,Off)



dy = zeros(2,1);

u = 1;

if t < Off
    u = 1;
    %u = sin(2*pi/10*t);
else
    u = 0;
end

dy(1) = 1/Vm * [u - k12 * (y(1) - y(2))];
dy(2) = 1/Vc * [k12 * (y(1) - y(2)) - k2 * y(2)];