Index: /issm/trunk-jpl/src/m/classes/materials.m =================================================================== --- /issm/trunk-jpl/src/m/classes/materials.m (revision 24755) +++ /issm/trunk-jpl/src/m/classes/materials.m (revision 24756) @@ -5,5 +5,5 @@ classdef materials < dynamicprops - properties (SetAccess=public) + properties (SetAccess=public) nature={}; %all properties are dynamic. @@ -13,18 +13,18 @@ if nargin==0 self.nature={'ice'}; - else + else self.nature=varargin; end - - %check this is acceptable: - for i=1:length(self.nature), - if ~(strcmpi(self.nature{i},'litho') | strcmpi(self.nature{i},'ice') | strcmpi(self.nature{i},'hydro')), - error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')'); - end - end - - %start filling in the dynamic fields: - for i=1:length(self.nature), - nat=self.nature{i}; + + %check this is acceptable: + for i=1:length(self.nature), + if ~(strcmpi(self.nature{i},'litho') | strcmpi(self.nature{i},'ice') | strcmpi(self.nature{i},'hydro')), + error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')'); + end + end + + %start filling in the dynamic fields: + for i=1:length(self.nature), + nat=self.nature{i}; switch nat case 'ice' @@ -37,5 +37,5 @@ self.addprop('thermalconductivity'); self.addprop('temperateiceconductivity'); - self.addprop('meltingpoint'); + self.addprop('meltingpoint'); self.addprop('beta'); self.addprop('mixed_layer_capacity'); @@ -46,5 +46,5 @@ case 'litho' self.addprop('numlayers'); - self.addprop('radius'); + self.addprop('radius'); self.addprop('viscosity'); self.addprop('lame_lambda'); @@ -70,5 +70,5 @@ for i=1:length(self.nature), - nat=self.nature{i}; + nat=self.nature{i}; switch nat case 'ice' @@ -83,5 +83,5 @@ %water viscosity (N.s/m^2) - self.mu_water=0.001787; + self.mu_water=0.001787; %ice heat capacity cp (J/kg/K) @@ -93,5 +93,5 @@ %ice thermal conductivity (W/m/K) self.thermalconductivity=2.4; - + %wet ice thermal conductivity (W/m/K) self.temperateiceconductivity=.24; @@ -114,9 +114,9 @@ case 'litho' - %we default to a configuration that enables running GIA solutions using giacaron and/or giaivins. + %we default to a configuration that enables running GIA solutions using giacaron and/or giaivins. self.numlayers=2; %center of the earth (approximation, must not be 0), then the lab (lithosphere/asthenosphere boundary) then the surface - %(with 1d3 to avoid numerical singularities) + %(with 1d3 to avoid numerical singularities) self.radius=[1e3;6278*1e3;6378*1e3]; @@ -136,6 +136,6 @@ %ocean water density (kg/m^3) self.rho_water=1023.; - %SLR - self.earth_density= 5512; % average density of the Earth, (kg/m^3) + % average density of the Earth (kg/m^3) + self.earth_density=5512; otherwise @@ -148,5 +148,5 @@ for i=1:length(self.nature), - nat=self.nature{i}; + nat=self.nature{i}; switch nat case 'ice' @@ -183,5 +183,5 @@ fielddisplay(self,'rho_ice','ice density [kg/m^3]'); fielddisplay(self,'rho_water','ocean water density [kg/m^3]'); - fielddisplay(self,'earth_density','Mantle density [kg/m^-3]'); + fielddisplay(self,'earth_density','mantle density [kg/m^3]'); otherwise @@ -193,5 +193,5 @@ for i=1:length(self.nature), - nat=self.nature{i}; + nat=self.nature{i}; switch nat case 'ice' @@ -217,9 +217,9 @@ for i=1:md.materials.numlayers, - if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))), + if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))), error('materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with isburgers choice'); end end - if md.materials.issolid(1)==0 | md.materials.lame_mu(1)==0 + if md.materials.issolid(1)==0 | md.materials.lame_mu(1)==0 error('First layer must be solid (issolid(1) > 0 AND lame_mu(1) > 0). Add a weak inner core if necessary.'); end @@ -242,10 +242,10 @@ function marshall(self,prefix,md,fid) % {{{ - %1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum + %1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum WriteData(fid,prefix,'name','md.materials.nature','data',naturetointeger(self.nature),'format','IntMat','mattype',3); - WriteData(fid,prefix,'name','md.materials.type','data',5,'format','Integer'); %DANGER, this can evolve if you have classes. - - for i=1:length(self.nature), - nat=self.nature{i}; + WriteData(fid,prefix,'name','md.materials.type','data',5,'format','Integer'); %DANGER, this can evolve if you have classes. + + for i=1:length(self.nature), + nat=self.nature{i}; switch nat case 'ice' @@ -271,9 +271,9 @@ WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_lambda','format','DoubleMat','mattype',3); WriteData(fid,prefix,'object',self,'class','materials','fieldname','issolid','format','DoubleMat','mattype',3); - WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3); - WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3); - WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3); - WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3); - WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3); + WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3); + WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3); + WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3); + WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3); + WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3); case 'hydro' WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double'); @@ -288,5 +288,5 @@ function self = extrude(self,md) % {{{ for i=1:length(self.nature), - nat=self.nature{i}; + nat=self.nature{i}; switch nat case 'ice' @@ -297,7 +297,7 @@ end % }}} function savemodeljs(self,fid,modelname) % {{{ - - for i=1:length(self.nature), - nat=self.nature{i}; + + for i=1:length(self.nature), + nat=self.nature{i}; switch nat case 'ice' @@ -324,9 +324,9 @@ writejsdouble(fid,[modelname '.materials.lame_lambda'],self.lame_lambda); writejsdouble(fid,[modelname '.materials.issolid'],self.issolid); - writejsdouble(fid,[modelname '.materials.density'],self.density); - writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity); - writejsdouble(fid,[modelname '.materials.isburgers'],self.isburgers); - writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity); - writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu); + writejsdouble(fid,[modelname '.materials.density'],self.density); + writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity); + writejsdouble(fid,[modelname '.materials.isburgers'],self.isburgers); + writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity); + writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu); case 'hydro' @@ -345,25 +345,26 @@ end end - function intnat = naturetointeger(strnat) % {{{ - intnat=zeros(length(strnat),1); - for i=1:length(strnat), - switch strnat{i}, - case 'damageice' - intnat(i)=1; - case 'estar' - intnat(i)=2; - case 'ice' - intnat(i)=3; - case 'enhancedice' - intnat(i)=4; - %case 'materials' %this case will never happen, kept to ensure equivalent of codes between IoCodeToMaterialsEnumm and IoCodeToNatureEnum - % intnat(i)=5; - case 'litho' - intnat(i)=6; - case 'hydro' - intnat(i)=7; - otherwise - error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')'); - end - end - end % }}} + +function intnat = naturetointeger(strnat) % {{{ + intnat=zeros(length(strnat),1); + for i=1:length(strnat), + switch strnat{i}, + case 'damageice' + intnat(i)=1; + case 'estar' + intnat(i)=2; + case 'ice' + intnat(i)=3; + case 'enhancedice' + intnat(i)=4; + %case 'materials' %this case will never happen, kept to ensure equivalent of codes between IoCodeToMaterialsEnum and IoCodeToNatureEnum + % intnat(i)=5; + case 'litho' + intnat(i)=6; + case 'hydro' + intnat(i)=7; + otherwise + error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')'); + end + end +end % }}} Index: /issm/trunk-jpl/src/m/classes/materials.py =================================================================== --- /issm/trunk-jpl/src/m/classes/materials.py (revision 24755) +++ /issm/trunk-jpl/src/m/classes/materials.py (revision 24756) @@ -4,28 +4,4 @@ from checkfield import checkfield from WriteData import WriteData - - -def naturetointeger(strnat): #{{{ - - intnat = np.zeros(len(strnat)) - for i in range(len(intnat)): - if strnat[i] == 'damageice': - intnat[i] = 1 - elif strnat[i] == 'estar': - intnat[i] = 2 - elif strnat[i] == 'ice': - intnat[i] = 3 - elif strnat[i] == 'enhancedice': - intnat[i] = 4 - elif strnat[i] == 'litho': - intnat[i] = 6 - elif strnat[i] == 'hydro': - intnat[i] = 7 - else: - raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')") - - return intnat - # }}} - class materials(object): @@ -45,6 +21,6 @@ for i in range(len(self.nature)): - if not(self.nature[i] == 'litho' or self.nature[i] == 'ice'): - raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')") + if not(self.nature[i] == 'litho' or self.nature[i] == 'ice' or self.nature[i] == 'hydro'): + raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") #start filling in the dynamic fields: @@ -52,5 +28,4 @@ nat = self.nature[i] if nat == 'ice': - setattr(self, 'rho_ice', 0) setattr(self, 'rho_ice', 0) setattr(self, 'rho_water', 0) @@ -79,57 +54,13 @@ setattr(self, 'density', 0) setattr(self, 'issolid', 0) - else: - raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')") - #set default parameters: + elif nat == 'hydro': + setattr(self, 'rho_ice', 0) + setattr(self, 'rho_water', 0) + setattr(self, 'earth_density', 0) + else: + raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") + + #set default parameters: self.setdefaultparameters() - #}}} - - def __repr__(self): # {{{ - string = " Materials:" - for i in range(len(self.nature)): - nat = self.nature[i] - if nat == 'ice': - string = "%s\n%s" % (string, 'Ice:') - string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "water density [kg / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]")) - string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K")) - string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]")) - string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', or 'NyeH2O'")) - elif nat == 'litho': - string = "%s\n%s" % (string, 'Litho:') - string = "%s\n%s" % (string, fielddisplay(self, 'numlayers', 'number of layers (default 2)')) - string = "%s\n%s" % (string, fielddisplay(self, 'radius', 'array describing the radius for each interface (numlayers + 1) [m]')) - string = "%s\n%s" % (string, fielddisplay(self, 'viscosity', 'array describing each layer''s viscosity (numlayers) [Pa.s]')) - string = "%s\n%s" % (string, fielddisplay(self, 'lame_lambda', 'array describing the lame lambda parameter (numlayers) [Pa]')) - string = "%s\n%s" % (string, fielddisplay(self, 'lame_mu', 'array describing the shear modulus for each layers (numlayers) [Pa]')) - string = "%s\n%s" % (string, fielddisplay(self, 'burgers_viscosity', 'array describing each layer''s transient viscosity, only for Burgers rheologies (numlayers) [Pa.s]')) - string = "%s\n%s" % (string, fielddisplay(self, 'burgers_mu', 'array describing each layer''s transient shear modulus, only for Burgers rheologies (numlayers) [Pa]')) - string = "%s\n%s" % (string, fielddisplay(self, 'isburgers', 'array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)')) - string = "%s\n%s" % (string, fielddisplay(self, 'density', 'array describing each layer''s density (numlayers) [kg / m^3]')) - string = "%s\n%s" % (string, fielddisplay(self, 'issolid', 'array describing whether the layer is solid or liquid (default 1) (numlayers)')) - - else: - raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')") - - return string - #}}} - - def extrude(self, md): # {{{ - for i in range(len(self.nature)): - nat = self.nature[i] - if nat == 'ice': - self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node') - self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element') - return self #}}} @@ -165,5 +96,4 @@ #available: none, paterson and arrhenius self.rheology_law = 'Paterson' - elif nat == 'litho': #we default to a configuration that enables running GIA solutions using giacaron and / or giaivins. @@ -180,9 +110,59 @@ self.density = [5.51 * 1e3, 5.50 * 1e3] # (Pa) #mantle and lithosphere density [kg / m^3] self.issolid = [1, 1] # is layer solid or liquid. - - else: - raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho')") + elif nat == 'hydro': + #ice density (kg / m^3) + self.rho_ice = 917. + #ocean water density (kg / m^3) + self.rho_water = 1023. + #average density of the Earth (kg / m^3) + self.earth_density = 5512 + else: + raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") return self + #}}} + + def __repr__(self): # {{{ + string = " Materials:" + for i in range(len(self.nature)): + nat = self.nature[i] + if nat == 'ice': + string = "%s\n%s" % (string, 'Ice:') + string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]")) + string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "ocean water density [kg / m^3]")) + string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]")) + string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]")) + string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]")) + string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]")) + string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]")) + string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K")) + string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]")) + string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]")) + string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]")) + string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]")) + string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]")) + string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent")) + string = "%s\n%s" % (string, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', or 'NyeH2O'")) + elif nat == 'litho': + string = "%s\n%s" % (string, 'Litho:') + string = "%s\n%s" % (string, fielddisplay(self, 'numlayers', 'number of layers (default 2)')) + string = "%s\n%s" % (string, fielddisplay(self, 'radius', 'array describing the radius for each interface (numlayers + 1) [m]')) + string = "%s\n%s" % (string, fielddisplay(self, 'viscosity', 'array describing each layer''s viscosity (numlayers) [Pa.s]')) + string = "%s\n%s" % (string, fielddisplay(self, 'lame_lambda', 'array describing the lame lambda parameter (numlayers) [Pa]')) + string = "%s\n%s" % (string, fielddisplay(self, 'lame_mu', 'array describing the shear modulus for each layers (numlayers) [Pa]')) + string = "%s\n%s" % (string, fielddisplay(self, 'burgers_viscosity', 'array describing each layer''s transient viscosity, only for Burgers rheologies (numlayers) [Pa.s]')) + string = "%s\n%s" % (string, fielddisplay(self, 'burgers_mu', 'array describing each layer''s transient shear modulus, only for Burgers rheologies (numlayers) [Pa]')) + string = "%s\n%s" % (string, fielddisplay(self, 'isburgers', 'array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)')) + string = "%s\n%s" % (string, fielddisplay(self, 'density', 'array describing each layer''s density (numlayers) [kg / m^3]')) + string = "%s\n%s" % (string, fielddisplay(self, 'issolid', 'array describing whether the layer is solid or liquid (default 1) (numlayers)')) + elif nat == 'hydro': + string = "%s\n%s" % (string, 'Hydro:') + string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]")) + string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "ocean water density [kg / m^3]")) + string = "%s\n%s" % (string, fielddisplay(self, "earth_density", "mantle density [kg / m^3]")) + else: + raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") + + return string #}}} @@ -198,5 +178,4 @@ md = checkfield(md, 'fieldname', 'materials.rheology_n', '>', 0, 'size', [md.mesh.numberofelements]) md = checkfield(md, 'fieldname', 'materials.rheology_law', 'values', ['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', 'NyeH2O']) - elif nat == 'litho': if 'LoveAnalysis' not in analyses: @@ -224,7 +203,10 @@ if (not md.materials.issolid[i]) and (not md.materials.issolid[i + 1]): #if there are at least two consecutive indices that contain issolid = 0 raise RuntimeError("%s%i%s" % ('2 or more adjacent fluid layers detected starting at layer ', i, '. This is not supported yet. Consider merging them.')) - - else: - raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho')") + elif nat == 'hydro': + md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0) + md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0) + md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', [1]) + else: + raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") return md @@ -233,6 +215,6 @@ def marshall(self, prefix, md, fid): # {{{ #1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum - WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 6, 'format', 'Integer') WriteData(fid, prefix, 'name', 'md.materials.nature', 'data', naturetointeger(self.nature), 'format', 'IntMat', 'mattype', 3) + WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 5, 'format', 'Integer') #DANGER, this can evolve if you have classes for i in range(len(self.nature)): @@ -255,5 +237,4 @@ WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_n', 'format', 'DoubleMat', 'mattype', 2) WriteData(fid, prefix, 'data', self.rheology_law, 'name', 'md.materials.rheology_law', 'format', 'String') - elif nat == 'litho': WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'numlayers', 'format', 'Integer') @@ -267,6 +248,42 @@ WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'burgers_viscosity', 'format', 'DoubleMat', 'mattype', 3) WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'burgers_mu', 'format', 'DoubleMat', 'mattype', 3) - - else: - raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')") + elif nat == 'hydro': + WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double') + WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_water', 'format', 'Double') + WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'earth_density', 'format', 'Double') + else: + raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") # }}} + + def extrude(self, md): # {{{ + for i in range(len(self.nature)): + nat = self.nature[i] + if nat == 'ice': + self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node') + self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element') + return self + #}}} + +def naturetointeger(strnat): #{{{ + intnat = np.zeros(len(strnat)) + + for i in range(len(intnat)): + if strnat[i] == 'damageice': + intnat[i] = 1 + elif strnat[i] == 'estar': + intnat[i] = 2 + elif strnat[i] == 'ice': + intnat[i] = 3 + elif strnat[i] == 'enhancedice': + intnat[i] = 4 + # elif strnat[i] == 'materials': + # intnat[i] = 5 #this case will never happen, kept to ensure equivalent of codes between IoCodeToMeterialsEnum and IoCodeToNatureEnum + elif strnat[i] == 'litho': + intnat[i] = 6 + elif strnat[i] == 'hydro': + intnat[i] = 7 + else: + raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')") + + return intnat +# }}}