source: issm/oecreview/Archive/24684-25833/ISSM-24755-24756.diff

Last change on this file was 25834, checked in by Mathieu Morlighem, 4 years ago

CHG: added 24684-25833
File size: 29.1 KB

  • ../trunk-jpl/src/m/classes/materials.py

     
    44from checkfield import checkfield
    55from WriteData import WriteData
    66
    7 
    8 def naturetointeger(strnat):  #{{{
    9 
    10     intnat = np.zeros(len(strnat))
    11     for i in range(len(intnat)):
    12         if strnat[i] == 'damageice':
    13             intnat[i] = 1
    14         elif strnat[i] == 'estar':
    15             intnat[i] = 2
    16         elif strnat[i] == 'ice':
    17             intnat[i] = 3
    18         elif strnat[i] == 'enhancedice':
    19             intnat[i] = 4
    20         elif strnat[i] == 'litho':
    21             intnat[i] = 6
    22         elif strnat[i] == 'hydro':
    23             intnat[i] = 7
    24         else:
    25             raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')")
    26 
    27         return intnat
    28     # }}}
    29 
    30 
    317class materials(object):
    328    """
    339    MATERIALS class definition
    … …  
    4420            self.nature = args
    4521
    4622        for i in range(len(self.nature)):
    47             if not(self.nature[i] == 'litho' or self.nature[i] == 'ice'):
    48                 raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
     23            if not(self.nature[i] == 'litho' or self.nature[i] == 'ice' or self.nature[i] == 'hydro'):
     24                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
    4925
    5026        #start filling in the dynamic fields:
    5127        for i in range(len(self.nature)):
    … …  
    5228            nat = self.nature[i]
    5329            if nat == 'ice':
    5430                setattr(self, 'rho_ice', 0)
    55                 setattr(self, 'rho_ice', 0)
    5631                setattr(self, 'rho_water', 0)
    5732                setattr(self, 'rho_freshwater', 0)
    5833                setattr(self, 'mu_water', 0)
    … …  
    7853                setattr(self, 'isburgers', 0)
    7954                setattr(self, 'density', 0)
    8055                setattr(self, 'issolid', 0)
     56            elif nat == 'hydro':
     57                setattr(self, 'rho_ice', 0)
     58                setattr(self, 'rho_water', 0)
     59                setattr(self, 'earth_density', 0)
    8160            else:
    82                 raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
    83     #set default parameters:
     61                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
     62
     63        #set default parameters:
    8464        self.setdefaultparameters()
    8565    #}}}
    8666
    87     def __repr__(self):  # {{{
    88         string = "   Materials:"
    89         for i in range(len(self.nature)):
    90             nat = self.nature[i]
    91             if nat == 'ice':
    92                 string = "%s\n%s" % (string, 'Ice:')
    93                 string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
    94                 string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "water density [kg / m^3]"))
    95                 string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]"))
    96                 string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]"))
    97                 string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]"))
    98                 string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]"))
    99                 string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]"))
    100                 string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K"))
    101                 string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]"))
    102                 string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]"))
    103                 string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]"))
    104                 string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]"))
    105                 string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]"))
    106                 string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent"))
    107                 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'"))
    108             elif nat == 'litho':
    109                 string = "%s\n%s" % (string, 'Litho:')
    110                 string = "%s\n%s" % (string, fielddisplay(self, 'numlayers', 'number of layers (default 2)'))
    111                 string = "%s\n%s" % (string, fielddisplay(self, 'radius', 'array describing the radius for each interface (numlayers + 1) [m]'))
    112                 string = "%s\n%s" % (string, fielddisplay(self, 'viscosity', 'array describing each layer''s viscosity (numlayers) [Pa.s]'))
    113                 string = "%s\n%s" % (string, fielddisplay(self, 'lame_lambda', 'array describing the lame lambda parameter (numlayers) [Pa]'))
    114                 string = "%s\n%s" % (string, fielddisplay(self, 'lame_mu', 'array describing the shear modulus for each layers (numlayers) [Pa]'))
    115                 string = "%s\n%s" % (string, fielddisplay(self, 'burgers_viscosity', 'array describing each layer''s transient viscosity, only for Burgers rheologies  (numlayers) [Pa.s]'))
    116                 string = "%s\n%s" % (string, fielddisplay(self, 'burgers_mu', 'array describing each layer''s transient shear modulus, only for Burgers rheologies  (numlayers) [Pa]'))
    117                 string = "%s\n%s" % (string, fielddisplay(self, 'isburgers', 'array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'))
    118                 string = "%s\n%s" % (string, fielddisplay(self, 'density', 'array describing each layer''s density (numlayers) [kg / m^3]'))
    119                 string = "%s\n%s" % (string, fielddisplay(self, 'issolid', 'array describing whether the layer is solid or liquid (default 1) (numlayers)'))
    120 
    121             else:
    122                 raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
    123 
    124         return string
    125     #}}}
    126 
    127     def extrude(self, md):  # {{{
    128         for i in range(len(self.nature)):
    129             nat = self.nature[i]
    130             if nat == 'ice':
    131                 self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node')
    132                 self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element')
    133             return self
    134     #}}}
    135 
    13667    def setdefaultparameters(self):  # {{{
    13768        for i in range(len(self.nature)):
    13869            nat = self.nature[i]
    … …  
    16495                #Rheology law: what is the temperature dependence of B with T
    16596                #available: none, paterson and arrhenius
    16697                self.rheology_law = 'Paterson'
    167 
    16898            elif nat == 'litho':
    16999                #we default to a configuration that enables running GIA solutions using giacaron and / or giaivins.
    170100                self.numlayers = 2
    … …  
    179109                self.isburgers = [0, 0]
    180110                self.density = [5.51 * 1e3, 5.50 * 1e3]  # (Pa)  #mantle and lithosphere density [kg / m^3]
    181111                self.issolid = [1, 1]  # is layer solid or liquid.
    182 
     112            elif nat == 'hydro':
     113                #ice density (kg / m^3)
     114                self.rho_ice = 917.
     115                #ocean water density (kg / m^3)
     116                self.rho_water = 1023.
     117                #average density of the Earth (kg / m^3)
     118                self.earth_density = 5512
    183119            else:
    184                 raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho')")
     120                raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
    185121
    186122        return self
    187123    #}}}
    188124
     125    def __repr__(self):  # {{{
     126        string = "   Materials:"
     127        for i in range(len(self.nature)):
     128            nat = self.nature[i]
     129            if nat == 'ice':
     130                string = "%s\n%s" % (string, 'Ice:')
     131                string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
     132                string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "ocean water density [kg / m^3]"))
     133                string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]"))
     134                string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]"))
     135                string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]"))
     136                string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]"))
     137                string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]"))
     138                string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K"))
     139                string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]"))
     140                string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]"))
     141                string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]"))
     142                string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]"))
     143                string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]"))
     144                string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent"))
     145                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'"))
     146            elif nat == 'litho':
     147                string = "%s\n%s" % (string, 'Litho:')
     148                string = "%s\n%s" % (string, fielddisplay(self, 'numlayers', 'number of layers (default 2)'))
     149                string = "%s\n%s" % (string, fielddisplay(self, 'radius', 'array describing the radius for each interface (numlayers + 1) [m]'))
     150                string = "%s\n%s" % (string, fielddisplay(self, 'viscosity', 'array describing each layer''s viscosity (numlayers) [Pa.s]'))
     151                string = "%s\n%s" % (string, fielddisplay(self, 'lame_lambda', 'array describing the lame lambda parameter (numlayers) [Pa]'))
     152                string = "%s\n%s" % (string, fielddisplay(self, 'lame_mu', 'array describing the shear modulus for each layers (numlayers) [Pa]'))
     153                string = "%s\n%s" % (string, fielddisplay(self, 'burgers_viscosity', 'array describing each layer''s transient viscosity, only for Burgers rheologies  (numlayers) [Pa.s]'))
     154                string = "%s\n%s" % (string, fielddisplay(self, 'burgers_mu', 'array describing each layer''s transient shear modulus, only for Burgers rheologies  (numlayers) [Pa]'))
     155                string = "%s\n%s" % (string, fielddisplay(self, 'isburgers', 'array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'))
     156                string = "%s\n%s" % (string, fielddisplay(self, 'density', 'array describing each layer''s density (numlayers) [kg / m^3]'))
     157                string = "%s\n%s" % (string, fielddisplay(self, 'issolid', 'array describing whether the layer is solid or liquid (default 1) (numlayers)'))
     158            elif nat == 'hydro':
     159                string = "%s\n%s" % (string, 'Hydro:')
     160                string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
     161                string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "ocean water density [kg / m^3]"))
     162                string = "%s\n%s" % (string, fielddisplay(self, "earth_density", "mantle density [kg / m^3]"))
     163            else:
     164                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
     165
     166        return string
     167    #}}}
     168
    189169    def checkconsistency(self, md, solution, analyses):  # {{{
    190170        for i in range(len(self.nature)):
    191171            nat = self.nature[i]
    … …  
    197177                md = checkfield(md, 'fieldname', 'materials.rheology_B', '>', 0, 'timeseries', 1, 'NaN', 1, 'Inf', 1)
    198178                md = checkfield(md, 'fieldname', 'materials.rheology_n', '>', 0, 'size', [md.mesh.numberofelements])
    199179                md = checkfield(md, 'fieldname', 'materials.rheology_law', 'values', ['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', 'NyeH2O'])
    200 
    201180            elif nat == 'litho':
    202181                if 'LoveAnalysis' not in analyses:
    203182                    return md
    … …  
    223202                    for i in range(md.materials.numlayers - 1):
    224203                        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
    225204                            raise RuntimeError("%s%i%s" % ('2 or more adjacent fluid layers detected starting at layer ', i, '. This is not supported yet. Consider merging them.'))
    226 
     205            elif nat == 'hydro':
     206                md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0)
     207                md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0)
     208                md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', [1])
    227209            else:
    228                 raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho')")
     210                raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
    229211
    230212        return md
    231213    # }}}
    … …  
    232214
    233215    def marshall(self, prefix, md, fid):  # {{{
    234216        #1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
    235         WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 6, 'format', 'Integer')
    236217        WriteData(fid, prefix, 'name', 'md.materials.nature', 'data', naturetointeger(self.nature), 'format', 'IntMat', 'mattype', 3)
     218        WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 5, 'format', 'Integer') #DANGER, this can evolve if you have classes
    237219
    238220        for i in range(len(self.nature)):
    239221            nat = self.nature[i]
    … …  
    254236                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_B', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts)
    255237                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_n', 'format', 'DoubleMat', 'mattype', 2)
    256238                WriteData(fid, prefix, 'data', self.rheology_law, 'name', 'md.materials.rheology_law', 'format', 'String')
    257 
    258239            elif nat == 'litho':
    259240                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'numlayers', 'format', 'Integer')
    260241                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'radius', 'format', 'DoubleMat', 'mattype', 3)
    … …  
    266247                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'isburgers', 'format', 'DoubleMat', 'mattype', 3)
    267248                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'burgers_viscosity', 'format', 'DoubleMat', 'mattype', 3)
    268249                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'burgers_mu', 'format', 'DoubleMat', 'mattype', 3)
    269 
     250            elif nat == 'hydro':
     251                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double')
     252                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_water', 'format', 'Double')
     253                WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'earth_density', 'format', 'Double')
    270254            else:
    271                 raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')")
     255                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
    272256    # }}}
     257
     258    def extrude(self, md):  # {{{
     259        for i in range(len(self.nature)):
     260            nat = self.nature[i]
     261            if nat == 'ice':
     262                self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node')
     263                self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element')
     264            return self
     265    #}}}
     266
     267def naturetointeger(strnat):  #{{{
     268    intnat = np.zeros(len(strnat))
     269
     270    for i in range(len(intnat)):
     271        if strnat[i] == 'damageice':
     272            intnat[i] = 1
     273        elif strnat[i] == 'estar':
     274            intnat[i] = 2
     275        elif strnat[i] == 'ice':
     276            intnat[i] = 3
     277        elif strnat[i] == 'enhancedice':
     278            intnat[i] = 4
     279        # elif strnat[i] == 'materials':
     280        #     intnat[i] = 5 #this case will never happen, kept to ensure equivalent of codes between IoCodeToMeterialsEnum and IoCodeToNatureEnum
     281        elif strnat[i] == 'litho':
     282            intnat[i] = 6
     283        elif strnat[i] == 'hydro':
     284            intnat[i] = 7
     285        else:
     286            raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
     287
     288    return intnat
     289# }}}

  • ../trunk-jpl/src/m/classes/materials.m

     
    44%      materials=materials();
    55
    66classdef materials < dynamicprops
    7         properties (SetAccess=public) 
     7        properties (SetAccess=public)
    88                nature={};
    99                %all properties are dynamic.
    1010        end
    … …  
    1212                function self = materials(varargin) % {{{
    1313                        if nargin==0
    1414                                self.nature={'ice'};
    15                         else 
     15                        else
    1616                                self.nature=varargin;
    1717                        end
    18                        
    19                         %check this is acceptable: 
     18
     19                        %check this is acceptable:
    2020                        for i=1:length(self.nature),
    21                                 if ~(strcmpi(self.nature{i},'litho') | strcmpi(self.nature{i},'ice') | strcmpi(self.nature{i},'hydro')), 
     21                                if ~(strcmpi(self.nature{i},'litho') | strcmpi(self.nature{i},'ice') | strcmpi(self.nature{i},'hydro')),
    2222                                        error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
    2323                                end
    2424                        end
    25                        
    26                         %start filling in the dynamic fields: 
     25
     26                        %start filling in the dynamic fields:
    2727                        for i=1:length(self.nature),
    28                                 nat=self.nature{i}; 
     28                                nat=self.nature{i};
    2929                                switch nat
    3030                                case 'ice'
    3131                                        self.addprop('rho_ice');
    … …  
    3636                                        self.addprop('latentheat');
    3737                                        self.addprop('thermalconductivity');
    3838                                        self.addprop('temperateiceconductivity');
    39                                     self.addprop('meltingpoint');
     39                                        self.addprop('meltingpoint');
    4040                                        self.addprop('beta');
    4141                                        self.addprop('mixed_layer_capacity');
    4242                                        self.addprop('thermal_exchange_velocity');
    … …  
    4545                                        self.addprop('rheology_law');
    4646                                case 'litho'
    4747                                        self.addprop('numlayers');
    48                                     self.addprop('radius');
     48                                        self.addprop('radius');
    4949                                        self.addprop('viscosity');
    5050                                        self.addprop('lame_lambda');
    5151                                        self.addprop('lame_mu');
    … …  
    6969                function self = setdefaultparameters(self) % {{{
    7070
    7171                        for i=1:length(self.nature),
    72                                 nat=self.nature{i}; 
     72                                nat=self.nature{i};
    7373                                switch nat
    7474                                case 'ice'
    7575                                        %ice density (kg/m^3)
    … …  
    8282                                        self.rho_freshwater=1000.;
    8383
    8484                                        %water viscosity (N.s/m^2)
    85                                         self.mu_water=0.001787; 
     85                                        self.mu_water=0.001787;
    8686
    8787                                        %ice heat capacity cp (J/kg/K)
    8888                                        self.heatcapacity=2093.;
    … …  
    9292
    9393                                        %ice thermal conductivity (W/m/K)
    9494                                        self.thermalconductivity=2.4;
    95                                        
     95
    9696                                        %wet ice thermal conductivity (W/m/K)
    9797                                        self.temperateiceconductivity=.24;
    9898
    … …  
    113113                                        self.rheology_law='Paterson';
    114114
    115115                                case 'litho'
    116                                         %we default to a configuration that enables running GIA solutions using giacaron and/or giaivins. 
     116                                        %we default to a configuration that enables running GIA solutions using giacaron and/or giaivins.
    117117                                        self.numlayers=2;
    118118
    119119                                        %center of the earth (approximation, must not be 0), then the lab (lithosphere/asthenosphere boundary) then the surface
    120                                         %(with 1d3 to avoid numerical singularities) 
     120                                        %(with 1d3 to avoid numerical singularities)
    121121                                        self.radius=[1e3;6278*1e3;6378*1e3];
    122122
    123123                                        self.viscosity=[1e21;1e40]; %mantle and lithosphere viscosity (respectively) [Pa.s]
    … …  
    135135
    136136                                        %ocean water density (kg/m^3)
    137137                                        self.rho_water=1023.;
    138                                         %SLR
    139                                         self.earth_density= 5512;  % average density of the Earth, (kg/m^3)
     138                                        % average density of the Earth (kg/m^3)
     139                                        self.earth_density=5512;
    140140
    141141                                otherwise
    142142                                        error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
    … …  
    147147                        disp(sprintf('   Materials:'));
    148148
    149149                        for i=1:length(self.nature),
    150                                 nat=self.nature{i}; 
     150                                nat=self.nature{i};
    151151                                switch nat
    152152                                case 'ice'
    153153                                        disp(sprintf('   \nIce:'));
    … …  
    182182                                        disp(sprintf('   \nHydro:'));
    183183                                        fielddisplay(self,'rho_ice','ice density [kg/m^3]');
    184184                                        fielddisplay(self,'rho_water','ocean water density [kg/m^3]');
    185                                         fielddisplay(self,'earth_density','Mantle density [kg/m^-3]');
     185                                        fielddisplay(self,'earth_density','mantle density [kg/m^3]');
    186186
    187187                                otherwise
    188188                                        error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
    … …  
    192192                function md = checkconsistency(self,md,solution,analyses) % {{{
    193193
    194194                        for i=1:length(self.nature),
    195                                 nat=self.nature{i}; 
     195                                nat=self.nature{i};
    196196                                switch nat
    197197                                case 'ice'
    198198                                        md = checkfield(md,'fieldname','materials.rho_ice','>',0);
    … …  
    216216                                        md = checkfield(md,'fieldname','materials.burgers_mu','Inf',1,'size',[md.materials.numlayers 1],'>=',0);
    217217
    218218                                        for i=1:md.materials.numlayers,
    219                                                 if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))), 
     219                                                if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))),
    220220                                                        error('materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with isburgers choice');
    221221                                                end
    222222                                        end
    223                                         if md.materials.issolid(1)==0 | md.materials.lame_mu(1)==0 
     223                                        if md.materials.issolid(1)==0 | md.materials.lame_mu(1)==0
    224224                                                        error('First layer must be solid (issolid(1) > 0 AND lame_mu(1) > 0). Add a weak inner core if necessary.');
    225225                                        end
    226226                                        ind=find(md.materials.issolid==0);
    … …  
    241241                end % }}}
    242242                function marshall(self,prefix,md,fid) % {{{
    243243
    244                         %1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum 
     244                        %1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
    245245                        WriteData(fid,prefix,'name','md.materials.nature','data',naturetointeger(self.nature),'format','IntMat','mattype',3);
    246                         WriteData(fid,prefix,'name','md.materials.type','data',5,'format','Integer'); %DANGER, this can evolve if you have classes. 
     246                        WriteData(fid,prefix,'name','md.materials.type','data',5,'format','Integer'); %DANGER, this can evolve if you have classes.
    247247
    248248                        for i=1:length(self.nature),
    249                                 nat=self.nature{i}; 
     249                                nat=self.nature{i};
    250250                                switch nat
    251251                                case 'ice'
    252252                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double');
    … …  
    270270                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_mu','format','DoubleMat','mattype',3);
    271271                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_lambda','format','DoubleMat','mattype',3);
    272272                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','issolid','format','DoubleMat','mattype',3);
    273                                         WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3); 
    274                                         WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3); 
    275                                         WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3); 
    276                                         WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3); 
    277                                         WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3); 
     273                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3);
     274                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3);
     275                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3);
     276                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3);
     277                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3);
    278278                                case 'hydro'
    279279                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double');
    280280                                        WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_water','format','Double');
    … …  
    287287                end % }}}
    288288                function self = extrude(self,md) % {{{
    289289                        for i=1:length(self.nature),
    290                                 nat=self.nature{i}; 
     290                                nat=self.nature{i};
    291291                                switch nat
    292292                                case 'ice'
    293293                                        self.rheology_B=project3d(md,'vector',self.rheology_B,'type','node');
    … …  
    296296                        end
    297297                end % }}}
    298298                function savemodeljs(self,fid,modelname) % {{{
    299        
     299
    300300                        for i=1:length(self.nature),
    301                                 nat=self.nature{i}; 
     301                                nat=self.nature{i};
    302302                                switch nat
    303303                                case 'ice'
    304304                                        writejsdouble(fid,[modelname '.materials.rho_ice'],self.rho_ice);
    … …  
    323323                                        writejsdouble(fid,[modelname '.materials.lame_mu'],self.lame_mu);
    324324                                        writejsdouble(fid,[modelname '.materials.lame_lambda'],self.lame_lambda);
    325325                                        writejsdouble(fid,[modelname '.materials.issolid'],self.issolid);
    326                                         writejsdouble(fid,[modelname '.materials.density'],self.density); 
    327                                         writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity); 
    328                                         writejsdouble(fid,[modelname '.materials.isburgers'],self.isburgers); 
    329                                         writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity); 
    330                                         writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu); 
     326                                        writejsdouble(fid,[modelname '.materials.density'],self.density);
     327                                        writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity);
     328                                        writejsdouble(fid,[modelname '.materials.isburgers'],self.isburgers);
     329                                        writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity);
     330                                        writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu);
    331331
    332332                                case 'hydro'
    333333                                        writejsdouble(fid,[modelname '.materials.rho_ice'],self.rho_ice);
    … …  
    344344                end % }}}
    345345        end
    346346end
    347                 function intnat = naturetointeger(strnat) % {{{
    348                         intnat=zeros(length(strnat),1);
    349                         for i=1:length(strnat),
    350                                 switch strnat{i},
    351                                         case 'damageice'
    352                                         intnat(i)=1;
    353                                 case 'estar'
    354                                         intnat(i)=2;
    355                                 case 'ice'
    356                                         intnat(i)=3;
    357                                 case 'enhancedice'
    358                                         intnat(i)=4;
    359                                 %case 'materials' %this case will never happen, kept to ensure equivalent of codes between IoCodeToMaterialsEnumm and IoCodeToNatureEnum
    360                                 %       intnat(i)=5;
    361                                 case 'litho'
    362                                         intnat(i)=6;
    363                                 case 'hydro'
    364                                         intnat(i)=7;
    365                                 otherwise
    366                                         error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
    367                                 end
    368                         end
    369                 end % }}}
     347
     348function intnat = naturetointeger(strnat) % {{{
     349        intnat=zeros(length(strnat),1);
     350        for i=1:length(strnat),
     351                switch strnat{i},
     352                        case 'damageice'
     353                        intnat(i)=1;
     354                case 'estar'
     355                        intnat(i)=2;
     356                case 'ice'
     357                        intnat(i)=3;
     358                case 'enhancedice'
     359                        intnat(i)=4;
     360                %case 'materials' %this case will never happen, kept to ensure equivalent of codes between IoCodeToMaterialsEnum and IoCodeToNatureEnum
     361                %       intnat(i)=5;
     362                case 'litho'
     363                        intnat(i)=6;
     364                case 'hydro'
     365                        intnat(i)=7;
     366                otherwise
     367                        error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
     368                end
     369        end
     370end % }}}
Note: See TracBrowser for help on using the repository browser.