Index: /issm/trunk-jpl/src/m/classes/SMB.m =================================================================== --- /issm/trunk-jpl/src/m/classes/SMB.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/SMB.m (revision 17079) @@ -15,4 +15,9 @@ error('constructor not supported'); end + end % }}} + function self = extrude(self,md) % {{{ + + self.mass_balance=project3d(md,'vector',self.mass_balance,'type','node'); + end % }}} function self = initialize(self,md) % {{{ Index: /issm/trunk-jpl/src/m/classes/SMB.py =================================================================== --- /issm/trunk-jpl/src/m/classes/SMB.py (revision 17078) +++ /issm/trunk-jpl/src/m/classes/SMB.py (revision 17079) @@ -2,4 +2,5 @@ from EnumDefinitions import * from checkfield import * +from project3d import * from WriteData import * @@ -20,4 +21,9 @@ return string #}}} + def extrude(self,md): # {{{ + + self.mass_balance=project3d(md,'vector',self.mass_balance,'type','node'); + return self + #}}} def initialize(self,md): # {{{ Index: /issm/trunk-jpl/src/m/classes/SMBgradients.m =================================================================== --- /issm/trunk-jpl/src/m/classes/SMBgradients.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/SMBgradients.m (revision 17079) @@ -19,4 +19,9 @@ error('constructor not supported'); end + end % }}} + function self = extrude(self,md) % {{{ + + %Nothing for now + end % }}} function self = initialize(self,md) % {{{ Index: /issm/trunk-jpl/src/m/classes/SMBgradients.py =================================================================== --- /issm/trunk-jpl/src/m/classes/SMBgradients.py (revision 17078) +++ /issm/trunk-jpl/src/m/classes/SMBgradients.py (revision 17079) @@ -2,4 +2,5 @@ from EnumDefinitions import * from checkfield import * +from project3d import * from WriteData import * @@ -29,4 +30,9 @@ return string #}}} + def extrude(self,md): # {{{ + + #Nothing for now + return self + #}}} def initialize(self,md): # {{{ Index: /issm/trunk-jpl/src/m/classes/SMBpdd.m =================================================================== --- /issm/trunk-jpl/src/m/classes/SMBpdd.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/SMBpdd.m (revision 17079) @@ -25,4 +25,14 @@ error('constructor not supported'); end + end % }}} + function self = extrude(self,md) % {{{ + + self.precipitation=project3d(md,'vector',self.precipitation,'type','node'); + self.monthlytemperatures=project3d(md,'vector',self.monthlytemperatures,'type','node'); + if(self.isdelta18o),self.temperatures_lgm=project3d(md,'vector',self.temperatures_lgm,'type','node');end + if(self.isdelta18o),self.temperatures_presentday=project3d(md,'vector',self.temperatures_presentday,'type','node');end + if(self.isdelta18o),self.precipitations_presentday=project3d(md,'vector',self.precipitations_presentday,'type','node');end + + end % }}} function self = initialize(self,md) % {{{ Index: /issm/trunk-jpl/src/m/classes/SMBpdd.py =================================================================== --- /issm/trunk-jpl/src/m/classes/SMBpdd.py (revision 17078) +++ /issm/trunk-jpl/src/m/classes/SMBpdd.py (revision 17079) @@ -2,4 +2,5 @@ from EnumDefinitions import * from checkfield import * +from project3d import * from WriteData import * @@ -43,4 +44,13 @@ return string #}}} + def extrude(self,md): # {{{ + + self.precipitation=project3d(md,'vector',self.precipitation,'type','node'); + self.monthlytemperatures=project3d(md,'vector',self.monthlytemperatures,'type','node'); + if self.isdelta18o: self.temperatures_lgm=project3d(md,'vector',self.temperatures_lgm,'type','node') + if self.isdelta18o: self.temperatures_presentday=project3d(md,'vector',self.temperatures_presentday,'type','node') + if self.isdelta18o: self.precipitations_presentday=project3d(md,'vector',self.precipitations_presentday,'type','node') + return self + #}}} def initialize(self,md): # {{{ Index: /issm/trunk-jpl/src/m/classes/model.m =================================================================== --- /issm/trunk-jpl/src/m/classes/model.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/model.m (revision 17079) @@ -103,4 +103,11 @@ if ~isa(md.outputdefinition,'outputdefinition'), md.outputdefinition=outputdefinition(); + end + %2014 January 9th + if isa(md.surfaceforcings,'surfaceforcings'), + disp('Recovering old surfaceforcings class'); + mass_balance=md.surfaceforcings.mass_balance; + md.surfaceforcings=SMB(); + md.surfaceforcings.mass_balance=mass_balance; end end% }}} @@ -673,8 +680,6 @@ md.inversion.vel_obs=project3d(md,'vector',md.inversion.vel_obs,'type','node'); md.inversion.thickness_obs=project3d(md,'vector',md.inversion.thickness_obs,'type','node'); - md.surfaceforcings.mass_balance=project3d(md,'vector',md.surfaceforcings.mass_balance,'type','node'); - md.surfaceforcings.precipitation=project3d(md,'vector',md.surfaceforcings.precipitation,'type','node'); + md.surfaceforcings = extrude(md.surfaceforcings,md); md.balancethickness.thickening_rate=project3d(md,'vector',md.balancethickness.thickening_rate,'type','node'); - md.surfaceforcings.monthlytemperatures=project3d(md,'vector',md.surfaceforcings.monthlytemperatures,'type','node'); %results @@ -760,7 +765,4 @@ if ~isnan(md.inversion.max_parameters),md.inversion.max_parameters=project3d(md,'vector',md.inversion.max_parameters,'type','node');end; if ~isnan(md.qmu.partition),md.qmu.partition=project3d(md,'vector',md.qmu.partition','type','node');end - if(md.surfaceforcings.isdelta18o),md.surfaceforcings.temperatures_lgm=project3d(md,'vector',md.surfaceforcings.temperatures_lgm,'type','node');end - if(md.surfaceforcings.isdelta18o),md.surfaceforcings.temperatures_presentday=project3d(md,'vector',md.surfaceforcings.temperatures_presentday,'type','node');end - if(md.surfaceforcings.isdelta18o),md.surfaceforcings.precipitations_presentday=project3d(md,'vector',md.surfaceforcings.precipitations_presentday,'type','node');end %Put lithostatic pressure if there is an existing pressure @@ -1064,5 +1066,5 @@ md.geometry = geometry(); md.initialization = initialization(); - md.surfaceforcings = surfaceforcings(); + md.surfaceforcings = SMB(); md.basalforcings = basalforcings(); md.friction = friction(); @@ -1079,7 +1081,7 @@ md.cluster = generic(); md.balancethickness = balancethickness(); - md.stressbalance = stressbalance(); + md.stressbalance = stressbalance(); md.hydrology = hydrologyshreve(); - md.masstransport = masstransport(); + md.masstransport = masstransport(); md.thermal = thermal(); md.steadystate = steadystate(); Index: /issm/trunk-jpl/src/m/classes/model.py =================================================================== --- /issm/trunk-jpl/src/m/classes/model.py (revision 17078) +++ /issm/trunk-jpl/src/m/classes/model.py (revision 17079) @@ -60,5 +60,5 @@ self.geometry = geometry() self.constants = constants() - self.surfaceforcings = surfaceforcings() + self.surfaceforcings = SMB() self.basalforcings = basalforcings() self.materials = matice() @@ -576,8 +576,6 @@ md.inversion.vel_obs=project3d(md,'vector',md.inversion.vel_obs,'type','node') md.inversion.thickness_obs=project3d(md,'vector',md.inversion.thickness_obs,'type','node') - md.surfaceforcings.mass_balance=project3d(md,'vector',md.surfaceforcings.mass_balance,'type','node') - md.surfaceforcings.precipitation=project3d(md,'vector',md.surfaceforcings.precipitation,'type','node') + md.surfaceforcings.extrude(md) md.balancethickness.thickening_rate=project3d(md,'vector',md.balancethickness.thickening_rate,'type','node') - md.surfaceforcings.monthlytemperatures=project3d(md,'vector',md.surfaceforcings.monthlytemperatures,'type','node') #results @@ -662,10 +660,4 @@ if not numpy.any(numpy.isnan(md.qmu.partition)): md.qmu.partition=project3d(md,'vector',numpy.transpose(md.qmu.partition),'type','node') - if(md.surfaceforcings.isdelta18o): - md.surfaceforcings.temperatures_lgm=project3d(md,'vector',md.surfaceforcings.temperatures_lgm,'type','node') - if(md.surfaceforcings.isdelta18o): - md.surfaceforcings.temperatures_presentday=project3d(md,'vector',md.surfaceforcings.temperatures_presentday,'type','node') - if(md.surfaceforcings.isdelta18o): - md.surfaceforcings.precipitations_presentday=project3d(md,'vector',md.surfaceforcings.precipitations_presentday,'type','node') #Put lithostatic pressure if there is an existing pressure Index: /issm/trunk-jpl/src/m/classes/oldclasses/README =================================================================== --- /issm/trunk-jpl/src/m/classes/oldclasses/README (revision 17078) +++ /issm/trunk-jpl/src/m/classes/oldclasses/README (revision 17079) @@ -1,3 +1,4 @@ We put here all old classes so that matlab can still load old models When an object is renamed or deleted, matlab is not able to load it -anymore. We keep it here so that matlab can load the object +anymore. We keep it here so that matlab can load the object (just remove all +methods) Index: /issm/trunk-jpl/src/m/classes/oldclasses/diagnostic.m =================================================================== --- /issm/trunk-jpl/src/m/classes/oldclasses/diagnostic.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/oldclasses/diagnostic.m (revision 17079) @@ -26,107 +26,3 @@ requested_outputs = NaN; end - methods - function obj = diagnostic(varargin) % {{{ - switch nargin - case 0 - obj=setdefaultparameters(obj); - otherwise - error('constructor not supported'); - end - end % }}} - function obj = setdefaultparameters(obj) % {{{ - - %maximum of non-linear iterations. - obj.maxiter=100; - - %Convergence criterion: absolute, relative and residual - obj.restol=10^-4; - obj.reltol=0.01; - obj.abstol=10; - - obj.FSreconditioning=10^13; - obj.shelf_dampening=0; - - %Penalty factor applied kappa=max(stiffness matrix)*10^penalty_factor - obj.penalty_factor=3; - - %coefficient to update the viscosity between each iteration of - %a diagnostic according to the following formula - %viscosity(n)=viscosity(n)+viscosity_overshoot(viscosity(n)-viscosity(n-1)) - obj.viscosity_overshoot=0; - - %Stop the iterations of rift if below a threshold - obj.rift_penalty_threshold=0; - - %in some solutions, it might be needed to stop a run when only - %a few constraints remain unstable. For thermal computation, this - %parameter is often used. - obj.rift_penalty_lock=10; - - end % }}} - function disp(obj) % {{{ - - disp(sprintf(' Stressbalance solution parameters:')); - - disp(sprintf('\n %s','Convergence criteria:')); - fielddisplay(obj,'restol','mechanical equilibrium residual convergence criterion'); - fielddisplay(obj,'reltol','velocity relative convergence criterion, NaN: not applied'); - fielddisplay(obj,'abstol','velocity absolute convergence criterion, NaN: not applied'); - fielddisplay(obj,'isnewton','0: Picard''s fixed point, 1: Newton''s method, 2: hybrid'); - fielddisplay(obj,'maxiter','maximum number of nonlinear iterations'); - fielddisplay(obj,'viscosity_overshoot','over-shooting constant new=new+C*(new-old)'); - - disp(sprintf('\n %s','boundary conditions:')); - fielddisplay(obj,'spcvx','x-axis velocity constraint (NaN means no constraint) [m/yr]'); - fielddisplay(obj,'spcvy','y-axis velocity constraint (NaN means no constraint) [m/yr]'); - fielddisplay(obj,'spcvz','z-axis velocity constraint (NaN means no constraint) [m/yr]'); - fielddisplay(obj,'icefront','segments on ice front list (last column 0: Air, 1: Water, 2: Ice)'); - - disp(sprintf('\n %s','Rift options:')); - fielddisplay(obj,'rift_penalty_threshold','threshold for instability of mechanical constraints'); - fielddisplay(obj,'rift_penalty_lock','number of iterations before rift penalties are locked'); - - disp(sprintf('\n %s','Penalty options:')); - fielddisplay(obj,'penalty_factor','offset used by penalties: penalty = Kmax*10^offset'); - fielddisplay(obj,'vertex_pairing','pairs of vertices that are penalized'); - - disp(sprintf('\n %s','Other:')); - fielddisplay(obj,'shelf_dampening','use dampening for floating ice ? Only for FS model'); - fielddisplay(obj,'FSreconditioning','multiplier for incompressibility equation. Only for FS model'); - fielddisplay(obj,'referential','local referential'); - fielddisplay(obj,'loadingforce','loading force applied on each point [N/m^3]'); - fielddisplay(obj,'requested_outputs','additional outputs requested'); - - end % }}} - - function marshall(obj,md,fid) % {{{ - WriteData(fid,'object',obj,'fieldname','spcvx','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','spcvy','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','spcvz','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','restol','format','Double'); - WriteData(fid,'object',obj,'fieldname','reltol','format','Double'); - WriteData(fid,'object',obj,'fieldname','abstol','format','Double'); - WriteData(fid,'object',obj,'fieldname','isnewton','format','Integer'); - WriteData(fid,'object',obj,'fieldname','FSreconditioning','format','Double'); - WriteData(fid,'object',obj,'fieldname','viscosity_overshoot','format','Double'); - WriteData(fid,'object',obj,'fieldname','maxiter','format','Integer'); - WriteData(fid,'object',obj,'fieldname','shelf_dampening','format','Integer'); - WriteData(fid,'object',obj,'fieldname','vertex_pairing','format','DoubleMat','mattype',3); - WriteData(fid,'object',obj,'fieldname','penalty_factor','format','Double'); - WriteData(fid,'object',obj,'fieldname','rift_penalty_lock','format','Integer'); - WriteData(fid,'object',obj,'fieldname','rift_penalty_threshold','format','Integer'); - WriteData(fid,'object',obj,'fieldname','referential','format','DoubleMat','mattype',1); - WriteData(fid,'data',obj.loadingforce(:,1),'format','DoubleMat','mattype',1,'enum',LoadingforceXEnum); - WriteData(fid,'data',obj.loadingforce(:,2),'format','DoubleMat','mattype',1,'enum',LoadingforceYEnum); - WriteData(fid,'data',obj.loadingforce(:,3),'format','DoubleMat','mattype',1,'enum',LoadingforceZEnum); - WriteData(fid,'object',obj,'fieldname','requested_outputs','format','DoubleMat','mattype',3); - - %marshall ice front - data=obj.icefront; - pos=find(data(:,end)==0); data(pos,end)=AirEnum(); - pos=find(data(:,end)==1); data(pos,end)=WaterEnum(); - pos=find(data(:,end)==2); data(pos,end)=IceEnum(); - WriteData(fid,'data',data,'enum',StressbalanceIcefrontEnum(),'format','DoubleMat','mattype',3); - end % }}} - end end Index: /issm/trunk-jpl/src/m/classes/oldclasses/hydrology.m =================================================================== --- /issm/trunk-jpl/src/m/classes/oldclasses/hydrology.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/oldclasses/hydrology.m (revision 17079) @@ -14,54 +14,3 @@ stabilization = 0; end - methods - function obj = hydrology(varargin) % {{{ - switch nargin - case 0 - obj=setdefaultparameters(obj); - otherwise - error('constructor not supported'); - end - end % }}} - function obj = setdefaultparameters(obj) % {{{ - - %Parameters from Johnson's 2002 thesis, section 3.5.4 - obj.n=.02; - obj.CR=.01; - obj.p=2; - obj.q=1; - obj.kn=0; - - %Type of stabilization to use 0:nothing 1:artificial_diffusivity - obj.stabilization=1; - end % }}} - function md = checkconsistency(obj,md,solution,analyses) % {{{ - - %Early return - if ~ismember(HydrologyAnalysisEnum(),analyses), return; end - - md = checkfield(md,'hydrology.spcwatercolumn','forcing',1); - md = checkfield(md,'hydrology.stabilization','>=',0); - end % }}} - function disp(obj) % {{{ - disp(sprintf(' hydrology solution parameters:')); - - fielddisplay(obj,'spcwatercolumn','water thickness constraints (NaN means no constraint)'); - fielddisplay(obj,'n','Manning roughness coefficient'); - fielddisplay(obj,'CR','tortuosity parameter'); - fielddisplay(obj,'p','dimensionless exponent in Manning velocity formula'); - fielddisplay(obj,'q','dimensionless exponent in Manning velocity formula'); - fielddisplay(obj,'kn','parameter in effective pressure formula'); - fielddisplay(obj,'stabilization','artificial diffusivity (default is 1). can be more than 1 to increase diffusivity.'); - - end % }}} - function marshall(obj,md,fid) % {{{ - WriteData(fid,'object',obj,'fieldname','spcwatercolumn','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','n','format','Double'); - WriteData(fid,'object',obj,'fieldname','CR','format','Double'); - WriteData(fid,'object',obj,'fieldname','p','format','Double'); - WriteData(fid,'object',obj,'fieldname','q','format','Double'); - WriteData(fid,'object',obj,'fieldname','kn','format','Double'); - WriteData(fid,'object',obj,'fieldname','stabilization','format','Double'); - end % }}} - end end Index: /issm/trunk-jpl/src/m/classes/oldclasses/prognostic.m =================================================================== --- /issm/trunk-jpl/src/m/classes/oldclasses/prognostic.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/oldclasses/prognostic.m (revision 17079) @@ -14,63 +14,3 @@ requested_outputs = NaN; end - methods - function obj = prognostic(varargin) % {{{ - switch nargin - case 0 - obj=setdefaultparameters(obj); - otherwise - error('constructor not supported'); - end - end % }}} - function obj = setdefaultparameters(obj) % {{{ - - %Type of stabilization to use 0:nothing 1:artificial_diffusivity 3:Discontinuous Galerkin - obj.stabilization=1; - - %Factor applied to compute the penalties kappa=max(stiffness matrix)*10^penalty_factor - obj.penalty_factor=3; - - %Minimum ice thickness that can be used - obj.min_thickness=1; - - %Hydrostatic adjustment - obj.hydrostatic_adjustment='Absolute'; - end % }}} - function md = checkconsistency(obj,md,solution,analyses) % {{{ - - %Early return - if ~ismember(PrognosticAnalysisEnum(),analyses) | (solution==TransientSolutionEnum() & md.transient.isprognostic==0), return; end - - md = checkfield(md,'prognostic.spcthickness','forcing',1); - md = checkfield(md,'prognostic.hydrostatic_adjustment','values',{'Absolute' 'Incremental'}); - md = checkfield(md,'prognostic.stabilization','values',[0 1 2 3]); - md = checkfield(md,'prognostic.min_thickness','>',0); - if ~isempty(md.prognostic.requested_outputs), - md = checkfield(md,'prognostic.requested_outputs','size',[NaN 1]); - end - - end % }}} - function disp(obj) % {{{ - disp(sprintf(' Prognostic solution parameters:')); - fielddisplay(obj,'spcthickness','thickness constraints (NaN means no constraint) [m]'); - fielddisplay(obj,'min_thickness','minimum ice thickness allowed [m]'); - fielddisplay(obj,'hydrostatic_adjustment','adjustment of ice shelves surface and bed elevations: ''Incremental'' or ''Absolute'' '); - fielddisplay(obj,'stabilization','0: no, 1: artificial_diffusivity, 2: streamline upwinding, 3: discontinuous Galerkin'); - - disp(sprintf('\n %s','Penalty options:')); - fielddisplay(obj,'penalty_factor','offset used by penalties: penalty = Kmax*10^offset'); - fielddisplay(obj,'vertex_pairing','pairs of vertices that are penalized'); - fielddisplay(obj,'requested_outputs','additional outputs requested'); - - end % }}} - function marshall(obj,md,fid) % {{{ - WriteData(fid,'object',obj,'fieldname','spcthickness','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','min_thickness','format','Double'); - WriteData(fid,'data',StringToEnum(obj.hydrostatic_adjustment),'format','Integer','enum',PrognosticHydrostaticAdjustmentEnum()); - WriteData(fid,'object',obj,'fieldname','stabilization','format','Integer'); - WriteData(fid,'object',obj,'fieldname','vertex_pairing','format','DoubleMat','mattype',3); - WriteData(fid,'object',obj,'fieldname','penalty_factor','format','Double'); - WriteData(fid,'object',obj,'fieldname','requested_outputs','format','DoubleMat','mattype',3); - end % }}} - end end Index: /issm/trunk-jpl/src/m/classes/oldclasses/surfaceforcings.m =================================================================== --- /issm/trunk-jpl/src/m/classes/oldclasses/surfaceforcings.m (revision 17079) +++ /issm/trunk-jpl/src/m/classes/oldclasses/surfaceforcings.m (revision 17079) @@ -0,0 +1,26 @@ +%SURFACEFORCING Class definition +% +% Usage: +% surfaceforcings=surfaceforcings(); + +classdef surfaceforcings + properties (SetAccess=public) + precipitation = NaN; + mass_balance = NaN; + ispdd = 0; + issmbgradients = 0; + isdelta18o = 0; + href = NaN; + smbref = NaN; + b_pos = NaN; + b_neg = NaN; + monthlytemperatures = NaN; + delta18o = NaN; + delta18o_surface = NaN; + temperatures_presentday = NaN; + temperatures_lgm = NaN; + precipitations_presentday = NaN; + desfac = 0.5; + s0p = 0; + end +end Index: sm/trunk-jpl/src/m/classes/surfaceforcings.m =================================================================== --- /issm/trunk-jpl/src/m/classes/surfaceforcings.m (revision 17078) +++ (revision ) @@ -1,147 +1,0 @@ -%SURFACEFORCING Class definition -% -% Usage: -% surfaceforcings=surfaceforcings(); - -classdef surfaceforcings - properties (SetAccess=public) - precipitation = NaN; - mass_balance = NaN; - ispdd = 0; - issmbgradients = 0; - isdelta18o = 0; - href = NaN; - smbref = NaN; - b_pos = NaN; - b_neg = NaN; - monthlytemperatures = NaN; - delta18o = NaN; - delta18o_surface = NaN; - temperatures_presentday = NaN; - temperatures_lgm = NaN; - precipitations_presentday = NaN; - desfac = 0.5; - s0p = 0; - end - methods - function obj = surfaceforcings(varargin) % {{{ - switch nargin - case 0 - obj=setdefaultparameters(obj); - otherwise - error('constructor not supported'); - end - end % }}} - function self = initialize(self,md) % {{{ - - if (isnan(self.precipitation) & (self.ispdd==1)), - self.precipitation=zeros(md.mesh.numberofvertices,1); - disp(' no surfaceforcings.precipitation specified: values set as zero'); - end - if (isnan(self.mass_balance) & (self.ispdd==0)), - self.mass_balance=zeros(md.mesh.numberofvertices,1); - disp(' no surfaceforcings.mass_balance specified: values set as zero'); - end - - end % }}} - function obj = setdefaultparameters(obj) % {{{ - - %pdd method not used in default mode - obj.ispdd=0; - obj.issmbgradients=0; - obj.isdelta18o=0; - obj.desfac=0.5; - obj.s0p=0; - end % }}} - function md = checkconsistency(obj,md,solution,analyses) % {{{ - - if ismember(MasstransportAnalysisEnum(),analyses), - md = checkfield(md,'fieldname','surfaceforcings.ispdd','numel',[1],'values',[0 1]); - md = checkfield(md,'fieldname','surfaceforcings.issmbgradients','numel',[1],'values',[0 1]); - if(obj.ispdd) - md = checkfield(md,'fieldname','surfaceforcings.desfac','<=',1,'numel',1); - md = checkfield(md,'fieldname','surfaceforcings.s0p','>=',0,'numel',1); - if(obj.isdelta18o==0) - md = checkfield(md,'fieldname','surfaceforcings.monthlytemperatures','forcing',1,'NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.precipitation','forcing',1,'NaN',1); - else - md = checkfield(md,'fieldname','surfaceforcings.delta18o','NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.delta18o_surface','NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.temperatures_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.temperatures_lgm','size',[md.mesh.numberofvertices+1 12],'NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.precipitations_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1); - end - elseif(obj.issmbgradients) - md = checkfield(md,'fieldname','surfaceforcings.href','forcing',1,'NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.smbref','forcing',1,'NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.b_pos','forcing',1,'NaN',1); - md = checkfield(md,'fieldname','surfaceforcings.b_neg','forcing',1,'NaN',1); - else - md = checkfield(md,'fieldname','surfaceforcings.mass_balance','forcing',1,'NaN',1); - end - end - if ismember(BalancethicknessAnalysisEnum(),analyses), - md = checkfield(md,'fieldname','surfaceforcings.mass_balance','size',[md.mesh.numberofvertices 1],'NaN',1); - end - end % }}} - function disp(obj) % {{{ - disp(sprintf(' surface forcings parameters:')); - - disp(sprintf('\n General parameters:')); - fielddisplay(obj,'mass_balance','surface mass balance [m/yr ice eq]'); - - disp(sprintf('\n PDD and deltaO18 parameters:')); - fielddisplay(obj,'ispdd','is pdd activated (0 or 1, default is 0)'); - fielddisplay(obj,'isdelta18o','is temperature and precipitation delta18o parametrisation activated (0 or 1, default is 0)'); - fielddisplay(obj,'desfac','desertification elevation factor (between 0 and 1, default is 0.5) [m]'); - fielddisplay(obj,'s0p','should be set to elevation from precip source (between 0 and a few 1000s m, default is 0) [m]'); - fielddisplay(obj,'monthlytemperatures','monthly surface temperatures [K], required if pdd is activated and delta18o not activated'); - fielddisplay(obj,'precipitation','surface precipitation [m/yr water eq]'); - fielddisplay(obj,'temperatures_presentday','monthly present day surface temperatures [K], required if pdd is activated and delta18o activated'); - fielddisplay(obj,'temperatures_lgm','monthly LGM surface temperatures [K], required if pdd is activated and delta18o activated'); - fielddisplay(obj,'precipitations_presentday','monthly surface precipitation [m/yr water eq], required if pdd is activated and delta18o activated'); - fielddisplay(obj,'delta18o','delta18o, required if pdd is activated and delta18o activated'); - fielddisplay(obj,'delta18o_surface','surface elevation of the delta18o site, required if pdd is activated and delta18o activated'); - - disp(sprintf('\n SMB gradients parameters:')); - fielddisplay(obj,'issmbgradients','is smb gradients method activated (0 or 1, default is 0)'); - fielddisplay(obj,'href',' reference elevation from which deviation is used to calculate SMB adjustment in smb gradients method'); - fielddisplay(obj,'smbref',' reference smb from which deviation is calculated in smb gradients method'); - fielddisplay(obj,'b_pos',' slope of hs - smb regression line for accumulation regime required if smb gradients is activated'); - fielddisplay(obj,'b_neg',' slope of hs - smb regression line for ablation regime required if smb gradients is activated'); - - end % }}} - function marshall(obj,md,fid) % {{{ - - yts=365.0*24.0*3600.0; - - WriteData(fid,'enum',SurfaceforcingsEnum(),'data',SurfaceforcingsEnum(),'format','Integer'); - - WriteData(fid,'object',obj,'fieldname','precipitation','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','mass_balance','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','ispdd','format','Boolean'); - WriteData(fid,'object',obj,'fieldname','isdelta18o','format','Boolean'); - if obj.ispdd, - WriteData(fid,'object',obj,'fieldname','desfac','format','Double'); - WriteData(fid,'object',obj,'fieldname','s0p','format','Double'); - if obj.isdelta18o - WriteData(fid,'object',obj,'fieldname','temperatures_presentday','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','temperatures_lgm','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','precipitations_presentday','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','delta18o_surface','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','delta18o','format','DoubleMat','mattype',1); - else - WriteData(fid,'object',obj,'fieldname','monthlytemperatures','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1); - end - end - WriteData(fid,'object',obj,'fieldname','issmbgradients','format','Boolean'); - if obj.issmbgradients, - WriteData(fid,'object',obj,'fieldname','href','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','smbref','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','b_pos','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1); - WriteData(fid,'object',obj,'fieldname','b_neg','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1); - end - - end % }}} - end -end Index: sm/trunk-jpl/src/m/classes/surfaceforcings.py =================================================================== --- /issm/trunk-jpl/src/m/classes/surfaceforcings.py (revision 17078) +++ (revision ) @@ -1,143 +1,0 @@ -from fielddisplay import fielddisplay -from EnumDefinitions import * -from checkfield import * -from WriteData import * - -class surfaceforcings(object): - """ - SURFACEFORCING Class definition - - Usage: - surfaceforcings=surfaceforcings(); - """ - - def __init__(self): # {{{ - self.precipitation = float('NaN') - self.mass_balance = float('NaN') - self.ispdd = 0 - self.issmbgradients = 0 - self.isdelta18o = 0 - self.href = float('NaN') - self.smbref = float('NaN') - self.b_pos = float('NaN') - self.b_neg = float('NaN') - self.monthlytemperatures = float('NaN') - self.delta18o = float('NaN') - self.delta18o_surface = float('NaN') - self.temperatures_presentday = float('NaN') - self.temperatures_lgm = float('NaN') - self.precipitations_presentday = float('NaN') - self.desfac = 0. - self.s0p = 0. - - #set defaults - self.setdefaultparameters() - - #}}} - def __repr__(self): # {{{ - string=" surface forcings parameters:" - - string="%s\n%s"%(string,fielddisplay(self,'precipitation','surface precipitation [m/yr water eq]')) - string="%s\n%s"%(string,fielddisplay(self,'mass_balance','surface mass balance [m/yr ice eq]')) - string="%s\n%s"%(string,fielddisplay(self,'ispdd','is pdd activated (0 or 1, default is 0)')) - string="%s\n%s"%(string,fielddisplay(self,'desfac','desertification elevation factor (between 0 and 1, default is 0.5) [m]')) - string="%s\n%s"%(string,fielddisplay(self,'isdelta18o','is temperature and precipitation delta18o parametrisation activated (0 or 1, default is 0)')) - string="%s\n%s"%(string,fielddisplay(self,'s0p','should be set to elevation from precip source (between 0 and a few 1000s m, default is 0) [m]')) - string="%s\n%s"%(string,fielddisplay(self,'monthlytemperatures','monthly surface temperatures [K], required if pdd is activated and delta18o not activated')) - string="%s\n%s"%(string,fielddisplay(self,'precipitation','surface precipitation [m/yr water eq]')) - string="%s\n%s"%(string,fielddisplay(self,'temperatures_presentday','monthly present day surface temperatures [K], required if pdd is activated and delta18o activated')) - string="%s\n%s"%(string,fielddisplay(self,'temperatures_lgm','monthly LGM surface temperatures [K], required if pdd is activated and delta18o activated')) - string="%s\n%s"%(string,fielddisplay(self,'precipitations_presentday','monthly surface precipitation [m/yr water eq], required if pdd is activated and delta18o activated')) - string="%s\n%s"%(string,fielddisplay(self,'delta18o','delta18o, required if pdd is activated and delta18o activated')) - string="%s\n%s"%(string,fielddisplay(self,'delta18o_surface','surface elevation of the delta18o site, required if pdd is activated and delta18o activated')) - string="%s\n%s"%(string,fielddisplay(self,'issmbgradients','is smb gradients method activated (0 or 1, default is 0)')) - string="%s\n%s"%(string,fielddisplay(self,'href',' reference elevation from which deviation is used to calculate SMB adjustment in smb gradients method')) - string="%s\n%s"%(string,fielddisplay(self,'smbref',' reference smb from which deviation is calculated in smb gradients method')) - string="%s\n%s"%(string,fielddisplay(self,'b_pos',' slope of hs - smb regression line for accumulation regime required if smb gradients is activated')) - string="%s\n%s"%(string,fielddisplay(self,'b_neg',' slope of hs - smb regression line for ablation regime required if smb gradients is activated')) - - return string - #}}} - def initialize(self,md): # {{{ - - if numpy.all(numpy.isnan(self.precipitation)) and (self.ispdd==1): - self.precipitation=numpy.zeros((md.mesh.numberofvertices,1)) - print " no surfaceforcings.precipitation specified: values set as zero" - if numpy.all(numpy.isnan(self.mass_balance)) and (self.ispdd==0): - self.mass_balance=numpy.zeros((md.mesh.numberofvertices,1)) - print " no surfaceforcings.mass_balance specified: values set as zero" - - return self - #}}} - def setdefaultparameters(self): # {{{ - - #pdd method not used in default mode - self.ispdd=0 - self.issmbgradients=0 - self.isdelta18o=0 - self.desfac=0.5 - self.s0p=0. - - return self - #}}} - def checkconsistency(self,md,solution,analyses): # {{{ - - if MasstransportAnalysisEnum() in analyses: - md = checkfield(md,'fieldname','surfaceforcings.ispdd','numel',[1],'values',[0,1]) - md = checkfield(md,'fieldname','surfaceforcings.issmbgradients','numel',[1],'values',[0,1]) - if self.ispdd: - md = checkfield(md,'fieldname','surfaceforcings.desfac','<=',1,'numel',[1]); - md = checkfield(md,'fieldname','surfaceforcings.s0p','>=',0,'numel',[1]); - if not self.isdelta18o: - md = checkfield(md,'fieldname','surfaceforcings.monthlytemperatures','forcing',1,'NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.precipitation','forcing',1,'NaN',1) - else: - md = checkfield(md,'fieldname','surfaceforcings.delta18o','NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.delta18o_surface','NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.temperatures_presentday','size',[md.mesh.numberofvertices+1,12],'NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.temperatures_lgm','size',[md.mesh.numberofvertices+1,12],'NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.precipitations_presentday','size',[md.mesh.numberofvertices+1,12],'NaN',1) - elif self.issmbgradients: - md = checkfield(md,'fieldname','surfaceforcings.href','forcing',1,'NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.smbref','forcing',1,'NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.b_pos','forcing',1,'NaN',1) - md = checkfield(md,'fieldname','surfaceforcings.b_neg','forcing',1,'NaN',1) - else: - md = checkfield(md,'fieldname','surfaceforcings.mass_balance','forcing',1,'NaN',1) - - if BalancethicknessAnalysisEnum() in analyses: - md = checkfield(md,'fieldname','surfaceforcings.mass_balance','size',[md.mesh.numberofvertices],'NaN',1) - - return md - # }}} - def marshall(self,md,fid): # {{{ - - yts=365.0*24.0*3600.0 - - WriteData(fid,'enum',SurfaceforcingsEnum(),'data',SurfaceforcingsEnum(),'format','Integer'); - - WriteData(fid,'object',self,'fieldname','precipitation','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1) - WriteData(fid,'object',self,'fieldname','mass_balance','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1) - WriteData(fid,'object',self,'fieldname','ispdd','format','Boolean') - WriteData(fid,'object',self,'fieldname','isdelta18o','format','Boolean') - - if self.ispdd: - WriteData(fid,'object',self,'fieldname','desfac','format','Double'); - WriteData(fid,'object',self,'fieldname','s0p','format','Double'); - if self.isdelta18o: - WriteData(fid,'object',self,'fieldname','temperatures_presentday','format','DoubleMat','mattype',1) - WriteData(fid,'object',self,'fieldname','temperatures_lgm','format','DoubleMat','mattype',1) - WriteData(fid,'object',self,'fieldname','precipitations_presentday','format','DoubleMat','mattype',1) - WriteData(fid,'object',self,'fieldname','delta18o_surface','format','DoubleMat','mattype',1) - WriteData(fid,'object',self,'fieldname','delta18o','format','DoubleMat','mattype',1) - else: - WriteData(fid,'object',self,'fieldname','monthlytemperatures','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1) - - WriteData(fid,'object',self,'fieldname','issmbgradients','format','Boolean') - - if self.issmbgradients: - WriteData(fid,'object',self,'fieldname','href','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1) - WriteData(fid,'object',self,'fieldname','smbref','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1) - WriteData(fid,'object',self,'fieldname','b_pos','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1) - WriteData(fid,'object',self,'fieldname','b_neg','format','DoubleMat','mattype',1,'scale',1./yts,'forcinglength',md.mesh.numberofvertices+1) - # }}} Index: /issm/trunk-jpl/src/m/classes/transient.m =================================================================== --- /issm/trunk-jpl/src/m/classes/transient.m (revision 17078) +++ /issm/trunk-jpl/src/m/classes/transient.m (revision 17079) @@ -6,11 +6,11 @@ classdef transient properties (SetAccess=public) - ismasstransport = 0; - isstressbalance = 0; + ismasstransport = 0; + isstressbalance = 0; isthermal = 0; isgroundingline = 0; isgia = 0; - isdamage = 0; - islevelset = 0; + isdamage = 0; + islevelset = 0; requested_outputs = {}; end Index: /issm/trunk-jpl/src/m/classes/transient.py =================================================================== --- /issm/trunk-jpl/src/m/classes/transient.py (revision 17078) +++ /issm/trunk-jpl/src/m/classes/transient.py (revision 17079) @@ -13,6 +13,6 @@ def __init__(self): # {{{ - self.ismasstransport = False - self.isstressbalance = False + self.ismasstransport = False + self.isstressbalance = False self.isthermal = False self.isgroundingline = False