%SLR class definition % % Usage: % slr=slr(); classdef slr properties (SetAccess=public) deltathickness = NaN; sealevel = NaN; maxiter = 0; reltol = 0; abstol = 0; love_h = 0; %provided by PREM model love_k = 0; %ideam tide_love_k = 0; %ideam tide_love_h = 0; %ideam rigid = 0; elastic = 0; rotation = 0; degacc = 0; requested_outputs = {}; transitions = {}; end methods function self = slr(varargin) % {{{ switch nargin case 0 self=setdefaultparameters(self); otherwise error('constructor not supported'); end end % }}} function self = setdefaultparameters(self) % {{{ %Convergence criterion: absolute, relative and residual self.reltol=NaN; %default self.abstol=0.001; %1 mm of sea level rise %maximum of non-linear iterations. self.maxiter=10; %computational flags: self.rigid=1; self.elastic=1; self.rotation=1; %tidal love numbers: self.tide_love_h=0.6149; %degree 2 self.tide_love_k=0.3055; % degree 2 %numerical discretization accuracy self.degacc=.01; %output default: self.requested_outputs={'default'}; %transitions should be a cell array of vectors: self.transitions={}; end % }}} function md = checkconsistency(self,md,solution,analyses) % {{{ if ~ismember(SealevelriseAnalysisEnum(),analyses), return; end md = checkfield(md,'fieldname','slr.deltathickness','NaN',1,'Inf',1,'size',[md.mesh.numberofelements 1]); %md = checkfield(md,'fieldname','slr.deltathickness','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]); md = checkfield(md,'fieldname','slr.sealevel','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]); md = checkfield(md,'fieldname','slr.love_h','NaN',1,'Inf',1); md = checkfield(md,'fieldname','slr.love_k','NaN',1,'Inf',1); md = checkfield(md,'fieldname','slr.tide_love_h','NaN',1,'Inf',1); md = checkfield(md,'fieldname','slr.tide_love_k','NaN',1,'Inf',1); md = checkfield(md,'fieldname','slr.love_k','NaN',1,'Inf',1); md = checkfield(md,'fieldname','slr.reltol','size',[1 1]); md = checkfield(md,'fieldname','slr.abstol','size',[1 1]); md = checkfield(md,'fieldname','slr.maxiter','size',[1 1],'>=',1); md = checkfield(md,'fieldname','slr.degacc','size',[1 1],'>=',1e-10); md = checkfield(md,'fieldname','slr.requested_outputs','stringrow',1); %check that love numbers are provided at the same level of accuracy: if (size(self.love_h,1) ~= size(self.love_k,1)), error('slr error message: love numbers should be provided at the same level of accuracy'); end %cross check that whereever we have an ice load, the mask is <0 on each vertex: pos=find(self.deltathickness); maskpos=md.mask.ice_levelset(md.mesh.elements(pos,:)); [els,vertices]=find(maskpos>0); if length(els), error('slr checkconsistency fail: there are elements with ice loads where some vertices are not on the ice!'); end end % }}} function list=defaultoutputs(self,md) % {{{ list = {'Sealevel'}; end % }}} function disp(self) % {{{ disp(sprintf(' slr parameters:')); fielddisplay(self,'deltathickness','thickness change (main loading of the slr solution core [m]'); fielddisplay(self,'sealevel','current sea level (prior to computation) [m]'); fielddisplay(self,'reltol','sea level rise relative convergence criterion, (default, NaN: not applied)'); fielddisplay(self,'abstol','sea level rise absolute convergence criterion, NaN: not applied'); fielddisplay(self,'maxiter','maximum number of nonlinear iterations'); fielddisplay(self,'love_h','love load number for radial displacement'); fielddisplay(self,'love_k','love load number for gravitational potential perturbation'); fielddisplay(self,'tide_love_k','tidal love number (deg 2)'); fielddisplay(self,'tide_love_h','tidal love number (deg 2)'); fielddisplay(self,'rotation','earth rotational potential perturbation'); fielddisplay(self,'rigid','rigid earth graviational potential perturbation'); fielddisplay(self,'elastic','elastic earth graviational potential perturbation'); fielddisplay(self,'degacc','accuracy (default .01 deg) for numerical discretization of the Green''s functions'); fielddisplay(self,'transitions','indices into parts of the mesh that will be icecaps'); fielddisplay(self,'requested_outputs','additional outputs requested'); end % }}} function marshall(self,md,fid) % {{{ WriteData(fid,'object',self,'class','sealevelrise','fieldname','deltathickness','format','DoubleMat','mattype',2); %WriteData(fid,'object',self,'class','sealevelrise','fieldname','deltathickness','format','DoubleMat','mattype',1); WriteData(fid,'data',self.sealevel,'mattype',1,'format','DoubleMat','enum',SealevelEnum(),'timeserieslength',md.mesh.numberofvertices+1); WriteData(fid,'object',self,'class','sealevelrise','fieldname','reltol','format','Double'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','abstol','format','Double'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','maxiter','format','Integer'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','love_h','format','DoubleMat','mattype',1); WriteData(fid,'object',self,'class','sealevelrise','fieldname','love_k','format','DoubleMat','mattype',1); WriteData(fid,'object',self,'class','sealevelrise','fieldname','tide_love_k','format','Double'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','tide_love_h','format','Double'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','rigid','format','Boolean'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','elastic','format','Boolean'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','rotation','format','Boolean'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','degacc','format','Double'); WriteData(fid,'object',self,'class','sealevelrise','fieldname','transitions','format','MatArray'); %process requested outputs outputs = self.requested_outputs; pos = find(ismember(outputs,'default')); if ~isempty(pos), outputs(pos) = []; %remove 'default' from outputs outputs = [outputs defaultoutputs(self,md)]; %add defaults end WriteData(fid,'data',outputs,'enum',SealevelriseRequestedOutputsEnum,'format','StringArray'); end % }}} function savemodeljs(self,fid,modelname) % {{{ writejs1Darray(fid,[modelname '.slr.deltathickness'],self.deltathickness); writejs1Darray(fid,[modelname '.slr.sealevel'],self.sealevel); writejsdouble(fid,[modelname '.slr.maxiter'],self.maxiter); writejsdouble(fid,[modelname '.slr.reltol'],self.reltol); writejsdouble(fid,[modelname '.slr.abstol'],self.abstol); writejs1Darray(fid,[modelname '.slr.love_h'],self.love_h); writejs1Darray(fid,[modelname '.slr.love_k'],self.love_k); writejsdouble(fid,[modelname '.slr.tide_love_k'],self.tide_love_k); writejsdouble(fid,[modelname '.slr.tide_love_h'],self.tide_love_h); writejsdouble(fid,[modelname '.slr.rigid'],self.rigid); writejsdouble(fid,[modelname '.slr.rotation'],self.rotation); writejsdouble(fid,[modelname '.slr.elastic'],self.elastic); writejsdouble(fid,[modelname '.slr.degacc'],self.degacc); writejscellstring(fid,[modelname '.slr.requested_outputs'],self.requested_outputs); writejscellarray(fid,[modelname '.slr.transitions'],self.transitions); end % }}} end end