Index: /issm/trunk-jpl/src/m/Makefile.am =================================================================== --- /issm/trunk-jpl/src/m/Makefile.am (revision 20180) +++ /issm/trunk-jpl/src/m/Makefile.am (revision 20181) @@ -119,4 +119,5 @@ ${ISSM_DIR}/src/m/classes/radaroverlay.js \ ${ISSM_DIR}/src/m/classes/rifts.js \ + ${ISSM_DIR}/src/m/classes/slr.js \ ${ISSM_DIR}/src/m/classes/settings.js \ ${ISSM_DIR}/src/m/classes/SMBforcing.js \ Index: /issm/trunk-jpl/src/m/boundaryconditions/love_numbers.m =================================================================== --- /issm/trunk-jpl/src/m/boundaryconditions/love_numbers.m (revision 20180) +++ /issm/trunk-jpl/src/m/boundaryconditions/love_numbers.m (revision 20181) @@ -17,5 +17,6 @@ love_numbers=[... - -0.28740059 -0.00000000 0.10141481 1.28740059 0.42519882 0.89858519 0.00000000; + 0 0 0 0 0 0 0 + -1.28740059 -1.00000000 -0.89858519 1.28740059 0.42519882 0.89858519 0.00000000 -1.00025365 -0.30922675 0.02060926 1.69102690 0.46358648 0.67016399 0.61829668; -1.06243501 -0.19927948 0.06801636 1.86315553 0.55741597 0.73270416 0.56270589; Index: /issm/trunk-jpl/src/m/classes/maskpsl.m =================================================================== --- /issm/trunk-jpl/src/m/classes/maskpsl.m (revision 20180) +++ /issm/trunk-jpl/src/m/classes/maskpsl.m (revision 20181) @@ -81,4 +81,12 @@ end % }}} + function savemodeljs(self,fid,modelname) % {{{ + + writejs1Darray(fid,[modelname '.mask.groundedice_levelset'],self.groundedice_levelset); + writejs2Darray(fid,[modelname '.mask.ice_levelset'],self.ice_levelset); + writejs2Darray(fid,[modelname '.mask.ocean_levelset'],self.ocean_levelset); + writejs2Darray(fid,[modelname '.mask.land_levelset'],self.land_levelset); + + end % }}} end end Index: /issm/trunk-jpl/src/m/classes/mesh3dsurface.m =================================================================== --- /issm/trunk-jpl/src/m/classes/mesh3dsurface.m (revision 20180) +++ /issm/trunk-jpl/src/m/classes/mesh3dsurface.m (revision 20181) @@ -166,4 +166,27 @@ z = self.z; end % }}} + function savemodeljs(self,fid,modelname) % {{{ + + writejs1Darray(fid,[modelname '.mesh.x'],self.x); + writejs1Darray(fid,[modelname '.mesh.y'],self.y); + writejs1Darray(fid,[modelname '.mesh.z'],self.z); + writejs2Darray(fid,[modelname '.mesh.elements'],self.elements); + writejsdouble(fid,[modelname '.mesh.numberofelements'],self.numberofelements); + writejsdouble(fid,[modelname '.mesh.numberofvertices'],self.numberofvertices); + writejsdouble(fid,[modelname '.mesh.numberofedges'],self.numberofedges); + writejs1Darray(fid,[modelname '.mesh.lat'],self.lat); + writejs1Darray(fid,[modelname '.mesh.long'],self.long); + writejs1Darray(fid,[modelname '.mesh.r'],self.r); + writejs1Darray(fid,[modelname '.mesh.vertexonboundary'],self.vertexonboundary); + writejs2Darray(fid,[modelname '.mesh.edges'],self.edges); + writejs2Darray(fid,[modelname '.mesh.segments'],self.segments); + writejs2Darray(fid,[modelname '.mesh.segmentmarkers'],self.segmentmarkers); + writejs2Darray(fid,[modelname '.mesh.vertexconnectivity'],self.vertexconnectivity); + writejs2Darray(fid,[modelname '.mesh.elementconnectivity'],self.elementconnectivity); + writejsdouble(fid,[modelname '.mesh.average_vertex_connectivity'],self.average_vertex_connectivity); + writejs1Darray(fid,[modelname '.mesh.extractedvertices'],self.extractedvertices); + writejs1Darray(fid,[modelname '.mesh.extractedelements'],self.extractedelements); + + end % }}} end end Index: /issm/trunk-jpl/src/m/classes/model.js =================================================================== --- /issm/trunk-jpl/src/m/classes/model.js (revision 20180) +++ /issm/trunk-jpl/src/m/classes/model.js (revision 20181) @@ -21,4 +21,5 @@ console.log(sprintf("%19s: %-22s -- %s","initialization" ,"[1x1 " + typeof(this.initialization) + "]","initial guess/state")); console.log(sprintf("%19s: %-22s -- %s","rifts" ,"[1x1 " + typeof(this.rifts) + "]","rifts properties")); + console.log(sprintf("%19s: %-22s -- %s","slr" ,"[1x1 " + typeof(this.slr) + "]","slr forcings")); console.log(sprintf("%19s: %-22s -- %s","debug" ,"[1x1 " + typeof(this.debug) + "]","debugging tools (valgrind, gprof)")); console.log(sprintf("%19s: %-22s -- %s","verbose" ,"[1x1 " + typeof(this.verbose) + "]","verbosity level in solve")); @@ -57,4 +58,5 @@ this.friction = new friction(); this.rifts = new rifts(); + this.slr = new slr(); this.timestepping = new timestepping(); this.groundingline = new groundingline(); @@ -124,4 +126,5 @@ this.initialization = 0; this.rifts = 0; + this.slr = 0; this.debug = 0; Index: /issm/trunk-jpl/src/m/classes/slr.js =================================================================== --- /issm/trunk-jpl/src/m/classes/slr.js (revision 20181) +++ /issm/trunk-jpl/src/m/classes/slr.js (revision 20181) @@ -0,0 +1,129 @@ +//SLR class definition +// +// Usage: +// slr=slr(); + +function slr(){ + //methods + this.setdefaultparameters = function (){ //{{{ + + //Convergence criterion: absolute, relative and residual + this.reltol=NaN; //default + this.abstol=0.001; //1 mm of sea level rise + + //maximum of non-linear iterations. + this.maxiter=10; + + //computational flags: + this.rigid=1; + this.elastic=1; + this.eustatic=1; + + //numerical discretization accuracy + this.degacc=.01; + + //output default: + this.requested_outputs=['default']; + + //transitions should be a cell array of vectors: + this.transitions=[]; + + }// }}} + this.checkconsistency = function(md,solution,analyses) { //{{{ + + //Early return + if(ArrayAnyEqual(ArrayIsMember(SealevelriseAnalysisEnum(),analyses),0))return; + + 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.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 (this.love_h.length != this.love_k.length){ + throw Error('slr error message: love numbers should be provided at the same level of accuracy'); + } + + } // }}} + this.defaultoutputs = function(md){ // {{{ + return ['Sealevel']; + }//}}} + this.disp= function(){// {{{ + + console.log(sprintf(' Sealevelrise solution parameters:')); + + fielddisplay(this,'deltathickness','thickness change (main loading of the slr solution core [m]'); + fielddisplay(this,'sealevel','current sea level (prior to computation) [m]'); + fielddisplay(this,'reltol','sea level rise relative convergence criterion, (default, NaN: not applied)'); + fielddisplay(this,'abstol','sea level rise absolute convergence criterion, NaN: not applied'); + fielddisplay(this,'maxiter','maximum number of nonlinear iterations'); + fielddisplay(this,'love_h','love load number for radial displacement'); + fielddisplay(this,'love_k','love load number for gravitational potential perturbation'); + fielddisplay(this,'rigid','rigid earth graviational potential perturbation'); + fielddisplay(this,'elastic','elastic earth graviational potential perturbation'); + fielddisplay(this,'eustatic','eustatic sea level rise'); + fielddisplay(this,'degacc',"accuracy (default .01 deg) for numerical discretization of the Green's functions"); + fielddisplay(this,'transitions','indices into parts of the mesh that will be icecaps'); + fielddisplay(this,'requested_outputs','additional outputs requested'); + } //}}} + this.marshall=function(md,fid) { //{{{ + + WriteData(fid,'object',this,'class','sealevelrise','fieldname','deltathickness','format','DoubleMat','mattype',1); + WriteData(fid,'data',this.sealevel,'mattype',1,'format','DoubleMat','enum',SealevelEnum(),'timeserieslength',md.mesh.numberofvertices+1); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','reltol','format','Double'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','abstol','format','Double'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','maxiter','format','Integer'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','love_h','format','DoubleMat','mattype',1); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','love_k','format','DoubleMat','mattype',1); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','rigid','format','Boolean'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','elastic','format','Boolean'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','eustatic','format','Boolean'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','degacc','format','Double'); + WriteData(fid,'object',this,'class','sealevelrise','fieldname','transitions','format','MatArray'); + + //process requested outputs + var outputs = this.requested_outputs; + for (var i=0;i