![(please configure the [header_logo] section in trac.ini)](/trac/issm/chrome/common/trac_banner.png)
Changeset 26233
- Timestamp:
- 05/04/21 10:00:04 (4 years ago)
- Location:
- issm/trunk-jpl/src
- Files:
-
- 4 edited
-
c/shared/Enum/EnumDefinitions.h (modified) (2 diffs)
-
m/classes/fourierlove.m (modified) (7 diffs)
-
m/classes/lovenumbers.m (modified) (2 diffs)
-
m/classes/materials.m (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
r26222 r26233 250 250 LockFileNameEnum, 251 251 LoveAllowLayerDeletionEnum, 252 LoveCoreMantleBoundaryEnum, 253 LoveEarthMassEnum, 252 254 LoveForcingTypeEnum, 253 255 LoveFrequenciesEnum, 256 LoveIsTemporalEnum, 254 257 LoveG0Enum, 258 LoveGravitationalConstantEnum, 259 LoveInnerCoreBoundaryEnum, 260 LoveIntStepsPerLayerEnum, 255 261 LoveKernelsEnum, 256 262 LoveMu0Enum, 257 263 LoveNfreqEnum, 264 LoveNTemporalIterationsEnum, 265 LoveNYiEquationsEnum, 258 266 LoveR0Enum, 259 267 LoveShNmaxEnum, 260 268 LoveShNminEnum, 269 LoveStartingLayerEnum, 270 LoveUnderflowTolEnum, 261 271 MassFluxSegmentsEnum, 262 272 MassFluxSegmentsPresentEnum, … … 779 789 SealevelGrotm2Enum, 780 790 SealevelGrotm3Enum, 791 SealevelGUrotm1Enum, 792 SealevelGUrotm2Enum, 793 SealevelGUrotm3Enum, 794 SealevelGNrotm1Enum, 795 SealevelGNrotm2Enum, 796 SealevelGNrotm3Enum, 797 SealevelGErotm1Enum, 798 SealevelGErotm2Enum, 799 SealevelGErotm3Enum, 781 800 SealevelRSLBarystaticEnum, 782 801 SealevelRSLRateEnum, -
issm/trunk-jpl/src/m/classes/fourierlove.m
r26059 r26233 6 6 classdef fourierlove 7 7 properties (SetAccess=public) 8 nfreq = 0; 9 frequencies = 0; 10 sh_nmax = 0; 11 sh_nmin = 0; 12 g0 = 0; 13 r0 = 0; 14 mu0 = 0; 15 allow_layer_deletion = 0; 16 love_kernels = 0; 17 forcing_type = 0; 8 nfreq = 0; 9 frequencies = 0; 10 sh_nmax = 0; 11 sh_nmin = 0; 12 g0 = 0; 13 r0 = 0; 14 mu0 = 0; 15 Gravitational_Constant = 0; 16 allow_layer_deletion = 0; 17 underflow_tol = 0; 18 integration_steps_per_layer= 0; 19 istemporal = 0; 20 n_temporal_iterations = 0; 21 time = 0; 22 love_kernels = 0; 23 forcing_type = 0; 24 inner_core_boundary = 0; 25 core_mantle_boundary = 0; 18 26 end 19 27 methods (Static) … … 44 52 self.r0=6371*1e3; %m; 45 53 self.mu0=10^11; % Pa 46 self.allow_layer_deletion=1; 54 self.Gravitational_Constant=6.67259e-11; % m^3 kg^-1 s^-2 55 self.allow_layer_deletion=1; 56 self.underflow_tol=1e-16; %threshold of deep to surface love number ratio to trigger the deletion of layer 57 self.integration_steps_per_layer=100; 58 self.istemporal=0; 59 self.n_temporal_iterations=8; 60 self.time=[0]; %s 47 61 self.love_kernels=0; 48 self.forcing_type = 11; 62 self.forcing_type = 11; % surface loading 63 self.inner_core_boundary=1; 64 self.core_mantle_boundary=2; 49 65 end % }}} 50 66 function disp(self) % {{{ … … 54 70 fielddisplay(self,'sh_nmin','minimum spherical harmonic degree (default 1)'); 55 71 fielddisplay(self,'g0','adimensioning constant for gravity (default 10) [m/s^2]'); 56 fielddisplay(self,'r0','adimensioning constant for radius (default 637 8*10^3) [m]');72 fielddisplay(self,'r0','adimensioning constant for radius (default 6371*10^3) [m]'); 57 73 fielddisplay(self,'mu0','adimensioning constant for stress (default 10^11) [Pa]'); 74 fielddisplay(self,'Gravitational_Constant','Newtonian constant of gravitation (default 6.67259e-11 [m^3 kg^-1 s^-2])'); 58 75 fielddisplay(self,'allow_layer_deletion','allow for migration of the integration boundary with increasing spherical harmonics degree (default 1)'); 76 fielddisplay(self,'underflow_tol','threshold of deep to surface love number ratio to trigger the deletion of layers (default 2.2204460492503131E-016)'); 77 fielddisplay(self,'integration_steps_per_layer','number of radial steps to propagate the yi system from the bottom to the top of each layer (default 100)'); 78 fielddisplay(self,'istemporal',{'1 for time-dependent love numbers, 0 for frequency-dependent or elastic love numbers (default 0)', 'If 1: use fourierlove function build_frequencies_from_time to meet consistency'}); 79 fielddisplay(self,'n_temporal_iterations','max number of iterations in the inverse Laplace transform. Also the number of spectral samples per time step requested (default 8)'); 80 fielddisplay(self,'time','time vector for deformation if istemporal (default 0) [s]'); 59 81 fielddisplay(self,'love_kernels','compute love numbers at depth? (default 0)'); 60 fielddisplay(self,'forcing_type',{'integer indicating the nature and depth of the forcing for the Love number calculation (default 11) :','1: Inner core boundary -- Volumic Potential','2: Inner core boundary -- Pressure','3: Inner core boundary -- Loading','4: Inner core boundary -- Tangential traction','5: Core mantle boundary -- Volumic Potential','6: Core mantle boundary -- Pressure','7: Core mantle boundary -- Loading','8: Core mantle boundary -- Tangential traction','9: Surface -- Volumic Potential','10: Surface -- Pressure','11: Surface -- Loading','12: Surface -- Tangential traction '}); 82 fielddisplay(self,'forcing_type',{'integer indicating the nature and depth of the forcing for the Love number calculation (default 11) :','1: Inner core boundary -- Volumic Potential','2: Inner core boundary -- Pressure','3: Inner core boundary -- Loading','4: Inner core boundary -- Tangential traction','5: Core mantle boundary -- Volumic Potential','6: Core mantle boundary -- Pressure','7: Core mantle boundary -- Loading','8: Core mantle boundary -- Tangential traction','9: Surface -- Volumic Potential','10: Surface -- Pressure','11: Surface -- Loading','12: Surface -- Tangential traction '}); 83 fielddisplay(self,'inner_core_boundary','interface index in materials.radius locating forcing. Only used for forcing_type 1--4 (default 1)'); 84 fielddisplay(self,'core_mantle_boundary','interface index in materials.radius locating forcing. Only used for forcing_type 5--8 (default 2)'); 61 85 62 86 end % }}} … … 72 96 md = checkfield(md,'fieldname','love.r0','NaN',1,'Inf',1,'numel',1,'>',0); 73 97 md = checkfield(md,'fieldname','love.mu0','NaN',1,'Inf',1,'numel',1,'>',0); 98 md = checkfield(md,'fieldname','love.Gravitational_Constant','NaN',1,'Inf',1,'numel',1,'>',0); 74 99 md = checkfield(md,'fieldname','love.allow_layer_deletion','values',[0 1]); 100 md = checkfield(md,'fieldname','love.underflow_tol','NaN',1,'Inf',1,'numel',1,'>',0); 101 md = checkfield(md,'fieldname','love.integration_steps_per_layer','NaN',1,'Inf',1,'numel',1,'>',0); 75 102 md = checkfield(md,'fieldname','love.love_kernels','values',[0 1]); 76 103 md = checkfield(md,'fieldname','love.forcing_type','NaN',1,'Inf',1,'numel',1,'>',0, '<=', 12); 77 if md.love.sh_nmin<=1 & md.love.forcing_type==9 104 105 md = checkfield(md,'fieldname','love.istemporal','values',[0 1]); 106 if md.love.istemporal==1 107 md = checkfield(md,'fieldname','love.n_temporal_iterations','NaN',1,'Inf',1,'numel',1,'>',0); 108 md = checkfield(md,'fieldname','love.time','NaN',1,'Inf',1,'numel',md.love.nfreq/2/md.love.n_temporal_iterations); 109 end 110 if md.love.sh_nmin<=1 & (md.love.forcing_type==9 || md.love.forcing_type==5 || md.love.forcing_type==1) 78 111 error('Degree 1 not supported for Volumetric Potential forcing. Use sh_min>=2 for this kind of calculation.') 79 112 end … … 82 115 if ~isa(md.materials,'materials') | ~sum(strcmpi(md.materials.nature,'litho')) 83 116 error('Need a ''litho'' material to run a Fourier Love number analysis'); 117 end 118 119 mat=find(strcmpi(md.materials.nature,'litho')); 120 if md.love.forcing_type<=4 121 md = checkfield(md,'fieldname','love.inner_core_boundary','NaN',1,'Inf',1,'numel',1,'>',0, '<=', md.materials(mat).numlayers); 122 elseif md.love.forcing_type<=8 123 md = checkfield(md,'fieldname','love.core_mantle_boundary','NaN',1,'Inf',1,'numel',1,'>',0, '<=', md.materials(mat).numlayers); 84 124 end 85 125 … … 94 134 WriteData(fid,prefix,'object',self,'fieldname','r0','format','Double'); 95 135 WriteData(fid,prefix,'object',self,'fieldname','mu0','format','Double'); 136 WriteData(fid,prefix,'object',self,'fieldname','Gravitational_Constant','format','Double'); 96 137 WriteData(fid,prefix,'object',self,'fieldname','allow_layer_deletion','format','Boolean'); 138 WriteData(fid,prefix,'object',self,'fieldname','underflow_tol','format','Double'); 139 WriteData(fid,prefix,'object',self,'fieldname','integration_steps_per_layer','format','Integer'); 140 WriteData(fid,prefix,'object',self,'fieldname','istemporal','format','Boolean'); 141 WriteData(fid,prefix,'object',self,'fieldname','n_temporal_iterations','format','Integer'); 142 %note: no need to marshall the time vector, we have frequencies 97 143 WriteData(fid,prefix,'object',self,'fieldname','love_kernels','format','Boolean'); 98 144 WriteData(fid,prefix,'object',self,'fieldname','forcing_type','format','Integer'); 145 WriteData(fid,prefix,'object',self,'fieldname','inner_core_boundary','format','Integer'); 146 WriteData(fid,prefix,'object',self,'fieldname','core_mantle_boundary','format','Integer'); 99 147 100 148 end % }}} … … 102 150 error('not implemented yet!'); 103 151 end % }}} 152 function self=build_frequencies_from_time(self) % {{{ 153 if ~self.istemporal 154 error('cannot build frequencies for temporal love numbers if love.istemporal==0') 155 end 156 disp('Temporal love numbers: Overriding md.love.nfreq and md.love.frequencies'); 157 self.nfreq=length(self.time)*2*self.n_temporal_iterations; 158 self.frequencies=zeros(self.nfreq,1); 159 for i=1:length(self.time) 160 for j=1:2*self.n_temporal_iterations 161 self.frequencies((i-1)*2*self.n_temporal_iterations +j) = j*log(2)/self.time(i)/2/pi; 162 end 163 end 164 end % }}} 104 165 end 105 166 end -
issm/trunk-jpl/src/m/classes/lovenumbers.m
r26047 r26233 57 57 md = checkfield(md,'fieldname','solidearth.lovenumbers.th','NaN',1,'Inf',1); 58 58 md = checkfield(md,'fieldname','solidearth.lovenumbers.tk','NaN',1,'Inf',1); 59 md = checkfield(md,'fieldname','solidearth.lovenumbers.tl','NaN',1,'Inf',1); 59 60 md = checkfield(md,'fieldname','solidearth.lovenumbers.tk2secular','NaN',1,'Inf',1); 60 61 … … 78 79 fielddisplay(self,'th','tidal load Love number (deg 2)'); 79 80 fielddisplay(self,'tk','tidal load Love number (deg 2)'); 81 fielddisplay(self,'tl','tidal load Love number (deg 2)'); 80 82 fielddisplay(self,'tk2secular','secular fluid Love number'); 81 83 -
issm/trunk-jpl/src/m/classes/materials.m
r26059 r26233 53 53 self.addprop('burgers_viscosity'); 54 54 self.addprop('burgers_mu'); 55 self.addprop('isburgers'); 55 self.addprop('ebm_alpha'); 56 self.addprop('ebm_delta'); 57 self.addprop('ebm_taul'); 58 self.addprop('ebm_tauh'); 59 self.addprop('rheologymodel'); 56 60 self.addprop('density'); 57 61 self.addprop('issolid'); … … 136 140 self.burgers_viscosity=[NaN;NaN]; 137 141 self.burgers_mu=[NaN;NaN]; 138 self.isburgers=[0;0]; 142 143 self.ebm_alpha=[NaN;NaN]; 144 self.ebm_delta=[NaN;NaN]; 145 self.ebm_taul=[NaN;NaN]; 146 self.ebm_tauh=[NaN;NaN]; 147 self.rheologymodel=[0;0]; 139 148 self.density=[5.51*1e3;5.50*1e3]; % (Pa) %mantle and lithosphere density [kg/m^3] 140 149 self.issolid=[1;1]; % is layer solid or liquid. … … 191 200 fielddisplay(self,'burgers_viscosity','array describing each layer''s transient viscosity, only for Burgers rheologies (numlayers) [Pa.s]'); 192 201 fielddisplay(self,'burgers_mu','array describing each layer''s transient shear modulus, only for Burgers rheologies (numlayers) [Pa]'); 193 fielddisplay(self,'isburgers','array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'); 202 203 fielddisplay(self,'ebm_alpha','array describing each layer''s exponent parameter controlling the shape of shear modulus curve between taul and tauh, only for EBM rheology (numlayers)'); 204 fielddisplay(self,'ebm_delta','array describing each layer''s amplitude of the transient relaxation (ratio between elastic rigity to pre-maxwell relaxation rigity), only for EBM rheology (numlayers)'); 205 fielddisplay(self,'ebm_taul','array describing each layer''s starting period for transient relaxation, only for EBM rheology (numlayers) [s]'); 206 fielddisplay(self,'ebm_tauh','array describing each layer''s array describing each layer''s end period for transient relaxation, only for Burgers rheology (numlayers) [s]'); 207 208 209 fielddisplay(self,'rheologymodel','array describing whether we adopt a MaxWell (0), Burgers (1) or EBM (2) rheology (default 0)'); 194 210 fielddisplay(self,'density','array describing each layer''s density (numlayers) [kg/m^3]'); 195 211 fielddisplay(self,'issolid','array describing whether the layer is solid or liquid (default 1) (numlayers)'); … … 228 244 md = checkfield(md,'fieldname','materials.density','NaN',1,'Inf',1,'size',[md.materials.numlayers 1],'>',0); 229 245 md = checkfield(md,'fieldname','materials.viscosity','NaN',1,'Inf',1,'size',[md.materials.numlayers 1],'>=',0); 230 md = checkfield(md,'fieldname','materials. isburgers','NaN',1,'Inf',1,'size',[md.materials.numlayers 1],'>=',0,'<=',1);246 md = checkfield(md,'fieldname','materials.rheologymodel','NaN',1,'Inf',1,'size',[md.materials.numlayers 1],'>=',0,'<=',2); 231 247 md = checkfield(md,'fieldname','materials.burgers_viscosity','Inf',1,'size',[md.materials.numlayers 1],'>=',0); 232 248 md = checkfield(md,'fieldname','materials.burgers_mu','Inf',1,'size',[md.materials.numlayers 1],'>=',0); 249 md = checkfield(md,'fieldname','materials.ebm_alpha','Inf',1,'size',[md.materials.numlayers 1],'>=',0); 250 md = checkfield(md,'fieldname','materials.ebm_delta','Inf',1,'size',[md.materials.numlayers 1],'>=',0); 251 md = checkfield(md,'fieldname','materials.ebm_taul','Inf',1,'size',[md.materials.numlayers 1],'>=',0); 252 md = checkfield(md,'fieldname','materials.ebm_tauh','Inf',1,'size',[md.materials.numlayers 1],'>=',0); 233 253 234 254 for i=1:md.materials.numlayers, 235 if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))), 236 error('materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with isburgers choice'); 255 if md.materials.rheologymodel(i)==1 & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))), 256 error('materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with rheologymodel choice'); 257 end 258 if md.materials.rheologymodel(i)==2 & (isnan(md.materials.ebm_alpha(i)) | isnan(md.materials.ebm_delta(i)) | isnan(md.materials.ebm_taul(i)) | isnan(md.materials.ebm_tauh(i)) ), 259 error('materials checkconsistency error message: Litho ebm_alpha, ebm_delta, ebm_taul or ebm_tauh has NaN values, inconsistent with rheologymodel choice'); 237 260 end 238 261 end … … 290 313 WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3); 291 314 WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3); 292 WriteData(fid,prefix,'object',self,'class','materials','fieldname',' isburgers','format','DoubleMat','mattype',3);315 WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheologymodel','format','DoubleMat','mattype',3); 293 316 WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3); 294 317 WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3); 318 WriteData(fid,prefix,'object',self,'class','materials','fieldname','ebm_alpha','format','DoubleMat','mattype',3); 319 WriteData(fid,prefix,'object',self,'class','materials','fieldname','ebm_delta','format','DoubleMat','mattype',3); 320 WriteData(fid,prefix,'object',self,'class','materials','fieldname','ebm_taul','format','DoubleMat','mattype',3); 321 WriteData(fid,prefix,'object',self,'class','materials','fieldname','ebm_tauh','format','DoubleMat','mattype',3); 295 322 %compute earth density compatible with our layer density distribution: 296 323 earth_density=0; … … 351 378 writejsdouble(fid,[modelname '.materials.density'],self.density); 352 379 writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity); 353 writejsdouble(fid,[modelname '.materials. isburgers'],self.isburgers);380 writejsdouble(fid,[modelname '.materials.rheologymodel'],self.rheologymodel); 354 381 writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity); 355 382 writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu); 383 writejsdouble(fid,[modelname '.materials.ebm_alpha'],self.ebm_alpha); 384 writejsdouble(fid,[modelname '.materials.ebm_delta'],self.ebm_delta); 385 writejsdouble(fid,[modelname '.materials.ebm_taul'],self.ebm_taul); 386 writejsdouble(fid,[modelname '.materials.ebm_tauh'],self.ebm_tauh); 356 387 357 388 case 'hydro'
Note:
See TracChangeset
for help on using the changeset viewer.
