Index: /issm/trunk-jpl/etc/environment.sh =================================================================== --- /issm/trunk-jpl/etc/environment.sh (revision 27429) +++ /issm/trunk-jpl/etc/environment.sh (revision 27430) @@ -675,6 +675,9 @@ 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" Index: /issm/trunk-jpl/src/m/classes/SMBpddSicopolis.m =================================================================== --- /issm/trunk-jpl/src/m/classes/SMBpddSicopolis.m (revision 27429) +++ /issm/trunk-jpl/src/m/classes/SMBpddSicopolis.m (revision 27430) @@ -72,5 +72,5 @@ self.desfac = -log(2.0)/1000; self.rlaps = 7.4; - self.requested_outputs={'default'}; + self.requested_outputs={'default'}; end % }}} Index: /issm/trunk-jpl/src/m/classes/initialization.m =================================================================== --- /issm/trunk-jpl/src/m/classes/initialization.m (revision 27429) +++ /issm/trunk-jpl/src/m/classes/initialization.m (revision 27430) @@ -25,5 +25,5 @@ str = NaN; sample = NaN; - debris = NaN; + debris = NaN; end methods @@ -157,9 +157,9 @@ 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); Index: /issm/trunk-jpl/src/m/classes/initialization.py =================================================================== --- /issm/trunk-jpl/src/m/classes/initialization.py (revision 27429) +++ /issm/trunk-jpl/src/m/classes/initialization.py (revision 27430) @@ -34,4 +34,5 @@ self.str = np.nan self.sample = np.nan + self.debris = np.nan self.setdefaultparameters() @@ -55,4 +56,5 @@ 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,6 +67,7 @@ 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,7 +82,6 @@ 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,8 +101,4 @@ 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,7 +112,19 @@ 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,6 +133,6 @@ 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,4 +150,5 @@ 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,4 +179,5 @@ 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 Index: /issm/trunk-jpl/src/m/classes/model.m =================================================================== --- /issm/trunk-jpl/src/m/classes/model.m (revision 27429) +++ /issm/trunk-jpl/src/m/classes/model.m (revision 27430) @@ -543,5 +543,5 @@ end - %Initialize the 2d mesh + %Initialize 2d mesh mesh=mesh2d(); mesh.x=md.mesh.x2d; @@ -721,5 +721,5 @@ end - %Initial 2d mesh + %Initial 2d mesh if isa(md1.mesh,'mesh3dprisms'), flag_elem_2d=flag_elem(1:md1.mesh.numberofelements2d); @@ -1091,5 +1091,5 @@ end - %Initialize with the 2d mesh + %Initialize with 2d mesh mesh2d = md.mesh; md.mesh=mesh3dprisms(); Index: /issm/trunk-jpl/src/m/classes/model.py =================================================================== --- /issm/trunk-jpl/src/m/classes/model.py (revision 27429) +++ /issm/trunk-jpl/src/m/classes/model.py (revision 27430) @@ -633,5 +633,5 @@ 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,4 +853,6 @@ 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,5 +962,5 @@ setattr(fieldr, solutionsubfield, project2d(md, subfield, 1)) - # Initialize he 2d mesh + # Initialize 2d mesh mesh = mesh2d() mesh.x = md.mesh.x2d