%LogisticReaction2D

n=100; % Number of spatial points
m=10000; % Number of temporal point % Final time point
D1=.01; % Diffusion constant
D2=.01;
b=3;
a=1;

%Construcstion of spatial and temporal domains and spacing.
x=linspace(0,1,n);
y=linspace(0,1,n);
dx=x(2)-x(1);
dt=dx^2/(4*D1); %This satisfies the CFL condition
t=linspace(0,(m-1)*dt,m);

[X,Y]=meshgrid(x,y);

%Construction of second derivative operator
DXX=-2*diag(ones(1,n),0)+diag(ones(1,n-1),1)+diag(ones(1,n-1),-1);
DXX(1,:)=0;
DXX(n,:)=0;
DXX=DXX./dx^2;

DXX=sparse(DXX);

u1=rand;
u2=rand;

v1=rand;
v2=rand;

% Construction of initial distribution
u0=sin(pi*X).*sin(pi*Y).*exp(-4*((X-u1).^2+(Y-u2).^2));
v0=3*sin(pi*X).*sin(pi*Y).*exp(-10*((X-v1).^2+(Y-v2).^2));
u0=u0';
v0=v0';


%% Method of lines
%Initialization of distribution
U=u0;
V=v0;
clinesu=linspace(0,3,10);
clinesv=linspace(0,7,10);

vidObj = VideoWriter('Example3.avi');
    open(vidObj);
    
%Plotting Initial guess
figure('units', 'normalized', 'outerposition', [0,0,1,1]);
hold on;
subplot(1,2,1)
contourf(X,Y,U,clinesu,'edgecolor','none');
ylabel('$y$','Interpreter','latex','Fontsize',18);
xlabel('$x$','Interpreter','latex','Fontsize',18);
title(strcat('Sheep, ','$\,t=$',num2str(t(1))),'Interpreter','latex','Fontsize',18 );
axis([0,1,0,1]);
caxis([0,3]);
colorbar;

subplot(1,2,2)
contourf(X,Y,V,clinesv,'edgecolor','none');
ylabel('$y$','Interpreter','latex','Fontsize',18);
xlabel('$x$','Interpreter','latex','Fontsize',18);
title(strcat('Wolves, ','$\,t=$',num2str(t(1))),'Interpreter','latex','Fontsize',18 );
axis([0,1,0,1]);
caxis([0,7]);
colorbar;
drawnow;

currFrame = getframe(gcf);
       writeVideo(vidObj,currFrame);
       
%Looping through Changes in time
for i=1:m-1
    Unew=U+(D1*(DXX*U+(DXX*U')')+U.*(b-U-V./(1+U)))*dt;
    Vnew=V+(D2*(DXX*V+(DXX*V')')+V.*(U./(1+U)-a*V))*dt;
    
    Unew(1,:)=0; %boundary conditions
    Unew(end,:)=0;%boundary conditions
    Unew(:,1)=0;
    Unew(:,end)=0;
    
    Vnew(1,:)=0; %boundary conditions
    Vnew(end,:)=0;%boundary conditions
    Vnew(:,1)=0;
    Vnew(:,end)=0;
    
    U=Unew;
    V=Vnew;

    
    if mod(i,10)==0, %Plot every 10th iterate
            clf;
            subplot(1,2,1);
        contourf(X,Y,U,clinesu,'edgecolor','none');
ylabel('$y$','Interpreter','latex','Fontsize',18);
xlabel('$x$','Interpreter','latex','Fontsize',18);
title(strcat('Sheep, ','$\,t=$',num2str(t(i))),'Interpreter','latex','Fontsize',18 );
axis([0,1,0,1]);
colorbar;
caxis([0,3]);

subplot(1,2,2);
contourf(X,Y,V,clinesv,'edgecolor','none');
ylabel('$y$','Interpreter','latex','Fontsize',18);
xlabel('$x$','Interpreter','latex','Fontsize',18);
title(strcat('Wolves, ','$\,t=$',num2str(t(i))),'Interpreter','latex','Fontsize',18 );
axis([0,1,0,1]);
caxis([0,7]);
colorbar;
drawnow;

currFrame = getframe(gcf);
       writeVideo(vidObj,currFrame);
    end
    
end

close(vidObj);
