Changeset 22955 for issm/trunk-jpl/src/m/classes/slr.m
- Timestamp:
- 07/16/18 15:39:26 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/m/classes/slr.m
r22808 r22955 8 8 deltathickness = NaN; 9 9 sealevel = NaN; 10 spcthickness = NaN; 10 11 maxiter = 0; 11 12 reltol = 0; … … 25 26 ocean_area_scaling = 0; 26 27 steric_rate = 0; %rate of ocean expansion from steric effects. 28 geodetic_run_frequency = 1; %how many time steps we skip before we run the geodetic part of the solver during transient 29 geodetic = 0; %compute geodetic SLR? (in addition to steric?) 27 30 degacc = 0; 31 loop_increment = 0; 32 horiz = 0; 33 Ngia = 0; 34 Ugia = 0; 28 35 requested_outputs = {}; 29 36 transitions = {}; … … 46 53 %maximum of non-linear iterations. 47 54 self.maxiter=5; 55 self.loop_increment=200; 48 56 49 57 %computational flags: 58 self.geodetic=0; 50 59 self.rigid=1; 51 60 self.elastic=1; 52 self.rotation=0;53 61 self.ocean_area_scaling=0; 62 self.rotation=1; 54 63 55 64 %tidal love numbers: … … 73 82 self.steric_rate=0; 74 83 84 %how many time steps we skip before we run SLR solver during transient 85 self.geodetic_run_frequency=1; 75 86 76 87 %output default: … … 79 90 %transitions should be a cell array of vectors: 80 91 self.transitions={}; 92 93 %horizontal displacement? (not by default) 94 self.horiz=0; 81 95 82 96 end % }}} … … 86 100 md = checkfield(md,'fieldname','slr.deltathickness','NaN',1,'Inf',1,'size',[md.mesh.numberofelements 1]); 87 101 md = checkfield(md,'fieldname','slr.sealevel','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]); 102 md = checkfield(md,'fieldname','slr.spcthickness','Inf',1,'timeseries',1); 88 103 md = checkfield(md,'fieldname','slr.love_h','NaN',1,'Inf',1); 89 104 md = checkfield(md,'fieldname','slr.love_k','NaN',1,'Inf',1); … … 98 113 md = checkfield(md,'fieldname','slr.abstol','size',[1 1]); 99 114 md = checkfield(md,'fieldname','slr.maxiter','size',[1 1],'>=',1); 115 md = checkfield(md,'fieldname','slr.geodetic_run_frequency','size',[1 1],'>=',1); 100 116 md = checkfield(md,'fieldname','slr.steric_rate','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]); 101 117 md = checkfield(md,'fieldname','slr.degacc','size',[1 1],'>=',1e-10); 102 118 md = checkfield(md,'fieldname','slr.requested_outputs','stringrow',1); 119 md = checkfield(md,'fieldname','slr.loop_increment','NaN',1,'Inf',1,'>=',1); 120 md = checkfield(md,'fieldname','slr.horiz','NaN',1,'Inf',1,'values',[0 1]); 121 md = checkfield(md,'fieldname','slr.Ngia','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]); 122 md = checkfield(md,'fieldname','slr.Ugia','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]); 103 123 104 124 %check that love numbers are provided at the same level of accuracy: … … 112 132 [els,vertices]=find(maskpos>0); 113 133 if length(els), 114 error('slr checkconsistency fail: there are elements with ice loads where some vertices are not on the ice!'); 115 end 134 warning('slr checkconsistency fail: there are elements with ice loads where some vertices are not on the ice!'); 135 end 136 137 %check that if geodetic is requested, we are a mesh3dsurface model (planet), or if we are not, 138 %a coupler to a planet model is provided. 139 if self.geodetic, 140 if md.transient.iscoupler, 141 %we are good; 142 else 143 if strcmpi(class(md.mesh),'mesh3dsurface'), 144 %we are good 145 else 146 error('model is requesting geodetic computations without being a mesh3dsurface, or being coupled to one!'); 147 end 148 end 149 end 150 116 151 117 152 end % }}} … … 124 159 fielddisplay(self,'deltathickness','thickness change: ice height equivalent [m]'); 125 160 fielddisplay(self,'sealevel','current sea level (prior to computation) [m]'); 126 fielddisplay(self,'reltol','sea level rise relative convergence criterion, (NaN: not applied)'); 127 fielddisplay(self,'abstol','sea level rise absolute convergence criterion, (default, NaN: not applied)'); 161 fielddisplay(self,'spcthickness','thickness constraints (NaN means no constraint) [m]'); 162 fielddisplay(self,'reltol','sea level rise relative convergence criterion, (default, NaN: not applied)'); 163 fielddisplay(self,'abstol','sea level rise absolute convergence criterion, NaN: not applied'); 128 164 fielddisplay(self,'maxiter','maximum number of nonlinear iterations'); 129 165 fielddisplay(self,'love_h','load Love number for radial displacement'); … … 136 172 fielddisplay(self,'polar_moi','polar moment of inertia [kg m^2]'); 137 173 fielddisplay(self,'angular_velocity','mean rotational velocity of earth [per second]'); 174 fielddisplay(self,'ocean_area_scaling','correction for model representation of ocean area [default: No correction]'); 175 fielddisplay(self,'steric_rate','rate of steric ocean expansion (in mm/yr)'); 176 fielddisplay(self,'Ngia','rate of viscous (GIA) geoid expansion (in mm/yr)'); 177 fielddisplay(self,'Ugia','rate of viscous (GIA) bedrock uplift (in mm/yr)'); 178 fielddisplay(self,'loop_increment','vector assembly (in the convolution) framentation'); 179 fielddisplay(self,'geodetic','compute geodetic SLR? ( in addition to steric?) default 0'); 180 fielddisplay(self,'geodetic_run_frequency','how many time steps we skip before we run SLR solver during transient (default: 1)'); 138 181 fielddisplay(self,'rigid','rigid earth graviational potential perturbation'); 139 182 fielddisplay(self,'elastic','elastic earth graviational potential perturbation'); 140 183 fielddisplay(self,'rotation','earth rotational potential perturbation'); 141 fielddisplay(self,'ocean_area_scaling','correction for model representation of ocean area (default: No correction)');142 fielddisplay(self,'steric_rate','rate of steric ocean expansion [mm/yr]');143 184 fielddisplay(self,'degacc','accuracy (default .01 deg) for numerical discretization of the Green''s functions'); 144 185 fielddisplay(self,'transitions','indices into parts of the mesh that will be icecaps'); … … 149 190 %WriteData(fid,prefix,'object',self,'fieldname','deltathickness','format','DoubleMat','mattype',2); 150 191 WriteData(fid,prefix,'object',self,'fieldname','deltathickness','format','DoubleMat','mattype',2,'timeserieslength',md.mesh.numberofelements+1,'yts',md.constants.yts); 192 %WriteData(fid,prefix,'object',self,'fieldname','deltathickness','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofelements+1); 151 193 WriteData(fid,prefix,'object',self,'fieldname','sealevel','mattype',1,'format','DoubleMat','timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts); 194 WriteData(fid,prefix,'object',self,'fieldname','spcthickness','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts); 152 195 WriteData(fid,prefix,'object',self,'fieldname','reltol','format','Double'); 153 196 WriteData(fid,prefix,'object',self,'fieldname','abstol','format','Double'); … … 166 209 WriteData(fid,prefix,'object',self,'fieldname','rotation','format','Boolean'); 167 210 WriteData(fid,prefix,'object',self,'fieldname','ocean_area_scaling','format','Boolean'); 211 WriteData(fid,prefix,'object',self,'fieldname','geodetic_run_frequency','format','Integer'); 168 212 WriteData(fid,prefix,'object',self,'fieldname','steric_rate','format','DoubleMat','mattype',1,'scale',1e-3/md.constants.yts); 213 WriteData(fid,prefix,'object',self,'fieldname','Ngia','format','DoubleMat','mattype',1,'scale',1e-3/md.constants.yts); 214 WriteData(fid,prefix,'object',self,'fieldname','Ugia','format','DoubleMat','mattype',1,'scale',1e-3/md.constants.yts); 169 215 WriteData(fid,prefix,'object',self,'fieldname','degacc','format','Double'); 170 216 WriteData(fid,prefix,'object',self,'fieldname','transitions','format','MatArray'); 217 WriteData(fid,prefix,'object',self,'fieldname','loop_increment','format','Integer'); 218 WriteData(fid,prefix,'object',self,'fieldname','horiz','format','Integer'); 219 WriteData(fid,prefix,'object',self,'fieldname','geodetic','format','Integer'); 171 220 172 221 %process requested outputs … … 184 233 writejs1Darray(fid,[modelname '.slr.deltathickness'],self.deltathickness); 185 234 writejs1Darray(fid,[modelname '.slr.sealevel'],self.sealevel); 235 writejs1Darray(fid,[modelname '.slr.spcthickness'],self.spcthickness); 186 236 writejsdouble(fid,[modelname '.slr.maxiter'],self.maxiter); 187 237 writejsdouble(fid,[modelname '.slr.reltol'],self.reltol); … … 200 250 writejsdouble(fid,[modelname '.slr.rotation'],self.rotation); 201 251 writejsdouble(fid,[modelname '.slr.ocean_area_scaling'],self.ocean_area_scaling); 252 writejsdouble(fid,[modelname '.slr.geodetic_run_frequency'],self.geodetic_run_frequency); 202 253 writejs1Darray(fid,[modelname '.slr.steric_rate'],self.steric_rate); 203 254 writejsdouble(fid,[modelname '.slr.degacc'],self.degacc); … … 205 256 writejscellarray(fid,[modelname '.slr.transitions'],self.transitions); 206 257 end % }}} 258 function self = extrude(self,md) % {{{ 259 self.sealevel=project3d(md,'vector',self.sealevel,'type','node'); 260 end % }}} 207 261 end 208 262 end
Note:
See TracChangeset
for help on using the changeset viewer.