%Hannah Bourne
%Matlab Final Project
%Greenland Ice Sheet Melt
%December 3, 2013

clf

%I have commented out this part because I only wanted to create the
%"Greenland.txt" file once
% ocean= load('topography_180x360_grid.txt');
% size(ocean)
% Greenland = ocean(5:33,300:350);
% image(Greenland)
% size(Greenland)
% 
% Greenland = num2str(Greenland,'%8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f %8.1f\n')
% 
% file=sprintf('Greenland.txt');
% fid=fopen(file,'w');
% 
% 
% for i=1:28
%     fprintf(fid,'%s\n',Greenland(i,:));
% end
% 
Greenland= load('Greenland.txt');
size(Greenland);
Greenland(19:23,38:48)=-100; %remove Iceland
%Greenland is 28 by 51 pixels

%Set Initial Conditions
Greenland_Ice=zeros(28,51);
f=0;
for i=1:28;
    for j=1:50;
        if Greenland(i,j)>0;
            Greenland_Ice(i,j)=40;
        elseif Greenland(i,j)<0;
            Greenland_Ice(i,j)=-5;
        end
    end
end
Greenland_Ice(:,1)=-5;
Greenland_Ice(:,2)=-5;

Greenland_Temp=zeros(28,51);
for i=1:28;
    for j=1:51;
        if Greenland(i,j)>0;
            Greenland_Temp(i,j)=-2;
        elseif Greenland(i,j)<0;
            Greenland_Temp(i,j)=5;
        end
    end
end

for v=1:50000;
while max(max(Greenland_Ice))>=0;
    Greenland_Temp_new=Greenland_Temp;
    Greenland_Ice_new=Greenland_Ice;
    for i = 2:27;
        for j = 2:50;
            Greenland_Temp_new(i,j)=(Greenland_Temp(i+1,j)+Greenland_Temp(i-1,j)+Greenland_Temp(i,j+1)+Greenland_Temp(i,j-1))/4;
            
            if Greenland(i,j)>0 && Greenland_Temp(i,j)>=0;
                Greenland_Ice_new(i,j)=Greenland_Ice(i,j)-.1;
                if Greenland_Ice_new(i,j)<0;
                    Greenland_Ice_new(i,j)=0;
                end
            end
            Greenland_Temp=Greenland_Temp_new;
            Greenland_Ice=Greenland_Ice_new;
        end
        
    end
    

    surf(Greenland_Ice);
    axis([0 51 0 28 -5 40])
    view([-155 76]);
    pause(.001);
    if mod(v+9,10) == 00
        f=f+1;
        surf(Greenland_Ice)
        axis([0 51 0 28 -5 40])
        view([-155 76]);
        saveas(gcf, ['Ice_Melt',sprintf('%03d',f),'.png'])
    end
end
end

