Changeset 27430

Timestamp:

11/27/22 15:42:45 (2 years ago)

Author:

jdquinn

Message:

CHG: Missing translations to initialization.py; cleanup

Location:

issm/trunk-jpl

Files:

• etc/environment.sh (modified) (1 diff)
• src/m/classes/SMBpddSicopolis.m (modified) (1 diff)
• src/m/classes/initialization.m (modified) (2 diffs)
• src/m/classes/initialization.py (modified) (9 diffs)
• src/m/classes/model.m (modified) (3 diffs)
• src/m/classes/model.py (modified) (3 diffs)

issm/trunk-jpl/etc/environment.sh ¶

@@ -675 +675 @@
 fi
 
-VALGRIND_ROOT="${ISSM_EXT_DIR}/valgrind/install"
-path_prepend "${VALGRIND_ROOT}/bin"
+VALGRIND_ROOT="${ISSM_DIR}/externalpackages/valgrind/install"
+if [ -d "${VALGRIND_ROOT="${ISSM_DIR}/valgrind/install"
+}" ]; then
+        path_prepend "${VALGRIND_ROOT}/bin"
+fi
 
 DOXYGEN_ROOT="${ISSM_EXT_DIR}/doxygen/install"

issm/trunk-jpl/src/m/classes/SMBpddSicopolis.m ¶

@@ -72 +72 @@
                         self.desfac        = -log(2.0)/1000;
                         self.rlaps         = 7.4;
-         self.requested_outputs={'default'};
+                        self.requested_outputs={'default'};
 
                 end % }}}

issm/trunk-jpl/src/m/classes/initialization.m ¶

@@ -25 +25 @@
                 str                 = NaN;
                 sample              = NaN;
-                debris              = NaN;
+                debris              = NaN;
         end
         methods
@@ -157 +157 @@
                         yts=md.constants.yts;
 
-                        WriteData(fid,prefix,'object',self,'fieldname','vx','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
-                        WriteData(fid,prefix,'object',self,'fieldname','vy','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
+                        WriteData(fid,prefix,'object',self,'fieldname','vx','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',yts);
+                        WriteData(fid,prefix,'object',self,'fieldname','vy','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',yts);
                         WriteData(fid,prefix,'object',self,'fieldname','vz','format','DoubleMat','mattype',1,'scale',1./yts);
                         WriteData(fid,prefix,'object',self,'fieldname','pressure','format','DoubleMat','mattype',1);
-                        WriteData(fid,prefix,'object',self,'fieldname','sealevel','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
+                        WriteData(fid,prefix,'object',self,'fieldname','sealevel','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',yts);
                         WriteData(fid,prefix,'object',self,'fieldname','bottompressure','format','DoubleMat','mattype',1);
                         WriteData(fid,prefix,'object',self,'fieldname','str','format','DoubleMat','mattype',1);

issm/trunk-jpl/src/m/classes/initialization.py ¶

@@ -34 +34 @@
         self.str = np.nan
         self.sample = np.nan
+        self.debris = np.nan
 
         self.setdefaultparameters()
@@ -55 +56 @@
         s += '{}\n'.format(fielddisplay(self, 'channelarea', 'subglaciale water channel area (for GlaDS) [m2]'))
         s += '{}\n'.format(fielddisplay(self, 'sample', 'Realization of a Gaussian random field'))
+        s += '{}\n'.format(fielddisplay(self, 'debris', 'Surface debris layer [m]'))
         return s
     #}}}
@@ -65 +67 @@
         if 'StressbalanceAnalysis' in analyses and not solution == 'TransientSolution' and not md.transient.isstressbalance:
             if not np.any(np.logical_or(np.isnan(md.initialization.vx), np.isnan(md.initialization.vy))):
-                md = checkfield(md, 'fieldname', 'initialization.vx', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
-                md = checkfield(md, 'fieldname', 'initialization.vy', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
+                if np.size(md.initialization.vx) > 1 or np.size(md.initialization.vy) > 1:
+                    md = checkfield(md, 'fieldname', 'initialization.vx', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
+                    md = checkfield(md, 'fieldname', 'initialization.vy', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
         if 'MasstransportAnalysis' in analyses and not solution == 'TransientSolution' and not md.transient.ismasstransport:
             md = checkfield(md, 'fieldname', 'initialization.vx', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
@@ -79 +82 @@
             md = checkfield(md, 'fieldname', 'initialization.vy', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
             # Triangle with zero velocity
-            if np.any(np.logical_and(np.sum(np.abs(md.initialization.vx[md.mesh.elements - 1]), axis=1) == 0,
-                                     np.sum(np.abs(md.initialization.vy[md.mesh.elements - 1]), axis=1) == 0)):
-                md.checkmessage("at least one triangle has all its vertices with a zero velocity")
+            if np.any(np.logical_and(np.sum(np.abs(md.initialization.vx[md.mesh.elements - 1]), axis=1) == 0, np.sum(np.abs(md.initialization.vy[md.mesh.elements - 1]), axis=1) == 0, np.min(md.mask.ice_levelset(md.mesh.elements), [], 2) < 0)):
+                md.checkmessage('at least one triangle has all its vertices with a zero velocity')
         if 'ThermalAnalysis' in analyses and not solution == 'TransientSolution' and not md.transient.isthermal:
             md = checkfield(md, 'fieldname', 'initialization.vx', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
@@ -99 +101 @@
                 if (solution == 'TransientSolution' and md.transient.ishydrology) or solution == 'HydrologySolution':
                     md = checkfield(md, 'fieldname', 'initialization.watercolumn', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
-        if 'HydrologyDCAnalysis' in analyses:
-            if type(md.hydrology).__name__ == 'hydrologydc':
-                if (solution == 'TransientSolution' and md.transient.ishydrology) or solution == 'HydrologySolution':
-                    md = checkfield(md, 'fieldname', 'initialization.sediment_head', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
         if 'HydrologyTwsAnalysis' in analyses:
             if type(md.hydrology).__name__ == 'hydrologytws':
@@ -114 +112 @@
                 md = checkfield(md, 'fieldname', 'initialization.hydraulic_potential', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
                 md = checkfield(md, 'fieldname', 'initialization.channelarea', 'NaN', 1, 'Inf', 1, '>=', 0, 'size', [md.mesh.numberofelements])
+        if 'HydrologyDCInefficientAnalysis' in analyses:
+            if type(md.hydrology).__name__ == 'hydrologydc':
+                md = checkfield(md, 'fieldname', 'initialization.sediment_head', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
+        if 'HydrologyDCEfficientAnalysis' in analyses:
+            if type(md.hydrology).__name__ == 'hydrologydc':
+                if md.hydrology.isefficientlayer:
+                    md = checkfield(md, 'fieldname', 'initialization.epl_head', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
+                    md = checkfield(md, 'fieldname', 'initialization.epl_thickness', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
         if 'SamplingAnalysis' in analyses and not solution == 'TransientSolution' and not md.transient.issampling:
             if np.any(np.isnan(md.initialization.sample)):
                 md = checkfield(md, 'fieldname', 'initialization.sample', 'NaN', 1,'Inf', 1, 'size', [md.mesh.numberofvertices])
+        if 'DebrisAnalysis' in analyses:
+            if not np.isnan(md.initialization.debris):
+                if (solution == 'TransientSolution' and md.transient.ishydrology) or solution == 'HydrologySolution':
+                    md = checkfield(md, 'fieldname', 'initialization.debris', 'NaN', 1,'Inf', 1, 'size', [md.mesh.numberofvertices])
         return md
     #}}}
@@ -123 +133 @@
         yts = md.constants.yts
 
-        WriteData(fid, prefix, 'object', self, 'fieldname', 'vx', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1 / yts)
-        WriteData(fid, prefix, 'object', self, 'fieldname', 'vy', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1 / yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'vx', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1 / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'vy', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1 / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
         WriteData(fid, prefix, 'object', self, 'fieldname', 'vz', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1 / yts)
         WriteData(fid, prefix, 'object', self, 'fieldname', 'pressure', 'format', 'DoubleMat', 'mattype', 1)
@@ -140 +150 @@
         WriteData(fid, prefix, 'object', self, 'fieldname', 'hydraulic_potential', 'format', 'DoubleMat', 'mattype', 1)
         WriteData(fid, prefix, 'object', self, 'fieldname', 'sample', 'format', 'DoubleMat', 'mattype', 1)
+        WriteData(fid, prefix, 'object', self, 'fieldname', 'debris', 'format', 'DoubleMat', 'mattype', 1)
 
         if md.thermal.isenthalpy:
@@ -168 +179 @@
         self.dsl = project3d(md, 'vector', self.dsl, 'type', 'node', 'layer', 1)
         self.str = project3d(md, 'vector', self.str, 'type', 'node', 'layer', 1)
+        self.debris = project3d(md, 'vector', self.debris, 'type', 'node', 'layer', 1)
 
         # Lithostatic pressure by default

issm/trunk-jpl/src/m/classes/model.m ¶

@@ -543 +543 @@
                         end
 
-                        %Initialize the 2d mesh
+                        %Initialize 2d mesh
                         mesh=mesh2d();
                         mesh.x=md.mesh.x2d;
@@ -721 +721 @@
                         end
 
-                        %Initial 2d mesh 
+                        %Initial 2d mesh
                         if isa(md1.mesh,'mesh3dprisms'),
                                 flag_elem_2d=flag_elem(1:md1.mesh.numberofelements2d);
@@ -1091 +1091 @@
                         end
 
-                        %Initialize with the 2d mesh
+                        %Initialize with 2d mesh
                         mesh2d = md.mesh;
                         md.mesh=mesh3dprisms();

issm/trunk-jpl/src/m/classes/model.py ¶

@@ -633 +633 @@
             raise TypeError('Cannot extrude a 3d mesh (extrude cannot be called more than once)')
 
-        #Initialize with the 2d mesh
+        #Initialize with 2d mesh
         mesh2d = md.mesh
         md.mesh = mesh3dprisms()
@@ -853 +853 @@
         if not np.isnan(md.initialization.watercolumn).all():
             md.initialization.watercolumn = project2d(md, md.initialization.watercolumn, 1)
+        if not np.isnan(md.initialization.debris).all():
+            md.initialization.debris = project2d(md, md.initialization.debris, 1)
 
         # elementstype
@@ -960 +962 @@
                                     setattr(fieldr, solutionsubfield, project2d(md, subfield, 1))
 
-        # Initialize he 2d mesh
+        # Initialize 2d mesh
         mesh = mesh2d()
         mesh.x = md.mesh.x2d