[19984] | 1 | %SLR class definition
|
---|
| 2 | %
|
---|
| 3 | % Usage:
|
---|
| 4 | % slr=slr();
|
---|
| 5 |
|
---|
| 6 | classdef slr
|
---|
| 7 | properties (SetAccess=public)
|
---|
| 8 | deltathickness = NaN;
|
---|
[20154] | 9 | sealevel = NaN;
|
---|
[19984] | 10 | maxiter = 0;
|
---|
| 11 | reltol = 0;
|
---|
| 12 | abstol = 0;
|
---|
| 13 | love_h = 0; %provided by PREM model
|
---|
| 14 | love_k = 0; %ideam
|
---|
[21000] | 15 | love_l = 0; %ideam
|
---|
[20349] | 16 | tide_love_k = 0; %ideam
|
---|
| 17 | tide_love_h = 0; %ideam
|
---|
| 18 | rigid = 0;
|
---|
| 19 | elastic = 0;
|
---|
| 20 | rotation = 0;
|
---|
[21295] | 21 | ocean_area_scaling = 0;
|
---|
[20033] | 22 | degacc = 0;
|
---|
[20036] | 23 | requested_outputs = {};
|
---|
[20137] | 24 | transitions = {};
|
---|
[19984] | 25 | end
|
---|
| 26 | methods
|
---|
| 27 | function self = slr(varargin) % {{{
|
---|
| 28 | switch nargin
|
---|
| 29 | case 0
|
---|
| 30 | self=setdefaultparameters(self);
|
---|
| 31 | otherwise
|
---|
| 32 | error('constructor not supported');
|
---|
| 33 | end
|
---|
| 34 | end % }}}
|
---|
| 35 | function self = setdefaultparameters(self) % {{{
|
---|
| 36 |
|
---|
| 37 | %Convergence criterion: absolute, relative and residual
|
---|
[20074] | 38 | self.reltol=NaN; %default
|
---|
[19984] | 39 | self.abstol=0.001; %1 mm of sea level rise
|
---|
| 40 |
|
---|
| 41 | %maximum of non-linear iterations.
|
---|
[20074] | 42 | self.maxiter=10;
|
---|
[19984] | 43 |
|
---|
| 44 | %computational flags:
|
---|
| 45 | self.rigid=1;
|
---|
| 46 | self.elastic=1;
|
---|
[20349] | 47 | self.rotation=1;
|
---|
[21309] | 48 | self.ocean_area_scaling=0;
|
---|
[20033] | 49 |
|
---|
[20349] | 50 | %tidal love numbers:
|
---|
| 51 | self.tide_love_h=0.6149; %degree 2
|
---|
| 52 | self.tide_love_k=0.3055; % degree 2
|
---|
| 53 |
|
---|
[20033] | 54 | %numerical discretization accuracy
|
---|
| 55 | self.degacc=.01;
|
---|
[19986] | 56 |
|
---|
[20036] | 57 | %output default:
|
---|
| 58 | self.requested_outputs={'default'};
|
---|
[20137] | 59 |
|
---|
| 60 | %transitions should be a cell array of vectors:
|
---|
| 61 | self.transitions={};
|
---|
[20036] | 62 |
|
---|
[19984] | 63 | end % }}}
|
---|
| 64 | function md = checkconsistency(self,md,solution,analyses) % {{{
|
---|
| 65 |
|
---|
[21049] | 66 | if ~ismember('SealevelriseAnalysis',analyses), return; end
|
---|
[20322] | 67 | md = checkfield(md,'fieldname','slr.deltathickness','NaN',1,'Inf',1,'size',[md.mesh.numberofelements 1]);
|
---|
[20154] | 68 | md = checkfield(md,'fieldname','slr.sealevel','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]);
|
---|
[19984] | 69 | md = checkfield(md,'fieldname','slr.love_h','NaN',1,'Inf',1);
|
---|
| 70 | md = checkfield(md,'fieldname','slr.love_k','NaN',1,'Inf',1);
|
---|
[21000] | 71 | md = checkfield(md,'fieldname','slr.love_l','NaN',1,'Inf',1);
|
---|
[20349] | 72 | md = checkfield(md,'fieldname','slr.tide_love_h','NaN',1,'Inf',1);
|
---|
| 73 | md = checkfield(md,'fieldname','slr.tide_love_k','NaN',1,'Inf',1);
|
---|
[19984] | 74 | md = checkfield(md,'fieldname','slr.reltol','size',[1 1]);
|
---|
| 75 | md = checkfield(md,'fieldname','slr.abstol','size',[1 1]);
|
---|
| 76 | md = checkfield(md,'fieldname','slr.maxiter','size',[1 1],'>=',1);
|
---|
[20033] | 77 | md = checkfield(md,'fieldname','slr.degacc','size',[1 1],'>=',1e-10);
|
---|
[20036] | 78 | md = checkfield(md,'fieldname','slr.requested_outputs','stringrow',1);
|
---|
[19984] | 79 |
|
---|
| 80 | %check that love numbers are provided at the same level of accuracy:
|
---|
[21000] | 81 | if (size(self.love_h,1)~=size(self.love_k,1) | size(self.love_h,1)~=size(self.love_l,1)),
|
---|
[19984] | 82 | error('slr error message: love numbers should be provided at the same level of accuracy');
|
---|
| 83 | end
|
---|
| 84 |
|
---|
[20322] | 85 | %cross check that whereever we have an ice load, the mask is <0 on each vertex:
|
---|
| 86 | pos=find(self.deltathickness);
|
---|
| 87 | maskpos=md.mask.ice_levelset(md.mesh.elements(pos,:));
|
---|
| 88 | [els,vertices]=find(maskpos>0);
|
---|
| 89 | if length(els),
|
---|
| 90 | error('slr checkconsistency fail: there are elements with ice loads where some vertices are not on the ice!');
|
---|
| 91 | end
|
---|
| 92 |
|
---|
[19984] | 93 | end % }}}
|
---|
[20036] | 94 | function list=defaultoutputs(self,md) % {{{
|
---|
[20154] | 95 | list = {'Sealevel'};
|
---|
[20036] | 96 | end % }}}
|
---|
[19984] | 97 | function disp(self) % {{{
|
---|
| 98 | disp(sprintf(' slr parameters:'));
|
---|
| 99 |
|
---|
| 100 | fielddisplay(self,'deltathickness','thickness change (main loading of the slr solution core [m]');
|
---|
[20154] | 101 | fielddisplay(self,'sealevel','current sea level (prior to computation) [m]');
|
---|
[20074] | 102 | fielddisplay(self,'reltol','sea level rise relative convergence criterion, (default, NaN: not applied)');
|
---|
[19984] | 103 | fielddisplay(self,'abstol','sea level rise absolute convergence criterion, NaN: not applied');
|
---|
| 104 | fielddisplay(self,'maxiter','maximum number of nonlinear iterations');
|
---|
[21000] | 105 | fielddisplay(self,'love_h','load Love number for radial displacement');
|
---|
| 106 | fielddisplay(self,'love_k','load Love number for gravitational potential perturbation');
|
---|
| 107 | fielddisplay(self,'love_l','load Love number for horizontal displacements');
|
---|
| 108 | fielddisplay(self,'tide_love_k','tidal load Love number (deg 2)');
|
---|
| 109 | fielddisplay(self,'tide_love_h','tidal load Love number (deg 2)');
|
---|
[19984] | 110 | fielddisplay(self,'rigid','rigid earth graviational potential perturbation');
|
---|
| 111 | fielddisplay(self,'elastic','elastic earth graviational potential perturbation');
|
---|
[21295] | 112 | fielddisplay(self,'rotation','earth rotational potential perturbation');
|
---|
[21309] | 113 | fielddisplay(self,'ocean_area_scaling','correction for model representation of ocean area [default: No correction]');
|
---|
[20033] | 114 | fielddisplay(self,'degacc','accuracy (default .01 deg) for numerical discretization of the Green''s functions');
|
---|
[20137] | 115 | fielddisplay(self,'transitions','indices into parts of the mesh that will be icecaps');
|
---|
[20036] | 116 | fielddisplay(self,'requested_outputs','additional outputs requested');
|
---|
[19984] | 117 |
|
---|
| 118 | end % }}}
|
---|
[20690] | 119 | function marshall(self,prefix,md,fid) % {{{
|
---|
| 120 | WriteData(fid,prefix,'object',self,'fieldname','deltathickness','format','DoubleMat','mattype',2);
|
---|
[20902] | 121 | WriteData(fid,prefix,'object',self,'fieldname','sealevel','mattype',1,'format','DoubleMat','timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
|
---|
[20690] | 122 | WriteData(fid,prefix,'object',self,'fieldname','reltol','format','Double');
|
---|
| 123 | WriteData(fid,prefix,'object',self,'fieldname','abstol','format','Double');
|
---|
| 124 | WriteData(fid,prefix,'object',self,'fieldname','maxiter','format','Integer');
|
---|
| 125 | WriteData(fid,prefix,'object',self,'fieldname','love_h','format','DoubleMat','mattype',1);
|
---|
| 126 | WriteData(fid,prefix,'object',self,'fieldname','love_k','format','DoubleMat','mattype',1);
|
---|
[21000] | 127 | WriteData(fid,prefix,'object',self,'fieldname','love_l','format','DoubleMat','mattype',1);
|
---|
[20891] | 128 | WriteData(fid,prefix,'object',self,'fieldname','tide_love_h','format','Double');
|
---|
[20690] | 129 | WriteData(fid,prefix,'object',self,'fieldname','tide_love_k','format','Double');
|
---|
| 130 | WriteData(fid,prefix,'object',self,'fieldname','rigid','format','Boolean');
|
---|
| 131 | WriteData(fid,prefix,'object',self,'fieldname','elastic','format','Boolean');
|
---|
| 132 | WriteData(fid,prefix,'object',self,'fieldname','rotation','format','Boolean');
|
---|
[21295] | 133 | WriteData(fid,prefix,'object',self,'fieldname','ocean_area_scaling','format','Boolean');
|
---|
[20690] | 134 | WriteData(fid,prefix,'object',self,'fieldname','degacc','format','Double');
|
---|
| 135 | WriteData(fid,prefix,'object',self,'fieldname','transitions','format','MatArray');
|
---|
[20036] | 136 |
|
---|
| 137 | %process requested outputs
|
---|
| 138 | outputs = self.requested_outputs;
|
---|
| 139 | pos = find(ismember(outputs,'default'));
|
---|
| 140 | if ~isempty(pos),
|
---|
| 141 | outputs(pos) = []; %remove 'default' from outputs
|
---|
| 142 | outputs = [outputs defaultoutputs(self,md)]; %add defaults
|
---|
| 143 | end
|
---|
[20690] | 144 | WriteData(fid,prefix,'data',outputs,'name','md.slr.requested_outputs','format','StringArray');
|
---|
[20036] | 145 |
|
---|
[19984] | 146 | end % }}}
|
---|
| 147 | function savemodeljs(self,fid,modelname) % {{{
|
---|
| 148 |
|
---|
[20181] | 149 | writejs1Darray(fid,[modelname '.slr.deltathickness'],self.deltathickness);
|
---|
| 150 | writejs1Darray(fid,[modelname '.slr.sealevel'],self.sealevel);
|
---|
[20262] | 151 | writejsdouble(fid,[modelname '.slr.maxiter'],self.maxiter);
|
---|
[19984] | 152 | writejsdouble(fid,[modelname '.slr.reltol'],self.reltol);
|
---|
| 153 | writejsdouble(fid,[modelname '.slr.abstol'],self.abstol);
|
---|
[20181] | 154 | writejs1Darray(fid,[modelname '.slr.love_h'],self.love_h);
|
---|
| 155 | writejs1Darray(fid,[modelname '.slr.love_k'],self.love_k);
|
---|
[21000] | 156 | writejs1Darray(fid,[modelname '.slr.love_l'],self.love_l);
|
---|
[20349] | 157 | writejsdouble(fid,[modelname '.slr.tide_love_k'],self.tide_love_k);
|
---|
| 158 | writejsdouble(fid,[modelname '.slr.tide_love_h'],self.tide_love_h);
|
---|
[19984] | 159 | writejsdouble(fid,[modelname '.slr.rigid'],self.rigid);
|
---|
[21295] | 160 | writejsdouble(fid,[modelname '.slr.elastic'],self.elastic);
|
---|
[20349] | 161 | writejsdouble(fid,[modelname '.slr.rotation'],self.rotation);
|
---|
[21295] | 162 | writejsdouble(fid,[modelname '.slr.ocean_area_scaling'],self.ocean_area_scaling);
|
---|
[20033] | 163 | writejsdouble(fid,[modelname '.slr.degacc'],self.degacc);
|
---|
[20036] | 164 | writejscellstring(fid,[modelname '.slr.requested_outputs'],self.requested_outputs);
|
---|
[20262] | 165 | writejscellarray(fid,[modelname '.slr.transitions'],self.transitions);
|
---|
[19984] | 166 | end % }}}
|
---|
| 167 | end
|
---|
| 168 | end
|
---|