Index: /issm/trunk-jpl/examples/SlrGRACE/runme.m
===================================================================
--- /issm/trunk-jpl/examples/SlrGRACE/runme.m	(revision 25936)
+++ /issm/trunk-jpl/examples/SlrGRACE/runme.m	(revision 25937)
@@ -2,8 +2,5 @@
 addpath('../Data','../Functions');
 
-steps=[1:7];
-
-try
- % [1:7]
+steps=[6];
 
 if any(steps==1) %{{{
@@ -11,5 +8,5 @@
 
 	numrefine=1;
-	resolution=150*1e3;			% inital resolution [m]
+	resolution=300*1e3;			% inital resolution [m]
 	radius = 6.371012*10^6;		% mean radius of Earth, m
 	mindistance_coast=150*1e3;	% coastal resolution [m]
@@ -22,14 +19,14 @@
 	for i=1:numrefine,
 
-		md.mask.ocean_levelset=gmtmask(md.mesh.lat,md.mesh.long);
+		ocean_mask=gmtmask(md.mesh.lat,md.mesh.long);
 
 		distance=zeros(md.mesh.numberofvertices,1);
 
-		pos=find(~md.mask.ocean_levelset);	coaste.lat=md.mesh.lat(pos);	coaste.long=md.mesh.long(pos); 
-		pos=find(md.mask.ocean_levelset);	coasto.lat=md.mesh.lat(pos);	coasto.long=md.mesh.long(pos); 
+		pos=find(~ocean_mask);	coaste.lat=md.mesh.lat(pos);	coaste.long=md.mesh.long(pos); 
+		pos=find(ocean_mask);	coasto.lat=md.mesh.lat(pos);	coasto.long=md.mesh.long(pos); 
 
 		for j=1:md.mesh.numberofvertices
 			phi1=md.mesh.lat(j)/180*pi; lambda1=md.mesh.long(j)/180*pi;
-			if md.mask.ocean_levelset(j),
+			if ocean_mask(j),
 				phi2=coaste.lat/180*pi; lambda2=coaste.long/180*pi;
 				deltaphi=abs(phi2-phi1); deltalambda=abs(lambda2-lambda1);
@@ -44,5 +41,5 @@
 		pos=find(distance<mindistance_coast); distance(pos)=mindistance_coast;
 		
-		pos2=find(md.mask.ocean_levelset~=1 & distance>mindistance_land);
+		pos2=find(ocean_mask~=1 & distance>mindistance_land);
 		distance(pos2)=mindistance_land;
 
@@ -51,5 +48,8 @@
 	end
 
-	md.mask.ocean_levelset=gmtmask(md.mesh.lat,md.mesh.long);
+	ocean_mask=gmtmask(md.mesh.lat,md.mesh.long);
+	pos = find(ocean_mask==0); 
+	md.mask.ocean_levelset=-ones(md.mesh.numberofvertices,1); 
+	md.mask.ocean_levelset(pos)=1; 
 
 	save ./Models/SlrGRACE_Mesh md;
@@ -67,9 +67,9 @@
 	water_load = grace(md.mesh.elements,md.mesh.lat,md.mesh.long,time_range(1),time_range(2));
 
-	md.slr.deltathickness = water_load*md.materials.rho_freshwater/md.materials.rho_ice;
+	md.solidearth.surfaceload.icethicknesschange = water_load*md.materials.rho_freshwater/md.materials.rho_ice;
 
 	save ./Models/SlrGRACE_Loads md;
 
-	plotmodel (md,'data',md.slr.deltathickness,'view',[90 -90],'caxis',[-.1 .1],'title','Ice height equivalent [m]');
+	plotmodel (md,'data',md.solidearth.surfaceload.icethicknesschange,'view',[90 -90],'caxis',[-.1 .1],'title','Ice height equivalent [m]');
 end %}}}
 
@@ -80,19 +80,18 @@
 	md.solidearth.lovenumbers = lovenumbers('maxdeg',10000);
 
-	md.mask.ocean_levelset = ones(md.mesh.numberofvertices,1);
-	md.mask.ice_levelset = ones(md.mesh.numberofvertices,1);
-	pos=find(md.slr.deltathickness~=0);
-	md.mask.ice_levelset(md.mesh.elements(pos,:))=-1;
-	md.mask.land_levelset = 1-md.mask.ocean_levelset;
-
-	md.slr.sealevel=zeros(md.mesh.numberofvertices,1);
-	md.slr.steric_rate=zeros(md.mesh.numberofvertices,1);
-	md.slr.ocean_area_scaling = 1;
-
-	di=md.materials.rho_ice/md.materials.rho_water;
-	md.geometry.thickness=ones(md.mesh.numberofvertices,1);
-	md.geometry.surface=(1-di)*zeros(md.mesh.numberofvertices,1);
-	md.geometry.base=md.geometry.surface-md.geometry.thickness;
-	md.geometry.bed=md.geometry.base;
+	md.mask.ice_levelset=-md.mask.ocean_levelset; 
+
+	md.solidearth.initialsealevel=zeros(md.mesh.numberofvertices,1);
+	md.dsl.global_average_thermosteric_sea_level_change=[0;0];
+	md.dsl.sea_surface_height_change_above_geoid=zeros(md.mesh.numberofvertices+1,1);
+	md.dsl.sea_water_pressure_change_at_sea_floor=zeros(md.mesh.numberofvertices+1,1);
+
+	md.solidearth.settings.ocean_area_scaling=1;
+
+	% arbitary to pass consistency check. 
+	md.geometry.bed=-ones(md.mesh.numberofvertices,1);
+	md.geometry.surface=ones(md.mesh.numberofvertices,1);
+	md.geometry.base=md.geometry.bed; 
+	md.geometry.thickness=md.geometry.surface-md.geometry.base; 
 
 	md.initialization.temperature=273.25*ones(md.mesh.numberofvertices,1);
@@ -109,7 +108,9 @@
 	md = loadmodel('./Models/SlrGRACE_Parameterization');
 
-	md.slr.rigid=1;
-	md.slr.elastic=1;
-	md.slr.rotation=1;
+	md.solidearth.settings.rigid=1;
+	md.solidearth.settings.elastic=1;
+	md.solidearth.settings.rotation=1;
+
+	md.solidearth.requested_outputs = {'Sealevel','SealevelRSL'};
 
 	md.cluster=generic('name',oshostname(),'np',3);
@@ -125,6 +126,6 @@
 	md = loadmodel('./Models/SlrGRACE_Solution');
 
-	sol1 = md.slr.deltathickness*100;						% [cm]
-	sol2 = md.results.SealevelriseSolution.Sealevel*1000;	% [mm]
+	sol1 = md.solidearth.surfaceload.icethicknesschange*100;						% [cm]
+	sol2 = md.results.SealevelriseSolution.SealevelRSL*1000;	% [mm]
 
 	sol_name={'Change in water equivalent height [cm]', 'Relative sea-level [mm]'};
@@ -190,7 +191,7 @@
 
 	[ne,nt]=size(water_load);
-	md.slr.deltathickness = zeros(ne+1,nt);
-	md.slr.deltathickness(1:ne,:) = water_load*md.materials.rho_freshwater/md.materials.rho_ice;
-	md.slr.deltathickness(ne+1,:)=[1:nt]; % times
+	md.solidearth.surfaceload.icethicknesschange = zeros(ne+1,nt);
+	md.solidearth.surfaceload.icethicknesschange(1:ne,:) = water_load*md.materials.rho_freshwater/md.materials.rho_ice;
+	md.solidearth.surfaceload.icethicknesschange(ne+1,:)=[1:nt]; % times
 
 	md.transient.issmb=0;