function A=Game_of_Life4(A0,tmax,m,n)
%Conway's Game of Life
%A0: seed array containing initial 'alive' cells
%tmax: number of generations
%m,n are size of game board
%GoL:
%1. Survival: an alive cell will stay alive if it has 2 or 3 live neighbors
%2. Birth: an unalive cell will 'give birth' if it has 3 live neighbors
%3. Death: an alive cell will die of loneliness if it has 0 or 1 alive
%neighbors, or die of overcrowding if it has 4 or more alive neighbors
%----------------------------------------
%default values for arguments
if ~exist('m','var') %if game board size is not given, default to 200
    m=200;
end
if ~exist('n','var') %if game board size is not given, default to 200
    n=200;
end
if ~exist('tmax','var') %if tmax is not given, do 2000 generations
    tmax=2000;
end
if ~exist('A0','var') %if seed is not given, generate random 100x100 arry of 0s and 1s
    A0=round(rand(100,100));
end

figure
colormap([0 1 1; 1 .5 .5]); %saturate colors so only black and white appear
A=zeros(m,n); %array to store values
[mm,nn]=size(A0); %define the size of the seed
i=floor((m-mm)/2); j=floor((n-nn)/2); %indices to put the seed in the approximate middle of the board
A(i:i+mm-1,j:j+nn-1)=A0; %put the seed into the game board
image(A+1); %display matrix in the figure
pause %pauses computation to show initial state

for t=1:tmax
    B=zeros(m,n); %auxilary array to store values
    %for this case, keep the boundary all zeros
    for i=2:m-1
        for j=2:n-2
            s=0; %initialize sum of neighbors to zero, then find how many alive neighbors
            s=s+A(i-1,j-1)+A(i-1,j)+A(i-1,j+1); %alive neighbors in row above
            s=s+A(i,j-1)+A(i,j+1); %alive neighbors in the same row
            s=s+A(i+1,j-1)+A(i+1,j)+A(i+1,j+1); %alive neighbors in row below
            if s==2 %if two alive neighbors, no change
                B(i,j)=A(i,j);
            elseif s==3 %if three neighbors, stay alive or give birth
                B(i,j)=1;
            else %0-1 alive neighbors dies of loneliness, 4+ alive dies of overcrowding
                B(i,j)=0;
            end
        end
    end
    A=B; %values in the auxilary matrix are copied into A
    image(A+1);
    drawnow %creates image
end 
end