%This test is a test from the EISMINT for Ice shelves 
% Vincent Rommelaere 1996

cluster='none';
np=2;

%The goal is to test the prognostic model
md=model;
md=mesh(md,'DomainOutline.exp',4550);
md=geography(md,'all','');
md=parameterize(md,'Square.par');
md=setelementstype(md,'macayeal','all');

%Evolution of the ice shelf
md.ndt=500*md.yts;
md.dt=5*md.yts;
md.artificial_diffusivity=1; %Better result with no artificial diffusivity

%launch transient solution
md.cluster=cluster;
md.np=np;
md=solve(md,'analysis_type','transient');

%plot results
plotmodel(md,'data',md.results.transient(end).thickness,'sectionvalue','CrossLine.exp')

%Don't forget to add these lines in icetransient2d.m, just before the computation of the thickness to change the thickness on the upper boundary condition
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%pos=find(m_p.ys~=0);
%m_p.ys(pos)=500+500/5*sin(2*pi*time/(200*md.yts))*ones(size(pos,1),1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%