Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/analyses/SealevelriseAnalysis.cpp =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/analyses/SealevelriseAnalysis.cpp (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/analyses/SealevelriseAnalysis.cpp (revision 21329) @@ -19,4 +19,5 @@ }/*}}}*/ void SealevelriseAnalysis::UpdateElements(Elements* elements,IoModel* iomodel,int analysis_counter,int analysis_type){/*{{{*/ + /*Update elements: */ @@ -38,5 +39,15 @@ iomodel->FetchDataToInput(elements,"md.slr.spcthickness",SealevelriseSpcthicknessEnum); iomodel->FetchDataToInput(elements,"md.slr.sealevel",SealevelEnum,0); - + + /*Initialize cumdeltalthickness input: unfortunately, we don't have femmodel, so we need to iterate on + *elements, otherwise we would have used InputUpdateFromConstantx: */ + for(int i=0;inumberofelements;i++){ + if(iomodel->my_elements[i]){ + Element* element=(Element*)elements->GetObjectByOffset(counter); + element->InputUpdateFromConstant(0,SealevelriseCumDeltathicknessEnum); + counter++; + } + } + }/*}}}*/ void SealevelriseAnalysis::UpdateParameters(Parameters* parameters,IoModel* iomodel,int solution_enum,int analysis_enum){/*{{{*/ @@ -74,4 +85,6 @@ parameters->AddObject(iomodel->CopyConstantObject("md.slr.ocean_area_scaling",SealevelriseOceanAreaScalingEnum)); parameters->AddObject(iomodel->CopyConstantObject("md.slr.steric_rate",SealevelriseStericRateEnum)); + parameters->AddObject(iomodel->CopyConstantObject("md.slr.run_frequency",SealevelriseRunFrequencyEnum)); + parameters->SetParam(1,SealevelriseRunCountEnum); iomodel->FetchData(&elastic,"md.slr.elastic"); Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/cores.h =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/cores.h (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/cores.h (revision 21329) @@ -65,4 +65,5 @@ void TransferSealevel(FemModel* femmodel,int forcingenum); void EarthMassTransport(FemModel* femmodel); +bool SlrAccumulateDeltaThickness(FemModel* femmodel); //solution configuration Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/sealevelrise_core.cpp =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/sealevelrise_core.cpp (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/sealevelrise_core.cpp (revision 21329) @@ -88,5 +88,6 @@ /*call sea-level rise sub cores:*/ - if(isslr){ + if(isslr && !SlrAccumulateDeltaThickness(femmodel)){ + Sg_eustatic=sealevelrise_core_eustatic(femmodel); //generalized eustatic (Farrel and Clark, Eq 4, 1st, 3rd and 4rd terms on the RHS. @@ -412,2 +413,51 @@ } /*}}}*/ +bool SlrAccumulateDeltaThickness(FemModel* femmodel){ /*{{{*/ + + Vector* deltathickness=NULL; + Vector* cumdeltathickness=NULL; + IssmDouble frequency=1; + int count; + + femmodel->parameters->FindParam(&frequency,SealevelriseRunFrequencyEnum); + femmodel->parameters->FindParam(&count,SealevelriseRunCountEnum); + + if (count==frequency){ + femmodel->parameters->SetParam(1,SealevelriseRunCountEnum); //reset counter. + + /*increment delta thickness load: */ + GetVectorFromInputsx(&cumdeltathickness,femmodel,SealevelriseCumDeltathicknessEnum,VertexSIdEnum); + GetVectorFromInputsx(&deltathickness,femmodel,SealevelriseDeltathicknessEnum,VertexSIdEnum); + + cumdeltathickness->AXPY(deltathickness,1); + + /*update deltathickness with the final cumulated thickness: */ + InputUpdateFromVectorx(femmodel,cumdeltathickness,SealevelriseDeltathicknessEnum,VertexSIdEnum); + + /*reset cumulative thickness to 0: */ + cumdeltathickness->Set(0); + InputUpdateFromVectorx(femmodel,cumdeltathickness,SealevelriseCumDeltathicknessEnum,VertexSIdEnum); + + delete deltathickness; + delete cumdeltathickness; + + return false; + } + else{ + count++; + femmodel->parameters->SetParam(count,SealevelriseRunCountEnum); //increment counter + + /*increment cumulative delta thickness load: */ + GetVectorFromInputsx(&cumdeltathickness,femmodel,SealevelriseCumDeltathicknessEnum,VertexSIdEnum); + GetVectorFromInputsx(&deltathickness,femmodel,SealevelriseDeltathicknessEnum,VertexSIdEnum); + + cumdeltathickness->AXPY(deltathickness,1); + InputUpdateFromVectorx(femmodel,cumdeltathickness,SealevelriseCumDeltathicknessEnum,VertexSIdEnum); + + delete deltathickness; + delete cumdeltathickness; + + return true; + } + +} /*}}}*/ Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/modules/ModelProcessorx/Transient/UpdateElementsTransient.cpp (revision 21329) @@ -17,3 +17,4 @@ iomodel->FetchDataToInput(elements,"md.geometry.bed",BedEnum); } + } Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/EnumDefinitions.h (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/EnumDefinitions.h (revision 21329) @@ -767,4 +767,5 @@ SealevelEustaticEnum, SealevelriseDeltathicknessEnum, + SealevelriseCumDeltathicknessEnum, SealevelriseSpcthicknessEnum, SealevelriseMaxiterEnum, @@ -775,4 +776,6 @@ SealevelriseOceanAreaScalingEnum, SealevelriseStericRateEnum, + SealevelriseRunFrequencyEnum, + SealevelriseRunCountEnum, SealevelriseGElasticEnum, SealevelriseUElasticEnum, Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/EnumToStringx.cpp (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/EnumToStringx.cpp (revision 21329) @@ -748,4 +748,5 @@ case SealevelEustaticEnum : return "SealevelEustatic"; case SealevelriseDeltathicknessEnum : return "SealevelriseDeltathickness"; + case SealevelriseCumDeltathicknessEnum : return "SealevelriseCumDeltathickness"; case SealevelriseSpcthicknessEnum : return "SealevelriseSpcthickness"; case SealevelriseMaxiterEnum : return "SealevelriseMaxiter"; @@ -756,4 +757,6 @@ case SealevelriseOceanAreaScalingEnum : return "SealevelriseOceanAreaScaling"; case SealevelriseStericRateEnum : return "SealevelriseStericRate"; + case SealevelriseRunFrequencyEnum : return "SealevelriseRunFrequency"; + case SealevelriseRunCountEnum : return "SealevelriseRunCount"; case SealevelriseGElasticEnum : return "SealevelriseGElastic"; case SealevelriseUElasticEnum : return "SealevelriseUElastic"; Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/StringToEnumx.cpp (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/shared/Enum/StringToEnumx.cpp (revision 21329) @@ -766,4 +766,5 @@ else if (strcmp(name,"SealevelEustatic")==0) return SealevelEustaticEnum; else if (strcmp(name,"SealevelriseDeltathickness")==0) return SealevelriseDeltathicknessEnum; + else if (strcmp(name,"SealevelriseCumDeltathickness")==0) return SealevelriseCumDeltathicknessEnum; else if (strcmp(name,"SealevelriseSpcthickness")==0) return SealevelriseSpcthicknessEnum; else if (strcmp(name,"SealevelriseMaxiter")==0) return SealevelriseMaxiterEnum; @@ -774,4 +775,6 @@ else if (strcmp(name,"SealevelriseOceanAreaScaling")==0) return SealevelriseOceanAreaScalingEnum; else if (strcmp(name,"SealevelriseStericRate")==0) return SealevelriseStericRateEnum; + else if (strcmp(name,"SealevelriseRunFrequency")==0) return SealevelriseRunFrequencyEnum; + else if (strcmp(name,"SealevelriseRunCount")==0) return SealevelriseRunCountEnum; else if (strcmp(name,"SealevelriseGElastic")==0) return SealevelriseGElasticEnum; else if (strcmp(name,"SealevelriseUElastic")==0) return SealevelriseUElasticEnum; @@ -872,11 +875,11 @@ else if (strcmp(name,"HydrologyDCEfficientAnalysis")==0) return HydrologyDCEfficientAnalysisEnum; else if (strcmp(name,"HydrologySommersAnalysis")==0) return HydrologySommersAnalysisEnum; - else if (strcmp(name,"HydrologySolution")==0) return HydrologySolutionEnum; - else if (strcmp(name,"MeltingAnalysis")==0) return MeltingAnalysisEnum; - else if (strcmp(name,"MasstransportAnalysis")==0) return MasstransportAnalysisEnum; else stage=8; } if(stage==8){ - if (strcmp(name,"MasstransportSolution")==0) return MasstransportSolutionEnum; + if (strcmp(name,"HydrologySolution")==0) return HydrologySolutionEnum; + else if (strcmp(name,"MeltingAnalysis")==0) return MeltingAnalysisEnum; + else if (strcmp(name,"MasstransportAnalysis")==0) return MasstransportAnalysisEnum; + else if (strcmp(name,"MasstransportSolution")==0) return MasstransportSolutionEnum; else if (strcmp(name,"FreeSurfaceBaseAnalysis")==0) return FreeSurfaceBaseAnalysisEnum; else if (strcmp(name,"FreeSurfaceTopAnalysis")==0) return FreeSurfaceTopAnalysisEnum; Index: /issm/branches/trunk-larour-NatGeoScience2016/src/m/classes/slr.m =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/m/classes/slr.m (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/m/classes/slr.m (revision 21329) @@ -22,4 +22,5 @@ ocean_area_scaling = 0; steric_rate = 0; %rate of ocean expansion from steric effects. + run_frequency = 1; %how many time steps we skip before we run SLR solver during transient degacc = 0; requested_outputs = {}; @@ -59,5 +60,8 @@ %steric: self.steric_rate=0; - + + %how many time steps we skip before we run SLR solver during transient + self.run_frequency=1; + %output default: self.requested_outputs={'default'}; @@ -82,4 +86,5 @@ md = checkfield(md,'fieldname','slr.maxiter','size',[1 1],'>=',1); md = checkfield(md,'fieldname','slr.steric_rate','size',[1 1]); + md = checkfield(md,'fieldname','slr.run_frequency','size',[1 1],'>=',1); md = checkfield(md,'fieldname','slr.degacc','size',[1 1],'>=',1e-10); md = checkfield(md,'fieldname','slr.requested_outputs','stringrow',1); @@ -119,4 +124,5 @@ fielddisplay(self,'ocean_area_scaling','correction for model representation of ocean area [default: No correction]'); fielddisplay(self,'steric_rate','rate of steric ocean expansion (in mm/yr)'); + fielddisplay(self,'run_frequency','how many time steps we skip before we run SLR solver during transient (default: 1)'); fielddisplay(self,'rigid','rigid earth graviational potential perturbation'); fielddisplay(self,'elastic','elastic earth graviational potential perturbation'); @@ -143,4 +149,5 @@ WriteData(fid,prefix,'object',self,'fieldname','ocean_area_scaling','format','Boolean'); WriteData(fid,prefix,'object',self,'fieldname','steric_rate','format','Double','scale',1e-3/md.constants.yts); + WriteData(fid,prefix,'object',self,'fieldname','run_frequency','format','Integer'); WriteData(fid,prefix,'object',self,'fieldname','degacc','format','Double'); WriteData(fid,prefix,'object',self,'fieldname','transitions','format','MatArray'); @@ -173,4 +180,5 @@ writejsdouble(fid,[modelname '.slr.ocean_area_scaling'],self.ocean_area_scaling); writejsdouble(fid,[modelname '.slr.steric_rate'],self.steric_rate); + writejsdouble(fid,[modelname '.slr.run_frequency'],self.run_frequency); writejsdouble(fid,[modelname '.slr.elastic'],self.elastic); writejsdouble(fid,[modelname '.slr.degacc'],self.degacc); Index: /issm/branches/trunk-larour-NatGeoScience2016/src/m/plot/plot_contour.m =================================================================== --- /issm/branches/trunk-larour-NatGeoScience2016/src/m/plot/plot_contour.m (revision 21328) +++ /issm/branches/trunk-larour-NatGeoScience2016/src/m/plot/plot_contour.m (revision 21329) @@ -247,5 +247,5 @@ else %clabel(c,h,'color',color,'FontSize',10,'labelspacing',20000); - clabel(c,h,'FontSize',10,'labelspacing',20000); - end -end + %clabel(c,h,'FontSize',10,'labelspacing',20000); + end +end