Changeset 23711


Ignore:
Timestamp:
02/11/19 08:43:03 (6 years ago)
Author:
bdef
Message:

BUG:fixing fieldname recognition in parseresults

Location:
issm/trunk-jpl/src/py3
Files:
7 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/py3/classes/matdamageice.py

    r23670 r23711  
    2121                self.thermalconductivity       = 0.
    2222                self.temperateiceconductivity  = 0.
     23    self.effectiveconductivity_averaging = 0.
    2324                self.meltingpoint              = 0.
    2425                self.beta                      = 0.
    … …  
    2930                self.rheology_law              = ''
    3031
    31                 #giaivins: 
     32                #giaivins:
    3233                self.lithosphere_shear_modulus  = 0.
    3334                self.lithosphere_density        = 0.
    3435                self.mantle_shear_modulus       = 0.
    3536                self.mantle_density             = 0.
    36                
     37
    3738                #SLR
    3839                self.earth_density= 5512;  # average density of the Earth, (kg/m^3)
    3940
    40 
    4141                self.setdefaultparameters()
    4242                #}}}
     43
    4344        def __repr__(self): # {{{
    4445                string="   Materials:"
    45 
    4646                string="%s\n%s"%(string,fielddisplay(self,"rho_ice","ice density [kg/m^3]"))
    4747                string="%s\n%s"%(string,fielddisplay(self,"rho_water","water density [kg/m^3]"))
    … …  
    5151                string="%s\n%s"%(string,fielddisplay(self,"thermalconductivity","ice thermal conductivity [W/m/K]"))
    5252                string="%s\n%s"%(string,fielddisplay(self,"temperateiceconductivity","temperate ice thermal conductivity [W/m/K]"))
     53                string="%s\n%s"%(string,fielddisplay(self,"effectiveconductivity_averaging","computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)"))
    5354                string="%s\n%s"%(string,fielddisplay(self,"meltingpoint","melting point of ice at 1atm in K"))
    5455                string="%s\n%s"%(string,fielddisplay(self,"latentheat","latent heat of fusion [J/m^3]"))
    … …  
    6465                string="%s\n%s"%(string,fielddisplay(self,"mantle_density","Mantle density [g/cm^-3]"))
    6566                string="%s\n%s"%(string,fielddisplay(self,"earth_density","Mantle density [kg/m^-3]"))
    66 
    67 
    6867                return string
    6968                #}}}
     69
    7070        def extrude(self,md): # {{{
    7171                self.rheology_B=project3d(md,'vector',self.rheology_B,'type','node')
    … …  
    7373                return self
    7474        #}}}
     75
    7576        def setdefaultparameters(self): # {{{
    7677                #ice density (kg/m^3)
    7778                self.rho_ice=917.
    78 
    7979                #ocean water density (kg/m^3)
    8080                self.rho_water=1023.
    81 
    8281                #fresh water density (kg/m^3)
    8382                self.rho_freshwater=1000.
    84 
    8583                #water viscosity (N.s/m^2)
    86                 self.mu_water=0.001787 
    87 
     84                self.mu_water=0.001787
    8885                #ice heat capacity cp (J/kg/K)
    8986                self.heatcapacity=2093.
    90 
    9187                #ice latent heat of fusion L (J/kg)
    9288                self.latentheat=3.34*10**5
    93 
    9489                #ice thermal conductivity (W/m/K)
    9590                self.thermalconductivity=2.4
    96 
    9791                #temperate ice thermal conductivity (W/m/K)
    9892                self.temperateiceconductivity=0.24
    99 
     93    #computation of effective conductivity
     94    self.effectiveconductivity_averaging=1
    10095                #the melting point of ice at 1 atmosphere of pressure in K
    10196                self.meltingpoint=273.15
    102 
    10397                #rate of change of melting point with pressure (K/Pa)
    10498                self.beta=9.8*10**-8
    105 
    10699                #mixed layer (ice-water interface) heat capacity (J/kg/K)
    107100                self.mixed_layer_capacity=3974.
    108 
    109101                #thermal exchange velocity (ice-water interface) (m/s)
    110102                self.thermal_exchange_velocity=1.00*10**-4
    111 
    112103                #Rheology law: what is the temperature dependence of B with T
    113104                #available: none, paterson and arrhenius
    … …  
    119110                self.mantle_shear_modulus       = 1.45*10**11 # (Pa)
    120111                self.mantle_density             = 3.34        # (g/cm^-3)
    121                
     112
    122113                #SLR
    123114                self.earth_density= 5512;  #average density of the Earth, (kg/m^3)
    124 
    125 
    126115                return self
    127116                #}}}
     117
    128118        def checkconsistency(self,md,solution,analyses):    # {{{
    129119                md = checkfield(md,'fieldname','materials.rho_ice','>',0)
    … …  
    134124                md = checkfield(md,'fieldname','materials.rheology_n','>',0,'size',[md.mesh.numberofelements])
    135125                md = checkfield(md,'fieldname','materials.rheology_law','values',['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson','Arrhenius','LliboutryDuval'])
     126                md = checkfield(md,'fieldname','materials.effectiveconductivity_averaging','numel',[1],'values',[0,1,2])
    136127                md = checkfield(md,'fieldname','materials.lithosphere_shear_modulus','>',0,'numel',[1]);
    137128                md = checkfield(md,'fieldname','materials.lithosphere_density','>',0,'numel',[1]);
    … …  
    139130                md = checkfield(md,'fieldname','materials.mantle_density','>',0,'numel',[1]);
    140131                md = checkfield(md,'fieldname','materials.earth_density','>',0,'numel',[1]);
    141 
    142132                return md
    143133        # }}}
     134
    144135        def marshall(self,prefix,md,fid):    # {{{
    145136                WriteData(fid,prefix,'name','md.materials.type','data',1,'format','Integer');
    … …  
    152143                WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermalconductivity','format','Double')
    153144                WriteData(fid,prefix,'object',self,'class','materials','fieldname','temperateiceconductivity','format','Double')
     145                WriteData(fid,prefix,'object',self,'class','materials','fieldname','effectiveconductivity_averaging','format','Integer')
    154146                WriteData(fid,prefix,'object',self,'class','materials','fieldname','meltingpoint','format','Double')
    155147                WriteData(fid,prefix,'object',self,'class','materials','fieldname','beta','format','Double')
    … …  
    159151                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2)
    160152                WriteData(fid,prefix,'data',self.rheology_law,'name','md.materials.rheology_law','format','String')
    161 
    162153                WriteData(fid,prefix,'object',self,'class','materials','fieldname','lithosphere_shear_modulus','format','Double');
    163154                WriteData(fid,prefix,'object',self,'class','materials','fieldname','lithosphere_density','format','Double','scale',10.**3.);

  • issm/trunk-jpl/src/py3/classes/matenhancedice.py

    r23677 r23711  
    2121                self.thermalconductivity       = 0.
    2222                self.temperateiceconductivity  = 0.
     23    self.effectiveconductivity_averaging = 0.
    2324                self.meltingpoint              = 0.
    2425                self.beta                      = 0.
    2526                self.mixed_layer_capacity      = 0.
    2627                self.thermal_exchange_velocity = 0.
    27                 self.rheology_E               = float('NaN')
     28                self.rheology_E                                                         = float('NaN')
    2829                self.rheology_B                = float('NaN')
    2930                self.rheology_n                = float('NaN')
    3031                self.rheology_law              = ''
    3132
    32                 #giaivins: 
     33                #giaivins:
    3334                self.lithosphere_shear_modulus  = 0.
    3435                self.lithosphere_density        = 0.
    3536                self.mantle_shear_modulus       = 0.
    3637                self.mantle_density             = 0.
    37                
     38
    3839                #SLR
    3940                self.earth_density= 0  # average density of the Earth, (kg/m^3)
    … …  
    4142                self.setdefaultparameters()
    4243                #}}}
     44
    4345        def __repr__(self): # {{{
    4446                string="   Materials:"
    45 
    4647                string="%s\n%s"%(string,fielddisplay(self,"rho_ice","ice density [kg/m^3]"))
    4748                string="%s\n%s"%(string,fielddisplay(self,"rho_water","water density [kg/m^3]"))
    … …  
    5152                string="%s\n%s"%(string,fielddisplay(self,"thermalconductivity","ice thermal conductivity [W/m/K]"))
    5253                string="%s\n%s"%(string,fielddisplay(self,"temperateiceconductivity","temperate ice thermal conductivity [W/m/K]"))
     54                string="%s\n%s"%(string,fielddisplay(self,"effectiveconductivity_averaging","computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)"))
    5355                string="%s\n%s"%(string,fielddisplay(self,"meltingpoint","melting point of ice at 1atm in K"))
    5456                string="%s\n%s"%(string,fielddisplay(self,"latentheat","latent heat of fusion [J/m^3]"))
    … …  
    6567                string="%s\n%s"%(string,fielddisplay(self,"mantle_density","Mantle density [g/cm^-3]"))
    6668                string="%s\n%s"%(string,fielddisplay(self,"earth_density","Mantle density [kg/m^-3]"))
    67 
    6869                return string
    6970                #}}}
     71
    7072        def extrude(self,md): # {{{
    7173                self.rheology_E=project3d(md,'vector',self.rheology_E,'type','node')
    … …  
    7476                return self
    7577        #}}}
     78
    7679        def setdefaultparameters(self): # {{{
    7780                #ice density (kg/m^3)
    7881                self.rho_ice=917.
    79 
    8082                #ocean water density (kg/m^3)
    8183                self.rho_water=1023.
    82 
    8384                #fresh water density (kg/m^3)
    8485                self.rho_freshwater=1000.
    85 
    8686                #water viscosity (N.s/m^2)
    87                 self.mu_water=0.001787 
    88 
     87                self.mu_water=0.001787
    8988                #ice heat capacity cp (J/kg/K)
    9089                self.heatcapacity=2093.
    91 
    9290                #ice latent heat of fusion L (J/kg)
    9391                self.latentheat=3.34*10**5
    94 
    9592                #ice thermal conductivity (W/m/K)
    9693                self.thermalconductivity=2.4
    97 
    9894                #temperate ice thermal conductivity (W/m/K)
    9995                self.temperateiceconductivity=0.24
    100 
     96    #computation of effective conductivity
     97    self.effectiveconductivity_averaging=1
    10198                #the melting point of ice at 1 atmosphere of pressure in K
    10299                self.meltingpoint=273.15
    103 
    104100                #rate of change of melting point with pressure (K/Pa)
    105101                self.beta=9.8*10**-8
    106 
    107102                #mixed layer (ice-water interface) heat capacity (J/kg/K)
    108103                self.mixed_layer_capacity=3974.
    109 
    110104                #thermal exchange velocity (ice-water interface) (m/s)
    111105                self.thermal_exchange_velocity=1.00*10**-4
    112 
    113106                #Rheology law: what is the temperature dependence of B with T
    114107                #available: none, paterson and arrhenius
    … …  
    120113                self.mantle_shear_modulus       = 1.45*10**11 # (Pa)
    121114                self.mantle_density             = 3.34        # (g/cm^-3)
    122                
     115
    123116                #SLR
    124117                self.earth_density= 5512  #average density of the Earth, (kg/m^3)
    … …  
    126119                return self
    127120                #}}}
     121
    128122        def checkconsistency(self,md,solution,analyses):    # {{{
    129123                md = checkfield(md,'fieldname','materials.rho_ice','>',0)
    … …  
    135129                md = checkfield(md,'fieldname','materials.rheology_n','>',0,'size',[md.mesh.numberofelements])
    136130                md = checkfield(md,'fieldname','materials.rheology_law','values',['None','BuddJacka','Cuffey','CuffeyTemperate','Paterson','Arrhenius','LliboutryDuval'])
     131                md = checkfield(md,'fieldname','materials.effectiveconductivity_averaging','numel',[1],'values',[0,1,2])
    137132
    138133                if 'GiaAnalysis' in analyses:
    … …  
    145140                return md
    146141        # }}}
     142
    147143        def marshall(self,prefix,md,fid):    # {{{
    148144                WriteData(fid,prefix,'name','md.materials.type','data',4,'format','Integer')
    … …  
    155151                WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermalconductivity','format','Double')
    156152                WriteData(fid,prefix,'object',self,'class','materials','fieldname','temperateiceconductivity','format','Double')
     153                WriteData(fid,prefix,'object',self,'class','materials','fieldname','effectiveconductivity_averaging','format','Integer')
    157154                WriteData(fid,prefix,'object',self,'class','materials','fieldname','meltingpoint','format','Double')
    158155                WriteData(fid,prefix,'object',self,'class','materials','fieldname','beta','format','Double')
    … …  
    163160                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2)
    164161                WriteData(fid,prefix,'data',self.rheology_law,'name','md.materials.rheology_law','format','String')
    165 
    166162                WriteData(fid,prefix,'object',self,'class','materials','fieldname','lithosphere_shear_modulus','format','Double')
    167163                WriteData(fid,prefix,'object',self,'class','materials','fieldname','lithosphere_density','format','Double','scale',10^3)

  • issm/trunk-jpl/src/py3/classes/materials.py

    r23677 r23711  
    44from checkfield import checkfield
    55from WriteData import WriteData
    6                
     6
    77def naturetointeger(strnat): #{{{
    8    
    9     intnat=np.zeros(len(strnat))
    10     for i in range(len(intnat)):
    11         if strnat[i]=='damageice':
    12             intnat[i]=1
    13         elif strnat[i]=='estar':
    14             intnat[i]=2
    15         elif strnat[i]=='ice':
    16             intnat[i]=3
    17         elif strnat[i]=='enhancedice':
    18             intnat[i]=4
    19         elif strnat[i]=='litho':
    20             intnat[i]=5
    21         else:
    22             raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')");
    23    
    24     return intnat
     8
     9        intnat=np.zeros(len(strnat))
     10        for i in range(len(intnat)):
     11                if strnat[i]=='damageice':
     12                        intnat[i]=1
     13                elif strnat[i]=='estar':
     14                        intnat[i]=2
     15                elif strnat[i]=='ice':
     16                        intnat[i]=3
     17                elif strnat[i]=='enhancedice':
     18                        intnat[i]=4
     19                elif strnat[i]=='litho':
     20                        intnat[i]=5
     21                else:
     22                        raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')");
     23
     24                return intnat
    2525# }}}
     26
    2627class materials(object):
    2728        """
    … …  
    3334
    3435        def __init__(self,*args): # {{{
    35                
    36                 self.nature                    = []
    37 
    38                 if not len(args):
    39                     self.nature=['ice']
    40                 else:
    41                     self.nature=args
    42 
    43                 for i in range(len(self.nature)):
    44                     if not(self.nature[i] == 'litho' or self.nature[i]=='ice'):
    45                         raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
    46                    
    47                 #start filling in the dynamic fields:
    48                 for i in range(len(self.nature)):
    49                     nat=self.nature[i];
    50                     if nat=='ice':
    51                         setattr(self,'rho_ice',0)
    52                         setattr(self,'rho_ice',0);
    53                         setattr(self,'rho_water',0);
    54                         setattr(self,'rho_freshwater',0);
    55                         setattr(self,'mu_water',0);
    56                         setattr(self,'heatcapacity',0);
    57                         setattr(self,'latentheat',0);
    58                         setattr(self,'thermalconductivity',0);
    59                         setattr(self,'temperateiceconductivity',0);
    60                         setattr(self,'meltingpoint',0);
    61                         setattr(self,'beta',0);
    62                         setattr(self,'mixed_layer_capacity',0);
    63                         setattr(self,'thermal_exchange_velocity',0);
    64                         setattr(self,'rheology_B',0);
    65                         setattr(self,'rheology_n',0);
    66                         setattr(self,'rheology_law',0);
    67                     elif nat=='litho':
    68                         setattr(self,'numlayers',0);
    69                         setattr(self,'radius',0);
    70                         setattr(self,'viscosity',0);
    71                         setattr(self,'lame_lambda',0);
    72                         setattr(self,'lame_mu',0);
    73                         setattr(self,'burgers_viscosity',0);
    74                         setattr(self,'burgers_mu',0);
    75                         setattr(self,'isburgers',0);
    76                         setattr(self,'density',0);
    77                         setattr(self,'issolid',0);
    78                     else:
    79                         raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')");
    80                 #set default parameters:
     36                self.nature                    = []
     37                if not len(args):
     38                        self.nature=['ice']
     39                else:
     40                        self.nature=args
     41
     42                for i in range(len(self.nature)):
     43                        if not(self.nature[i] == 'litho' or self.nature[i]=='ice'):
     44                                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
     45
     46                #start filling in the dynamic fields:
     47                for i in range(len(self.nature)):
     48                        nat=self.nature[i];
     49                        if nat=='ice':
     50                                setattr(self,'rho_ice',0)
     51                                setattr(self,'rho_ice',0);
     52                                setattr(self,'rho_water',0);
     53                                setattr(self,'rho_freshwater',0);
     54                                setattr(self,'mu_water',0);
     55                                setattr(self,'heatcapacity',0);
     56                                setattr(self,'latentheat',0);
     57                                setattr(self,'thermalconductivity',0);
     58                                setattr(self,'temperateiceconductivity',0);
     59                                setattr(self,'meltingpoint',0);
     60                                setattr(self,'beta',0);
     61                                setattr(self,'mixed_layer_capacity',0);
     62                                setattr(self,'thermal_exchange_velocity',0);
     63                                setattr(self,'rheology_B',0);
     64                                setattr(self,'rheology_n',0);
     65                                setattr(self,'rheology_law',0);
     66                        elif nat=='litho':
     67                                setattr(self,'numlayers',0);
     68                                setattr(self,'radius',0);
     69                                setattr(self,'viscosity',0);
     70                                setattr(self,'lame_lambda',0);
     71                                setattr(self,'lame_mu',0);
     72                                setattr(self,'burgers_viscosity',0);
     73                                setattr(self,'burgers_mu',0);
     74                                setattr(self,'isburgers',0);
     75                                setattr(self,'density',0);
     76                                setattr(self,'issolid',0);
     77                        else:
     78                                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')");
     79                        #set default parameters:
    8180                self.setdefaultparameters()
    8281                #}}}
     82
    8383        def __repr__(self): # {{{
    8484                string="   Materials:"
    85                 for i in range(len(self.nature)):
    86                     nat=self.nature[i];
    87                     if nat=='ice':
    88                         string="%s\n%s"%(string,'Ice:');
    89                         string="%s\n%s"%(string,fielddisplay(self,"rho_ice","ice density [kg/m^3]"))
    90                         string="%s\n%s"%(string,fielddisplay(self,"rho_water","water density [kg/m^3]"))
    91                         string="%s\n%s"%(string,fielddisplay(self,"rho_freshwater","fresh water density [kg/m^3]"))
    92                         string="%s\n%s"%(string,fielddisplay(self,"mu_water","water viscosity [N s/m^2]"))
    93                         string="%s\n%s"%(string,fielddisplay(self,"heatcapacity","heat capacity [J/kg/K]"))
    94                         string="%s\n%s"%(string,fielddisplay(self,"thermalconductivity","ice thermal conductivity [W/m/K]"))
    95                         string="%s\n%s"%(string,fielddisplay(self,"temperateiceconductivity","temperate ice thermal conductivity [W/m/K]"))
    96                         string="%s\n%s"%(string,fielddisplay(self,"meltingpoint","melting point of ice at 1atm in K"))
    97                         string="%s\n%s"%(string,fielddisplay(self,"latentheat","latent heat of fusion [J/m^3]"))
    98                         string="%s\n%s"%(string,fielddisplay(self,"beta","rate of change of melting point with pressure [K/Pa]"))
    99                         string="%s\n%s"%(string,fielddisplay(self,"mixed_layer_capacity","mixed layer capacity [W/kg/K]"))
    100                         string="%s\n%s"%(string,fielddisplay(self,"thermal_exchange_velocity","thermal exchange velocity [m/s]"))
    101                         string="%s\n%s"%(string,fielddisplay(self,"rheology_B","flow law parameter [Pa s^(1/n)]"))
    102                         string="%s\n%s"%(string,fielddisplay(self,"rheology_n","Glen's flow law exponent"))
    103                         string="%s\n%s"%(string,fielddisplay(self,"rheology_law","law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius' or 'LliboutryDuval'"))
    104                     elif nat=='litho':
    105                         string="%s\n%s"%(string,'Litho:');
    106                         string="%s\n%s"%(string,fielddisplay(self,'numlayers','number of layers (default 2)'))
    107                         string="%s\n%s"%(string,fielddisplay(self,'radius','array describing the radius for each interface (numlayers+1) [m]'))
    108                         string="%s\n%s"%(string,fielddisplay(self,'viscosity','array describing each layer''s viscosity (numlayers) [Pa.s]'))
    109                         string="%s\n%s"%(string,fielddisplay(self,'lame_lambda','array describing the lame lambda parameter (numlayers) [Pa]'))
    110                         string="%s\n%s"%(string,fielddisplay(self,'lame_mu','array describing the shear modulus for each layers (numlayers) [Pa]'))
    111                         string="%s\n%s"%(string,fielddisplay(self,'burgers_viscosity','array describing each layer''s transient viscosity, only for Burgers rheologies  (numlayers) [Pa.s]'))
    112                         string="%s\n%s"%(string,fielddisplay(self,'burgers_mu','array describing each layer''s transient shear modulus, only for Burgers rheologies  (numlayers) [Pa]'))
    113                         string="%s\n%s"%(string,fielddisplay(self,'isburgers','array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'))
    114                         string="%s\n%s"%(string,fielddisplay(self,'density','array describing each layer''s density (numlayers) [kg/m^3]'))
    115                         string="%s\n%s"%(string,fielddisplay(self,'issolid','array describing whether the layer is solid or liquid (default 1) (numlayers)'))
    116 
    117                     else:
    118                         raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')");
     85                for i in range(len(self.nature)):
     86                        nat=self.nature[i];
     87                        if nat=='ice':
     88                                string="%s\n%s"%(string,'Ice:');
     89                                string="%s\n%s"%(string,fielddisplay(self,"rho_ice","ice density [kg/m^3]"))
     90                                string="%s\n%s"%(string,fielddisplay(self,"rho_water","water density [kg/m^3]"))
     91                                string="%s\n%s"%(string,fielddisplay(self,"rho_freshwater","fresh water density [kg/m^3]"))
     92                                string="%s\n%s"%(string,fielddisplay(self,"mu_water","water viscosity [N s/m^2]"))
     93                                string="%s\n%s"%(string,fielddisplay(self,"heatcapacity","heat capacity [J/kg/K]"))
     94                                string="%s\n%s"%(string,fielddisplay(self,"thermalconductivity","ice thermal conductivity [W/m/K]"))
     95                                string="%s\n%s"%(string,fielddisplay(self,"temperateiceconductivity","temperate ice thermal conductivity [W/m/K]"))
     96                                string="%s\n%s"%(string,fielddisplay(self,"meltingpoint","melting point of ice at 1atm in K"))
     97                                string="%s\n%s"%(string,fielddisplay(self,"latentheat","latent heat of fusion [J/m^3]"))
     98                                string="%s\n%s"%(string,fielddisplay(self,"beta","rate of change of melting point with pressure [K/Pa]"))
     99                                string="%s\n%s"%(string,fielddisplay(self,"mixed_layer_capacity","mixed layer capacity [W/kg/K]"))
     100                                string="%s\n%s"%(string,fielddisplay(self,"thermal_exchange_velocity","thermal exchange velocity [m/s]"))
     101                                string="%s\n%s"%(string,fielddisplay(self,"rheology_B","flow law parameter [Pa s^(1/n)]"))
     102                                string="%s\n%s"%(string,fielddisplay(self,"rheology_n","Glen's flow law exponent"))
     103                                string="%s\n%s"%(string,fielddisplay(self,"rheology_law","law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius' or 'LliboutryDuval'"))
     104                        elif nat=='litho':
     105                                string="%s\n%s"%(string,'Litho:');
     106                                string="%s\n%s"%(string,fielddisplay(self,'numlayers','number of layers (default 2)'))
     107                                string="%s\n%s"%(string,fielddisplay(self,'radius','array describing the radius for each interface (numlayers+1) [m]'))
     108                                string="%s\n%s"%(string,fielddisplay(self,'viscosity','array describing each layer''s viscosity (numlayers) [Pa.s]'))
     109                                string="%s\n%s"%(string,fielddisplay(self,'lame_lambda','array describing the lame lambda parameter (numlayers) [Pa]'))
     110                                string="%s\n%s"%(string,fielddisplay(self,'lame_mu','array describing the shear modulus for each layers (numlayers) [Pa]'))
     111                                string="%s\n%s"%(string,fielddisplay(self,'burgers_viscosity','array describing each layer''s transient viscosity, only for Burgers rheologies  (numlayers) [Pa.s]'))
     112                                string="%s\n%s"%(string,fielddisplay(self,'burgers_mu','array describing each layer''s transient shear modulus, only for Burgers rheologies  (numlayers) [Pa]'))
     113                                string="%s\n%s"%(string,fielddisplay(self,'isburgers','array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'))
     114                                string="%s\n%s"%(string,fielddisplay(self,'density','array describing each layer''s density (numlayers) [kg/m^3]'))
     115                                string="%s\n%s"%(string,fielddisplay(self,'issolid','array describing whether the layer is solid or liquid (default 1) (numlayers)'))
     116
     117                        else:
     118                                raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')");
    119119
    120120                return string
    121121                #}}}
     122
    122123        def extrude(self,md): # {{{
    123             for i in range(len(self.nature)):
    124                 nat=self.nature[i];
    125                 if nat=='ice':
    126                     self.rheology_B=project3d(md,'vector',self.rheology_B,'type','node')
    127                     self.rheology_n=project3d(md,'vector',self.rheology_n,'type','element')
    128             return self
     124                for i in range(len(self.nature)):
     125                        nat=self.nature[i];
     126                        if nat=='ice':
     127                                self.rheology_B=project3d(md,'vector',self.rheology_B,'type','node')
     128                                self.rheology_n=project3d(md,'vector',self.rheology_n,'type','element')
     129                        return self
    129130        #}}}
     131
    130132        def setdefaultparameters(self): # {{{
    131             for i in range(len(self.nature)):
    132                 nat=self.nature[i];
    133                 if nat=='ice':
    134                     #ice density (kg/m^3)
    135                     self.rho_ice=917.
    136 
    137                     #ocean water density (kg/m^3)
    138                     self.rho_water=1023.
    139 
    140                     #fresh water density (kg/m^3)
    141                     self.rho_freshwater=1000.
    142 
    143                     #water viscosity (N.s/m^2)
    144                     self.mu_water=0.001787 
    145 
    146                     #ice heat capacity cp (J/kg/K)
    147                     self.heatcapacity=2093.
    148 
    149                     #ice latent heat of fusion L (J/kg)
    150                     self.latentheat=3.34*10^5
    151 
    152                     #ice thermal conductivity (W/m/K)
    153                     self.thermalconductivity=2.4
    154                    
    155                     #wet ice thermal conductivity (W/m/K)
    156                     self.temperateiceconductivity=.24
    157 
    158                     #the melting point of ice at 1 atmosphere of pressure in K
    159                     self.meltingpoint=273.15
    160 
    161                     #rate of change of melting point with pressure (K/Pa)
    162                     self.beta=9.8*10^-8
    163 
    164                     #mixed layer (ice-water interface) heat capacity (J/kg/K)
    165                     self.mixed_layer_capacity=3974.
    166 
    167                     #thermal exchange velocity (ice-water interface) (m/s)
    168                     self.thermal_exchange_velocity=1.00*10^-4
    169 
    170                     #Rheology law: what is the temperature dependence of B with T
    171                     #available: none, paterson and arrhenius
    172                     self.rheology_law='Paterson'
    173 
    174                 elif nat=='litho':
    175                     #we default to a configuration that enables running GIA solutions using giacaron and/or giaivins.
    176                     self.numlayers=2
    177 
    178                     #center of the earth (approximation, must not be 0), then the lab (lithosphere/asthenosphere boundary) then the surface
    179                     #(with 1d3 to avoid numerical singularities)
    180                     self.radius=[1e3,6278*1e3,6378*1e3]
    181 
    182                     self.viscosity=[1e21,1e40] #mantle and lithosphere viscosity (respectively) [Pa.s]
    183                     self.lame_mu=[1.45*1e11,6.7*1e10]  # (Pa) #lithosphere and mantle shear modulus (respectively) [Pa]
    184                     self.lame_lambda=self.lame_mu  # (Pa) #mantle and lithosphere lamba parameter (respectively) [Pa]
    185                     self.burgers_viscosity=[np.nan,np.nan]
    186                     self.burgers_mu=[np.nan,np.nan]
    187                     self.isburgers=[0,0]
    188                     self.density=[5.51*1e3,5.50*1e3]  # (Pa) #mantle and lithosphere density [kg/m^3]
    189                     self.issolid=[1,1] # is layer solid or liquid.
    190 
    191                 else:
    192                     raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho')");
     133                for i in range(len(self.nature)):
     134                        nat=self.nature[i];
     135                        if nat=='ice':
     136                                #ice density (kg/m^3)
     137                                self.rho_ice=917.
     138                                #ocean water density (kg/m^3)
     139                                self.rho_water=1023.
     140                                #fresh water density (kg/m^3)
     141                                self.rho_freshwater=1000.
     142                                #water viscosity (N.s/m^2)
     143                                self.mu_water=0.001787
     144                                #ice heat capacity cp (J/kg/K)
     145                                self.heatcapacity=2093.
     146                                #ice latent heat of fusion L (J/kg)
     147                                self.latentheat=3.34*10^5
     148                                #ice thermal conductivity (W/m/K)
     149                                self.thermalconductivity=2.4
     150                                #wet ice thermal conductivity (W/m/K)
     151                                self.temperateiceconductivity=.24
     152                                #the melting point of ice at 1 atmosphere of pressure in K
     153                                self.meltingpoint=273.15
     154                                #rate of change of melting point with pressure (K/Pa)
     155                                self.beta=9.8*10^-8
     156                                #mixed layer (ice-water interface) heat capacity (J/kg/K)
     157                                self.mixed_layer_capacity=3974.
     158                                #thermal exchange velocity (ice-water interface) (m/s)
     159                                self.thermal_exchange_velocity=1.00*10^-4
     160                                #Rheology law: what is the temperature dependence of B with T
     161                                #available: none, paterson and arrhenius
     162                                self.rheology_law='Paterson'
     163
     164                        elif nat=='litho':
     165                                #we default to a configuration that enables running GIA solutions using giacaron and/or giaivins.
     166                                self.numlayers=2
     167                                #center of the earth (approximation, must not be 0), then the lab (lithosphere/asthenosphere boundary) then the surface
     168                                #(with 1d3 to avoid numerical singularities)
     169                                self.radius=[1e3,6278*1e3,6378*1e3]
     170                                self.viscosity=[1e21,1e40] #mantle and lithosphere viscosity (respectively) [Pa.s]
     171                                self.lame_mu=[1.45*1e11,6.7*1e10]  # (Pa) #lithosphere and mantle shear modulus (respectively) [Pa]
     172                                self.lame_lambda=self.lame_mu  # (Pa) #mantle and lithosphere lamba parameter (respectively) [Pa]
     173                                self.burgers_viscosity=[np.nan,np.nan]
     174                                self.burgers_mu=[np.nan,np.nan]
     175                                self.isburgers=[0,0]
     176                                self.density=[5.51*1e3,5.50*1e3]  # (Pa) #mantle and lithosphere density [kg/m^3]
     177                                self.issolid=[1,1] # is layer solid or liquid.
     178
     179                        else:
     180                                raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho')");
    193181
    194182                return self
    195183                #}}}
    196184        def checkconsistency(self,md,solution,analyses):    # {{{
    197             for i in range(len(self.nature)):
    198                 nat=self.nature[i];
    199                 if nat=='ice':
    200                     md = checkfield(md,'fieldname','materials.rho_ice','>',0)
    201                     md = checkfield(md,'fieldname','materials.rho_water','>',0)
    202                     md = checkfield(md,'fieldname','materials.rho_freshwater','>',0)
    203                     md = checkfield(md,'fieldname','materials.mu_water','>',0)
    204                     md = checkfield(md,'fieldname','materials.rheology_B','>',0,'timeseries',1,'NaN',1,'Inf',1)
    205                     md = checkfield(md,'fieldname','materials.rheology_n','>',0,'size',[md.mesh.numberofelements])
    206                     md = checkfield(md,'fieldname','materials.rheology_law','values',['None','BuddJacka','Cuffey','CuffeyTemperate','Paterson','Arrhenius','LliboutryDuval'])
    207                 elif nat=='litho':
    208                     if 'LoveAnalysis' not in analyses:
    209                         return md
    210 
    211                     md = checkfield(md,'fieldname','materials.numlayers','NaN',1,'Inf',1,'>',0,'numel',[1])
    212                     md = checkfield(md,'fieldname','materials.radius','NaN',1,'Inf',1,'size',[md.materials.numlayers+1,1],'>',0)
    213                     md = checkfield(md,'fieldname','materials.lame_mu','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0)
    214                     md = checkfield(md,'fieldname','materials.lame_lambda','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0)
    215                     md = checkfield(md,'fieldname','materials.issolid','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0,'<',2)
    216                     md = checkfield(md,'fieldname','materials.density','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>',0)
    217                     md = checkfield(md,'fieldname','materials.viscosity','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0)
    218                     md = checkfield(md,'fieldname','materials.isburgers','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0,'<=',1)
    219                     md = checkfield(md,'fieldname','materials.burgers_viscosity','Inf',1,'size',[md.materials.numlayers,1],'>=',0)
    220                     md = checkfield(md,'fieldname','materials.burgers_mu','Inf',1,'size',[md.materials.numlayers,1],'>=',0)
    221 
    222                     for i in range(md.materials.numlayers):
    223                         if md.materials.isburgers[i] and (np.isnan(md.materials.burgers_viscosity[i] or np.isnan(md.materials.burgers_mu[i]))):
    224                             raise RuntimeError("materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with isburgers choice")
    225                        
    226                     if md.materials.issolid[0]==0 or md.materials.lame_mu[0]==0:
    227                         raise RuntimeError('First layer must be solid (issolid(1) > 0 AND lame_mu(1) > 0). Add a weak inner core if necessary.')
    228                    
    229                     for i in range(md.materials.numlayers-1):
    230                         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
    231                             raise RuntimeError("%s%i%s"%('2 or more adjacent fluid layers detected starting at layer ',i,'. This is not supported yet. Consider merging them.'))
    232 
    233                 else:
    234                     raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho')");
     185                for i in range(len(self.nature)):
     186                        nat=self.nature[i];
     187                        if nat=='ice':
     188                                md = checkfield(md,'fieldname','materials.rho_ice','>',0)
     189                                md = checkfield(md,'fieldname','materials.rho_water','>',0)
     190                                md = checkfield(md,'fieldname','materials.rho_freshwater','>',0)
     191                                md = checkfield(md,'fieldname','materials.mu_water','>',0)
     192                                md = checkfield(md,'fieldname','materials.rheology_B','>',0,'timeseries',1,'NaN',1,'Inf',1)
     193                                md = checkfield(md,'fieldname','materials.rheology_n','>',0,'size',[md.mesh.numberofelements])
     194                                md = checkfield(md,'fieldname','materials.rheology_law','values',['None','BuddJacka','Cuffey','CuffeyTemperate','Paterson','Arrhenius','LliboutryDuval'])
     195
     196                        elif nat=='litho':
     197                                if 'LoveAnalysis' not in analyses:
     198                                        return md
     199
     200                                md = checkfield(md,'fieldname','materials.numlayers','NaN',1,'Inf',1,'>',0,'numel',[1])
     201                                md = checkfield(md,'fieldname','materials.radius','NaN',1,'Inf',1,'size',[md.materials.numlayers+1,1],'>',0)
     202                                md = checkfield(md,'fieldname','materials.lame_mu','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0)
     203                                md = checkfield(md,'fieldname','materials.lame_lambda','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0)
     204                                md = checkfield(md,'fieldname','materials.issolid','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0,'<',2)
     205                                md = checkfield(md,'fieldname','materials.density','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>',0)
     206                                md = checkfield(md,'fieldname','materials.viscosity','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0)
     207                                md = checkfield(md,'fieldname','materials.isburgers','NaN',1,'Inf',1,'size',[md.materials.numlayers,1],'>=',0,'<=',1)
     208                                md = checkfield(md,'fieldname','materials.burgers_viscosity','Inf',1,'size',[md.materials.numlayers,1],'>=',0)
     209                                md = checkfield(md,'fieldname','materials.burgers_mu','Inf',1,'size',[md.materials.numlayers,1],'>=',0)
     210
     211                                for i in range(md.materials.numlayers):
     212                                        if md.materials.isburgers[i] and (np.isnan(md.materials.burgers_viscosity[i] or np.isnan(md.materials.burgers_mu[i]))):
     213                                                raise RuntimeError("materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with isburgers choice")
     214
     215                                        if md.materials.issolid[0]==0 or md.materials.lame_mu[0]==0:
     216                                                raise RuntimeError('First layer must be solid (issolid(1) > 0 AND lame_mu(1) > 0). Add a weak inner core if necessary.')
     217
     218                                        for i in range(md.materials.numlayers-1):
     219                                                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
     220                                                        raise RuntimeError("%s%i%s"%('2 or more adjacent fluid layers detected starting at layer ',i,'. This is not supported yet. Consider merging them.'))
     221
     222                                                else:
     223                                                        raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho')");
    235224
    236225                return md
    237226        # }}}
     227
    238228        def marshall(self,prefix,md,fid):    # {{{
    239 
    240             #1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
    241             WriteData(fid,prefix,'name','md.materials.type','data',6,'format','Integer')
    242             WriteData(fid,prefix,'name','md.materials.nature','data',naturetointeger(self.nature),'format','IntMat','mattype',3)
    243 
    244             for i in range(len(self.nature)):
    245                 nat=self.nature[i];
    246                 if nat=='ice':
    247 
    248                     WriteData(fid,prefix,'name','md.materials.type','data',3,'format','Integer')
    249                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double')
    250                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_water','format','Double')
    251                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_freshwater','format','Double')
    252                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','mu_water','format','Double')
    253                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','heatcapacity','format','Double')
    254                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','latentheat','format','Double')
    255                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermalconductivity','format','Double')
    256                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','temperateiceconductivity','format','Double')
    257                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','meltingpoint','format','Double')
    258                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','beta','format','Double')
    259                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','mixed_layer_capacity','format','Double')
    260                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermal_exchange_velocity','format','Double')
    261                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts)
    262                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2)
    263                     WriteData(fid,prefix,'data',self.rheology_law,'name','md.materials.rheology_law','format','String')
    264 
    265                 elif nat=='litho':
    266                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','numlayers','format','Integer')
    267                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','radius','format','DoubleMat','mattype',3)
    268                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_mu','format','DoubleMat','mattype',3)
    269                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_lambda','format','DoubleMat','mattype',3)
    270                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','issolid','format','DoubleMat','mattype',3)
    271                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3)
    272                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3)
    273                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3)
    274                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3)
    275                     WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3)
    276 
    277                 else:
    278                     raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')")
    279 
     229                #1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
     230                WriteData(fid,prefix,'name','md.materials.type','data',6,'format','Integer')
     231                WriteData(fid,prefix,'name','md.materials.nature','data',naturetointeger(self.nature),'format','IntMat','mattype',3)
     232
     233                for i in range(len(self.nature)):
     234                        nat=self.nature[i];
     235                        if nat=='ice':
     236                                WriteData(fid,prefix,'name','md.materials.type','data',3,'format','Integer')
     237                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double')
     238                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_water','format','Double')
     239                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_freshwater','format','Double')
     240                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','mu_water','format','Double')
     241                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','heatcapacity','format','Double')
     242                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','latentheat','format','Double')
     243                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermalconductivity','format','Double')
     244                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','temperateiceconductivity','format','Double')
     245                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','meltingpoint','format','Double')
     246                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','beta','format','Double')
     247                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','mixed_layer_capacity','format','Double')
     248                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermal_exchange_velocity','format','Double')
     249                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts)
     250                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2)
     251                                WriteData(fid,prefix,'data',self.rheology_law,'name','md.materials.rheology_law','format','String')
     252
     253                        elif nat=='litho':
     254                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','numlayers','format','Integer')
     255                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','radius','format','DoubleMat','mattype',3)
     256                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_mu','format','DoubleMat','mattype',3)
     257                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_lambda','format','DoubleMat','mattype',3)
     258                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','issolid','format','DoubleMat','mattype',3)
     259                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3)
     260                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3)
     261                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3)
     262                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3)
     263                                WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3)
     264
     265                        else:
     266                                raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')")
    280267        # }}}

  • issm/trunk-jpl/src/py3/classes/matestar.py

    r23709 r23711  
    2222                thermalconductivity       = 0.
    2323                temperateiceconductivity  = 0.
     24    self.effectiveconductivity_averaging = 0.
    2425                meltingpoint              = 0.
    2526                beta                      = 0.
    … …  
    5354                string="%s\n%s"%(string,fielddisplay(self,'thermalconductivity',['ice thermal conductivity [W/m/K]']))
    5455                string="%s\n%s"%(string,fielddisplay(self,'temperateiceconductivity','temperate ice thermal conductivity [W/m/K]'))
     56                string="%s\n%s"%(string,fielddisplay(self,"effectiveconductivity_averaging","computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)"))
    5557                string="%s\n%s"%(string,fielddisplay(self,'meltingpoint','melting point of ice at 1atm in K'))
    5658                string="%s\n%s"%(string,fielddisplay(self,'latentheat','latent heat of fusion [J/kg]'))
    … …  
    9496                #wet ice thermal conductivity (W/m/K)
    9597                self.temperateiceconductivity=.24
     98    #computation of effective conductivity
     99    self.effectiveconductivity_averaging=1
    96100                #the melting point of ice at 1 atmosphere of pressure in K
    97101                self.meltingpoint=273.15
    … …  
    125129                md = checkfield(md,'fieldname','materials.rheology_Es','>',0,'size',[md.mesh.numberofvertices],'NaN',1,'Inf',1)
    126130                md = checkfield(md,'fieldname','materials.rheology_law','values',['None','BuddJacka', 'Cuffey','CuffeyTemperate','Paterson','Arrhenius','LliboutryDuval'])
     131                md = checkfield(md,'fieldname','materials.effectiveconductivity_averaging','numel',[1],'values',[0,1,2])
    127132
    128133                if 'GiaAnalysis' in analyses:
    … …  
    148153                WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermalconductivity','format','Double')
    149154                WriteData(fid,prefix,'object',self,'class','materials','fieldname','temperateiceconductivity','format','Double')
     155                WriteData(fid,prefix,'object',self,'class','materials','fieldname','effectiveconductivity_averaging','format','Integer')
    150156                WriteData(fid,prefix,'object',self,'class','materials','fieldname','meltingpoint','format','Double')
    151157                WriteData(fid,prefix,'object',self,'class','materials','fieldname','beta','format','Double')

  • issm/trunk-jpl/src/py3/classes/matice.py

    r23670 r23711  
    2121                self.thermalconductivity       = 0.
    2222                self.temperateiceconductivity  = 0.
     23                self.effectiveconductivity_averaging = 0.
    2324                self.meltingpoint              = 0.
    2425                self.beta                      = 0.
    … …  
    2930                self.rheology_law              = ''
    3031
    31                 #giaivins: 
     32                #giaivins:
    3233                self.lithosphere_shear_modulus  = 0.
    3334                self.lithosphere_density        = 0.
    3435                self.mantle_shear_modulus       = 0.
    35                 self.mantle_density             = 0. 
    36                
     36                self.mantle_density             = 0.
     37
    3738                #SLR
    38                 self.earth_density= 5512; 
    39 
    40 
     39                self.earth_density= 5512;
    4140
    4241                self.setdefaultparameters()
    4342                #}}}
     43
    4444        def __repr__(self): # {{{
    4545                string="   Materials:"
    … …  
    5252                string="%s\n%s"%(string,fielddisplay(self,"thermalconductivity","ice thermal conductivity [W/m/K]"))
    5353                string="%s\n%s"%(string,fielddisplay(self,"temperateiceconductivity","temperate ice thermal conductivity [W/m/K]"))
     54                string="%s\n%s"%(string,fielddisplay(self,"effectiveconductivity_averaging","computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)"))
    5455                string="%s\n%s"%(string,fielddisplay(self,"meltingpoint","melting point of ice at 1atm in K"))
    5556                string="%s\n%s"%(string,fielddisplay(self,"latentheat","latent heat of fusion [J/m^3]"))
    … …  
    6566                string="%s\n%s"%(string,fielddisplay(self,"mantle_density","Mantle density [g/cm^-3]"))
    6667                string="%s\n%s"%(string,fielddisplay(self,"earth_density","Mantle density [kg/m^-3]"))
    67 
    68 
    6968                return string
    7069                #}}}
     70
    7171        def extrude(self,md): # {{{
    7272                self.rheology_B=project3d(md,'vector',self.rheology_B,'type','node')
    … …  
    7474                return self
    7575        #}}}
     76
    7677        def setdefaultparameters(self): # {{{
    7778                #ice density (kg/m^3)
    7879                self.rho_ice=917.
    79 
    8080                #ocean water density (kg/m^3)
    8181                self.rho_water=1023.
    82 
    8382                #fresh water density (kg/m^3)
    8483                self.rho_freshwater=1000.
    85 
    8684                #water viscosity (N.s/m^2)
    87                 self.mu_water=0.001787 
    88 
     85                self.mu_water=0.001787
    8986                #ice heat capacity cp (J/kg/K)
    9087                self.heatcapacity=2093.
    91 
    9288                #ice latent heat of fusion L (J/kg)
    9389                self.latentheat=3.34*10**5
    94 
    9590                #ice thermal conductivity (W/m/K)
    9691                self.thermalconductivity=2.4
    97 
     92    #computation of effective conductivity
     93                self.effectiveconductivity_averaging=1
    9894                #temperate ice thermal conductivity (W/m/K)
    9995                self.temperateiceconductivity=0.24
    100 
    10196                #the melting point of ice at 1 atmosphere of pressure in K
    10297                self.meltingpoint=273.15
    103 
    10498                #rate of change of melting point with pressure (K/Pa)
    10599                self.beta=9.8*10**-8
    106 
    107100                #mixed layer (ice-water interface) heat capacity (J/kg/K)
    108101                self.mixed_layer_capacity=3974.
    109 
    110102                #thermal exchange velocity (ice-water interface) (m/s)
    111103                self.thermal_exchange_velocity=1.00*10**-4
    112 
    113104                #Rheology law: what is the temperature dependence of B with T
    114105                #available: none, paterson and arrhenius
    … …  
    120111                self.mantle_shear_modulus       = 1.45*10**11 # (Pa)
    121112                self.mantle_density             = 3.34        # (g/cm^-3)
    122                
     113
    123114                #SLR
    124115                self.earth_density= 5512;  # average density of the Earth, (kg/m^3)
    125 
    126 
    127116                return self
    128117                #}}}
     118
    129119        def checkconsistency(self,md,solution,analyses):    # {{{
    130120                md = checkfield(md,'fieldname','materials.rho_ice','>',0)
    … …  
    135125                md = checkfield(md,'fieldname','materials.rheology_n','>',0,'size',[md.mesh.numberofelements])
    136126                md = checkfield(md,'fieldname','materials.rheology_law','values',['None','BuddJacka','Cuffey','CuffeyTemperate','Paterson','Arrhenius','LliboutryDuval'])
     127                md = checkfield(md,'fieldname','materials.effectiveconductivity_averaging','numel',[1],'values',[0,1,2])
    137128                md = checkfield(md,'fieldname','materials.lithosphere_shear_modulus','>',0,'numel',[1]);
    138129                md = checkfield(md,'fieldname','materials.lithosphere_density','>',0,'numel',[1]);
    … …  
    140131                md = checkfield(md,'fieldname','materials.mantle_density','>',0,'numel',[1]);
    141132                md = checkfield(md,'fieldname','materials.earth_density','>',0,'numel',[1]);
    142 
    143133                return md
    144134        # }}}
     135
    145136        def marshall(self,prefix,md,fid):    # {{{
    146137                WriteData(fid,prefix,'name','md.materials.type','data',3,'format','Integer');
    … …  
    153144                WriteData(fid,prefix,'object',self,'class','materials','fieldname','thermalconductivity','format','Double')
    154145                WriteData(fid,prefix,'object',self,'class','materials','fieldname','temperateiceconductivity','format','Double')
     146                WriteData(fid,prefix,'object',self,'class','materials','fieldname','effectiveconductivity_averaging','format','Integer')
    155147                WriteData(fid,prefix,'object',self,'class','materials','fieldname','meltingpoint','format','Double')
    156148                WriteData(fid,prefix,'object',self,'class','materials','fieldname','beta','format','Double')

  • issm/trunk-jpl/src/py3/qmu/rlev_write.py

    r23709 r23711  
    44from vector_write import *
    55from param_write import *
     6from response_function import *
    67#import relevent qmu classes
    78from MatlabArray import *
    … …  
    3839  function to write response levels
    3940'''
    40         from response_function import *
     41#       from response_function import *
    4142
    4243        if len(dresp) == 0 or len(fieldnames(dresp[0])) == 0:

  • issm/trunk-jpl/src/py3/solve/parseresultsfromdisk.py

    r23691 r23711  
    146146                length=struct.unpack('i',fid.read(struct.calcsize('i')))[0]
    147147                fieldname=struct.unpack('{}s'.format(length),fid.read(length))[0][:-1]
     148                fieldname=fieldname.decode() #strings are booleans when stored so need to be converted back
    148149                time=struct.unpack('d',fid.read(struct.calcsize('d')))[0]
    149150                step=struct.unpack('i',fid.read(struct.calcsize('i')))[0]
    … …  
    246247
    247248                saveres=OrderedDict()
    248                 saveres['fieldname']=fieldname.decode()
     249                saveres['fieldname']=fieldname
    249250                saveres['time']=time
    250251                saveres['step']=step
Note: See TracChangeset for help on using the changeset viewer.