Changeset 24249 for issm/trunk-jpl/test/NightlyRun/test4003.m

Timestamp:

10/17/19 16:23:25 (5 years ago)

Author:

jdquinn

Message:

CHG: Eased tolerances.

File:

• issm/trunk-jpl/test/NightlyRun/test4003.m (modified) (35 diffs)

issm/trunk-jpl/test/NightlyRun/test4003.m

@@ -16 +16 @@
 presentdirectory=pwd;
 
-% {{{ Parameters: 
+% {{{ Parameters:
 if perform(org,'Parameters'),
         Nx=3; %number of longitude cells
@@ -28 +28 @@
         dLat=1; %latitude grid spacing
         delZ=10; %thickness of vertical levels
-        icefront_position_ratio=.75; 
+        icefront_position_ratio=.75;
         ice_thickness=1000;
         rho_ice=917;
@@ -42 +42 @@
         mlSalt   = 33.4; % open boundary salinity (PSU)
         mlTheta  = -1.9; % open boundary potential temperature (deg C)
-        obcUvel  = -0.1; % open boundary velocity (m/s)        
+        obcUvel  = -0.1; % open boundary velocity (m/s)
 
         MITgcmDeltaT=600; % MITgcm time step in seconds
@@ -69 +69 @@
 end
 % }}}
-% {{{ Bathymetry: 
+% {{{ Bathymetry:
 if perform(org,'Bathymetry'),
 
@@ -94 +94 @@
 end
 % }}}
-% {{{ IceSheetGeometry: 
+% {{{ IceSheetGeometry:
 if perform(org,'IceSheetGeometry'),
-        
+
         loaddata(org,'Parameters');
         loaddata(org,'Bathymetry');
@@ -102 +102 @@
         latmax=max(lat(:));
 
-        %put ice_thickness constant layer of ice over the bathymetry, unless it floats: 
+        %put ice_thickness constant layer of ice over the bathymetry, unless it floats:
         s=size(bathymetry);
         thickness=ice_thickness*ones(s);
 
-        %figure out ice shelf: 
+        %figure out ice shelf:
         pos=find(-di*thickness>bathymetry);
         iceshelf_mask=zeros(s);
@@ -116 +116 @@
         iceshelf_mask(pos)=0;
 
-        %compute draft of ice shelf: 
+        %compute draft of ice shelf:
         draft=bathymetry;
         pos=find(iceshelf_mask);
         draft(pos)=-di*thickness(pos);
-        pos=find(~ice_mask); 
+        pos=find(~ice_mask);
         draft(pos)=0;
 
@@ -128 +128 @@
 
 %Configure MITgcm
-% {{{ GetMITgcm: 
+% {{{ GetMITgcm:
 if perform(org,'GetMITgcm'),
   system([pwd '/../MITgcm/get_mitgcm_dngoldberg.sh']);
 end
 % }}}
-% {{{ BuildMITgcm: 
+% {{{ BuildMITgcm:
 if perform(org,'BuildMITgcm'),
     system(['../MITgcm/build_4003.sh generic ' pwd '/../MITgcm']);
@@ -139 +139 @@
 % }}}
 addpath(recursivepath([pwd '/../MITgcm']));
-% {{{ RunUncoupledMITgcm: 
+% {{{ RunUncoupledMITgcm:
 if perform(org,'RunUncoupledMITgcm'),
     % rename previous run directory and create new one
@@ -161 +161 @@
 
 %Configure ISSM
-% {{{ CreateMesh: 
+% {{{ CreateMesh:
 if perform(org,'CreateMesh'),
-        
-        %create model: 
+
+        %create model:
         md=model();
-        
-        %Grab lat,long from MITgcm: 
+
+        %Grab lat,long from MITgcm:
         long=readbin('run/XG.data',[3 200]);
         long=[long long(:,end)]; long=[long; -105.1250*ones(1,size(long,2))];
@@ -173 +173 @@
         lat=[lat -73.8832*ones(size(lat,1),1)]; lat=[lat; lat(end,:)];
 
-        %project lat,long: 
+        %project lat,long:
         [x,y]=ll2xy(lat(:),long(:),-1);
 
@@ -179 +179 @@
         index=[];
         %  C  D
-        %  A  B 
+        %  A  B
         for j=1:Ny-1,
-                for i=1:Nx-1, 
+                for i=1:Nx-1,
                         A=(j-1)*Nx+i;
                         B=(j-1)*Nx+i+1;
@@ -191 +191 @@
         end
 
-        %fill mesh and model: 
+        %fill mesh and model:
         md=meshconvert(md,index,x,y);
         md.mesh.lat=lat(:);
@@ -200 +200 @@
 end
 % }}}
-% {{{ MeshGeometry: 
+% {{{ MeshGeometry:
 if perform(org,'MeshGeometry'),
-        
+
         loaddata(org,'CreateMesh');
 
-        %transfer to vertices: 
+        %transfer to vertices:
         bathymetry=readbin('run/bathy.box',[3 200],1,'real*8');
         bathymetry=[bathymetry bathymetry(:,end)]; bathymetry=[bathymetry(1,:); bathymetry];
@@ -214 +214 @@
         thickness=[thickness thickness(:,end)]; thickness=[thickness; thickness(end,:)];
 
-        %start filling some of the fields 
+        %start filling some of the fields
         md.geometry.bed=bathymetry(:);
         md.geometry.thickness=thickness(:);
-        md.geometry.base=-917/1028*md.geometry.thickness; 
+        md.geometry.base=-917/1028*md.geometry.thickness;
         md.geometry.surface=md.geometry.base+md.geometry.thickness;
 
-        %nothing passes icefront: 
+        %nothing passes icefront:
         pos=find((~ice_mask(:) & ice_mask(:)~=0) | thickness(:)==0);
         md.geometry.thickness(pos)=1;
@@ -226 +226 @@
         md.geometry.base(pos)=-di*md.geometry.thickness(pos);
 
-        %level sets: 
+        %level sets:
         md.mask.groundedice_levelset=iceshelf_mask(:);
         md.mask.ice_levelset=-ones(md.mesh.numberofvertices,1);
 
-        pos=find(~ice_mask(:) & thickness(:)==0); md.mask.ice_levelset(pos)=1; 
+        pos=find(~ice_mask(:) & thickness(:)==0); md.mask.ice_levelset(pos)=1;
 
         savemodel(org,md);
 end
 % }}}
-% {{{ ParameterizeIce: 
+% {{{ ParameterizeIce:
 if perform(org,'ParameterizeIce'),
-        
+
         loaddata(org,'MeshGeometry');
 
@@ -243 +243 @@
         md.miscellaneous.name='test4003';
 
-        %initial velocity: 
+        %initial velocity:
         md.initialization.vx=zeros(md.mesh.numberofvertices,1);
         md.initialization.vy=zeros(md.mesh.numberofvertices,1);
         md.initialization.vz=zeros(md.mesh.numberofvertices,1);
 
-        %friction: 
+        %friction:
         md.friction.coefficient=0*ones(md.mesh.numberofvertices,1);
         pos=find(md.mask.groundedice_levelset>0);
@@ -260 +260 @@
         md.smb.mass_balance = 0*ones(md.mesh.numberofvertices,1);
 
-        %Flow law 
+        %Flow law
         md.materials.rheology_B=paterson(md.initialization.temperature);
         md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
         md.damage.D=zeros(md.mesh.numberofvertices,1);
         md.damage.spcdamage=NaN*ones(md.mesh.numberofvertices,1);
-        
+
         %the spcs going
         md.stressbalance.spcvx=NaN*ones(md.mesh.numberofvertices,1);
@@ -272 +272 @@
         md.stressbalance.referential=NaN*ones(md.mesh.numberofvertices,6);
         md.stressbalance.loadingforce=0*ones(md.mesh.numberofvertices,3);
-        md.masstransport.spcthickness=NaN*ones(md.mesh.numberofvertices,1); 
-
-        %get some flux at the ice divide: 
+        md.masstransport.spcthickness=NaN*ones(md.mesh.numberofvertices,1);
+
+        %get some flux at the ice divide:
         pos=find(md.mesh.lat==min(md.mesh.lat));
         md.masstransport.spcthickness(pos)=md.geometry.thickness(pos);
@@ -280 +280 @@
         md.stressbalance.spcvy(pos)=0;
 
-        %deal with boundaries, excluding icefront: 
+        %deal with boundaries, excluding icefront:
         pos=find(md.mesh.long==min(md.mesh.long) | md.mesh.long==max(md.mesh.long));
         md.stressbalance.spcvy(pos)=0;
@@ -295 +295 @@
         md.basalforcings.geothermalflux=.064*ones(md.mesh.numberofvertices,1);
 
-        %flow equations: 
+        %flow equations:
         md=setflowequation(md,'SSA','all');
 
@@ -301 +301 @@
 end
 % }}}
-% {{{ RunUncoupledISSM: 
+% {{{ RunUncoupledISSM:
 if perform(org,'RunUncoupledISSM'),
-        
+
         loaddata(org,'Parameters');
         loaddata(org,'ParameterizeIce');
 
-        %timestepping: 
+        %timestepping:
         md.timestepping.final_time=100;
         md.timestepping.time_step=0.5;
@@ -327 +327 @@
 
 %Run MITgcm/ISSM
-% {{{ RunCoupledMITgcmISSM: 
+% {{{ RunCoupledMITgcmISSM:
 if perform(org,'RunCoupledMITgcmISSM'),
-        
+
         loaddata(org,'Parameters');
         loaddata(org,'ParameterizeIce');
@@ -352 +352 @@
 
         % }}}
-        % {{{ prepare MITgcm 
+        % {{{ prepare MITgcm
         % rename previous run directory and create new one
         if exist ('run.old')
@@ -379 +379 @@
 end
 % }}}
-% {{{ RunCoupledMITgcmISSMMatlabCoupling: 
+% {{{ RunCoupledMITgcmISSMMatlabCoupling:
 if perform(org,'RunCoupledMITgcmISSMMatlabCoupling'),
 
-        %load data: 
+        %load data:
         loaddata(org,'Parameters');
         loaddata(org,'ParameterizeIce');
@@ -390 +390 @@
          floor(time_step*y2s*async_step_MITgcm_multiplier/MITgcmDeltaT));
 
-        % {{{ prepare MITgcm 
+        % {{{ prepare MITgcm
         % rename previous run directory and create new one
         if exist ('run.old')
@@ -449 +449 @@
             tline = fgetl(fidi);
             if ~ischar(tline), break, end
-            %do the change here: 
+            %do the change here:
             if strcmpi(tline,' OB_Iwest = 40*1,'),
                 fprintf(fido,'%s%i%s\n',' OB_Iwest = ',Ny,'*1,');
@@ -467 +467 @@
         writebin('run/bathymetry.bin',bathymetry);
         % }}}
-            
+
         % {{{ ISSM settings:
 
-        setenv('DYLD_LIBRARY_PATH', '/usr/local/gfortran/lib') 
-        %timestepping: 
+        setenv('DYLD_LIBRARY_PATH', '/usr/local/gfortran/lib')
+        %timestepping:
         md.timestepping.start_time=start_time;
         md.timestepping.final_time=final_time;
@@ -499 +499 @@
             writebin('run/icetopo.bin',draft);
 
-            % {{{ generate MITgcm parameter file data 
+            % {{{ generate MITgcm parameter file data
             fidi=fopen('../MITgcm/input/data','r');
             fido=fopen('run/data','w');
@@ -507 +507 @@
                 tline = fgetl(fidi);
                 if ~ischar(tline), break, end
-                %do the change here: 
+                %do the change here:
                 if strcmpi(tline,' xgOrigin = 0.0,'),
                     fprintf(fido,'%s%i%s\n',' xgOrigin = ',xgOrigin,',');
@@ -578 +578 @@
                 tmp(find(tmp<0))=0;
                 [im jm]=find(tmp); % horizontal indices where there is melt
-                
+
                 % Extrapolate T/S to locations where ice shelf retreated
                 for i=1:length(im)
@@ -612 +612 @@
             % }}}
 
-            % {{{ system call to run MITgcm 
+            % {{{ system call to run MITgcm
             cd run
             eval(['!mpiexec -np ' int2str(nPx*nPy) ' ./mitgcmuv']);
@@ -638 +638 @@
                  md.basalforcings.floatingice_melting_rate=melting_rate;
 
-            % {{{ run ISSM and recover results 
+            % {{{ run ISSM and recover results
 
             md.timestepping.start_time=t;
@@ -655 +655 @@
                 md.mask.groundedice_levelset=md.results.TransientSolution(end).MaskGroundediceLevelset;
                 md.results.TransientSolution(end).FloatingiceMeltingRate=md.basalforcings.floatingice_melting_rate;
-                
+
                 %save these results in the model, otherwise, they'll be wiped out
                 results(end+1)=md.results;
@@ -681 +681 @@
         'Melting2','Vx3','Vy3','Thickness3','Base3','MaskGroundediceLevelset3','FloatingiceMeltingRate3',...
         'Melting3','Vx4','Vy4','Thickness4','Base4','MaskGroundediceLevelset4','FloatingiceMeltingRate4','Melting4'};
-field_tolerances={2e-13,1e-13,5e-13,5e-13,1e-13,1e-13,1e-13,1e-12,...
+field_tolerances={2e-13,1e-13,6e-13,6e-13,1e-13,1e-13,1e-13,1e-12,...
         1e-13, 1e-12, 1e-12, 1e-13, 1e-13, 1e-13, 1e-12,...
         1e-13, 1e-12, 1e-12, 1e-13, 1e-13, 1e-13, 1e-12, 1e-13 };