radius=1.e6;
shelfextent=2.e5;

md=roundmesh(model(),radius,50000.);
%fix center node to 0,0
rad=sqrt(md.mesh.x.^2+md.mesh.y.^2);
pos=find(rad==min(rad));
md.mesh.x(pos)=0.; md.mesh.y(pos)=0.; %the closest node to the center is changed to be exactly at the center
xelem=md.mesh.x(md.mesh.elements)*[1;1;1]/3.;
yelem=md.mesh.y(md.mesh.elements)*[1;1;1]/3.;
rad=sqrt(xelem.^2+yelem.^2);
flags=zeros(md.mesh.numberofelements,1);
pos=find(rad>=(radius-shelfextent));
flags(pos)=1;
md=setmask(md,flags,''); 
md=parameterize(md,'../Par/RoundSheetShelf.par');
md=setflowequation(md,'SSA','all');
md=extrude(md,3,1);
md.cluster=generic('name',oshostname(),'np',3);

md.transient.isthermal=0;
md.transient.ismasstransport=0;
md.transient.isstressbalance=0;
md.transient.isgroundingline=1;

%test different grounding line dynamics.
md.groundingline.migration='AgressiveMigration';
md=solve(md,TransientSolutionEnum());
element_on_iceshelf_agressive=(md.results.TransientSolution.MaskGroundediceLevelset);

md.groundingline.migration='SoftMigration';
md=solve(md,TransientSolutionEnum());
element_on_iceshelf_soft=(md.results.TransientSolution.MaskGroundediceLevelset);

md.groundingline.migration='SubelementMigration';
md=solve(md,TransientSolutionEnum());
element_on_iceshelf_subelement=(md.results.TransientSolution.MaskGroundediceLevelset);

md.groundingline.migration='SubelementMigration2';
md=solve(md,TransientSolutionEnum());
element_on_iceshelf_subelement2=(md.results.TransientSolution.MaskGroundediceLevelset);

%Fields and tolerances to track changes
field_names     ={'ElementOnIceShelfAgressive','ElementOnIceShelfSoft','ElementOnIceShelfSubelement','ElementOnIceShelfSubelement'};
field_tolerances={1e-13,1e-13,1e-13,1e-13};
field_values={element_on_iceshelf_agressive,element_on_iceshelf_soft,element_on_iceshelf_subelement,element_on_iceshelf_subelement2};
