%Start defining model parameters here

%Geometry and observation
load('./../Data/Pig.data','-mat');
md.vx_obs   =InterpFromMeshToMesh2d(index,x,y,vx_obs,md.x,md.y,0);
md.vy_obs   =InterpFromMeshToMesh2d(index,x,y,vy_obs,md.x,md.y,0);
md.surface  =InterpFromMeshToMesh2d(index,x,y,surface,md.x,md.y,0);
md.thickness=InterpFromMeshToMesh2d(index,x,y,thickness,md.x,md.y,0);
md.bed=md.surface-md.thickness;
clear surface thickness vx_obs vy_obs x y index;
md.vx=md.vx_obs;
md.vy=md.vy_obs;
md.vz=zeros(md.numberofgrids,1);
md.pressure=zeros(md.numberofgrids,1);

%Materials
md.observed_temperature=(273-20)*ones(md.numberofgrids,1);
md.rheology_B=paterson(md.observed_temperature);
md.rheology_n=3*ones(md.numberofelements,1);
md.temperature=md.observed_temperature;

%Friction
pos=find(md.elementoniceshelf);
md.drag_type=2;
md.drag_coefficient=50*ones(md.numberofgrids,1);
md.drag_coefficient(md.elements(pos,:))=0;
md.drag_p=ones(md.numberofelements,1);
md.drag_q=ones(md.numberofelements,1);

%Numerical parameters
md.viscosity_overshoot=0.3;
md.artificial_diffusivity=1;
md.verbose=0;
md.waitonlock=30;
md.np=3;
md.dt=1;
md.ndt=2;
md.eps_res=0.05;
md.eps_rel=1;
md.eps_abs=NaN;

%Boundary conditions:
md=SetMarineIceSheetBC(md);

%Change name so that no test have the same name
A=dbstack;
if (length(A)==3), md.name=A(3).file(1:end-2); end