Index: /issm/trunk-jpl/test/Par/GiaBenchmarksAB.par
===================================================================
--- /issm/trunk-jpl/test/Par/GiaBenchmarksAB.par	(revision 14847)
+++ /issm/trunk-jpl/test/Par/GiaBenchmarksAB.par	(revision 14847)
@@ -0,0 +1,72 @@
+%Geometry specific to Experiments A and B 
+rad=600000;
+nv=md.mesh.numberofvertices; 
+for i=1:nv,
+	dist=sqrt(md.mesh.x(i).^2+md.mesh.y(i)^2);
+	if (dist <= rad)
+		md.geometry.thickness(i)=2000.0;
+	else
+		md.geometry.thickness(i)=1.0; % non-zero thickness
+	end
+end
+md.geometry.thickness=md.geometry.thickness';
+md.geometry.bed=zeros(md.mesh.numberofvertices,1); 
+md.geometry.surface=md.geometry.thickness+md.geometry.bed; 
+
+%Ice density used for benchmarking, not 917 kg/m^3
+md.materials.rho_ice=1000; %kg m^3
+
+%GIA parameters specific to Experiments A  and B
+md.gia.lithosphere_thickness=100*ones(md.mesh.numberofvertices,1); %in km
+md.materials.lithosphere_shear_modulus=6.7*10^10;                  %in Pa
+md.materials.lithosphere_density=3.36;                             %in g/cm^3
+md.materials.mantle_shear_modulus=1.45*10^11;                      %in Pa
+md.materials.mantle_viscosity=10^21;                               %in Pa.s
+md.materials.mantle_density=3.38;                                  %in g/cm^3
+
+
+%Initial velocity 
+x     = transpose(ncread('../Data/SquareSheetConstrained.nc','x'));
+y     = transpose(ncread('../Data/SquareSheetConstrained.nc','y'));
+vx    = transpose(ncread('../Data/SquareSheetConstrained.nc','vx'));
+vy    = transpose(ncread('../Data/SquareSheetConstrained.nc','vy'));
+index = transpose(ncread('../Data/SquareSheetConstrained.nc','index'));
+
+md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y);
+md.initialization.vy=InterpFromMeshToMesh2d(index,x,y,vy,md.mesh.x,md.mesh.y);
+clear vx vy x y index;
+md.initialization.vz=zeros(md.mesh.numberofvertices,1);
+md.initialization.pressure=zeros(md.mesh.numberofvertices,1);
+
+%Materials
+md.initialization.temperature=(273.-20.)*ones(md.mesh.numberofvertices,1);
+md.materials.rheology_B=paterson(md.initialization.temperature);
+md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
+
+%Friction
+pos=find(md.mask.elementonfloatingice);
+md.friction.coefficient=20.*ones(md.mesh.numberofvertices,1);
+md.friction.coefficient(md.mesh.elements(pos,:))=0;
+md.friction.p=ones(md.mesh.numberofelements,1);
+md.friction.q=ones(md.mesh.numberofelements,1);
+
+%Numerical parameters
+md.diagnostic.viscosity_overshoot=0.0;
+md.prognostic.stabilization=1.;
+md.thermal.stabilization=1.;
+md.verbose=verbose(0);
+md.settings.waitonlock=30;
+md.diagnostic.restol=0.05;
+md.steadystate.reltol=0.05;
+md.diagnostic.reltol=0.05;
+md.diagnostic.abstol=NaN;
+md.timestepping.time_step=1.;
+md.timestepping.final_time=3.;
+
+%Boundary conditions:
+md=SetIceSheetBC(md);
+
+%Change name so that no test have the same name
+A=dbstack;
+if (length(A)>2), md.miscellaneous.name=A(3).file(1:end-2); end
+