Index: ../trunk-jpl/src/c/cores/masstransport_core.cpp =================================================================== --- ../trunk-jpl/src/c/cores/masstransport_core.cpp (revision 24096) +++ ../trunk-jpl/src/c/cores/masstransport_core.cpp (revision 24097) @@ -60,6 +60,7 @@ depthaverage_core(femmodel); } femmodel->SetCurrentConfiguration(MasstransportAnalysisEnum); + InputDuplicatex(femmodel,ThicknessEnum,ThicknessOldEnum); if(stabilization==4){ solutionsequence_fct(femmodel); } Index: ../trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- ../trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 24096) +++ ../trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 24097) @@ -1807,7 +1807,7 @@ /*Fetch dof list and allocate solution vector*/ GetDofListLocal(&doflist,NoneApproximationEnum,GsetEnum); - IssmDouble* values = xNew(numnodes); + IssmDouble* values = xNew(numnodes); /*Use the dof list to index into the solution vector: */ for(int i=0;iIsOnBase()) return; - /*Get basal element*/ - int domaintype; element->FindParam(&domaintype,DomainTypeEnum); - Element* basalelement=element; - if(domaintype!=Domain2DhorizontalEnum) basalelement = element->SpawnBasalElement(); - - /*Fetch number of nodes and dof for this finite element*/ - int numnodes = basalelement->GetNumberOfNodes(); - int numvertices = basalelement->GetNumberOfVertices(); - - /*Keep old thickness for later*/ - IssmDouble* oldthickness = xNew(numvertices); - basalelement->GetInputListOnNodes(&oldthickness[0],ThicknessEnum); - /*Fetch dof list and allocate solution vector*/ int *doflist = NULL; element->GetDofListLocal(&doflist,NoneApproximationEnum,GsetEnum); + + int numnodes = element->GetNumberOfNodes(); IssmDouble* newthickness = xNew(numnodes); /*Use the dof list to index into the solution vector: */ - IssmDouble minthickness = basalelement->FindParam(MasstransportMinThicknessEnum); + IssmDouble minthickness = element->FindParam(MasstransportMinThicknessEnum); for(int i=0;i(newthickness[i])) _error_("Inf found in solution vector"); if(newthickness[i]AddBasalInput(ThicknessEnum,newthickness,element->GetElementType()); - - /*Free ressources:*/ + xDelete(doflist); xDelete(newthickness); - xDelete(doflist); + /*Update bed and surface accordingly*/ + + /*Get basal element*/ + int domaintype; element->FindParam(&domaintype,DomainTypeEnum); + Element* basalelement=element; + if(domaintype!=Domain2DhorizontalEnum) basalelement = element->SpawnBasalElement(); + + /*Fetch number of nodes and dof for this finite element*/ + int numvertices = basalelement->GetNumberOfVertices(); + /*Now, we need to do some "processing"*/ newthickness = xNew(numvertices); + IssmDouble* oldthickness = xNew(numvertices); IssmDouble* cumdeltathickness = xNew(numvertices); IssmDouble* deltathickness = xNew(numvertices); IssmDouble* newbase = xNew(numvertices); @@ -855,6 +852,7 @@ /*Get previous base, thickness, surfac and current sealevel and bed:*/ basalelement->GetInputListOnVertices(&newthickness[0],ThicknessEnum); + basalelement->GetInputListOnVertices(&oldthickness[0],ThicknessOldEnum); basalelement->GetInputListOnVertices(&oldbase[0],BaseEnum); basalelement->GetInputListOnVertices(&oldsurface[0],SurfaceEnum); basalelement->GetInputListOnVertices(&phi[0],MaskGroundediceLevelsetEnum); @@ -868,8 +866,8 @@ /*What is the delta thickness forcing the sea-level rise core: cumulated over time, hence the +=:*/ for(int i=0;iFindParam(MaterialsRhoSeawaterEnum); for(int i=0;i0.){ //this is an ice sheet: just add thickness to base. - /*Update! actually, the bed has moved too!:*/ + if (phi[i]>0.){ //this is grounded ice: just add thickness to base. if(isgroundingline){ newsurface[i] = bed[i]+newthickness[i]; //surface = bed + newthickness newbase[i] = bed[i]; //new base at new bed @@ -887,7 +884,7 @@ } else{ newsurface[i] = oldbase[i]+newthickness[i]; //surface = oldbase + newthickness - newbase[i] = oldbase[i]; //same base: do nothing + newbase[i] = oldbase[i]; //same base: do nothing } } else{ //this is an ice shelf: hydrostatic equilibrium*/ Index: ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 24096) +++ ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 24097) @@ -776,6 +776,7 @@ ThicknessAcrossGradientEnum, ThicknessAlongGradientEnum, ThicknessEnum, + ThicknessOldEnum, ThicknessPositiveEnum, VelEnum, VxAverageEnum, Index: ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 24096) +++ ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 24097) @@ -782,6 +782,7 @@ case ThicknessAcrossGradientEnum : return "ThicknessAcrossGradient"; case ThicknessAlongGradientEnum : return "ThicknessAlongGradient"; case ThicknessEnum : return "Thickness"; + case ThicknessOldEnum : return "ThicknessOld"; case ThicknessPositiveEnum : return "ThicknessPositive"; case VelEnum : return "Vel"; case VxAverageEnum : return "VxAverage"; Index: ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 24096) +++ ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 24097) @@ -800,6 +800,7 @@ else if (strcmp(name,"ThicknessAcrossGradient")==0) return ThicknessAcrossGradientEnum; else if (strcmp(name,"ThicknessAlongGradient")==0) return ThicknessAlongGradientEnum; else if (strcmp(name,"Thickness")==0) return ThicknessEnum; + else if (strcmp(name,"ThicknessOld")==0) return ThicknessOldEnum; else if (strcmp(name,"ThicknessPositive")==0) return ThicknessPositiveEnum; else if (strcmp(name,"Vel")==0) return VelEnum; else if (strcmp(name,"VxAverage")==0) return VxAverageEnum; @@ -873,11 +874,11 @@ else if (strcmp(name,"Outputdefinition54")==0) return Outputdefinition54Enum; else if (strcmp(name,"Outputdefinition55")==0) return Outputdefinition55Enum; else if (strcmp(name,"Outputdefinition56")==0) return Outputdefinition56Enum; - else if (strcmp(name,"Outputdefinition57")==0) return Outputdefinition57Enum; else stage=8; } if(stage==8){ - if (strcmp(name,"Outputdefinition58")==0) return Outputdefinition58Enum; + if (strcmp(name,"Outputdefinition57")==0) return Outputdefinition57Enum; + else if (strcmp(name,"Outputdefinition58")==0) return Outputdefinition58Enum; else if (strcmp(name,"Outputdefinition59")==0) return Outputdefinition59Enum; else if (strcmp(name,"Outputdefinition5")==0) return Outputdefinition5Enum; else if (strcmp(name,"Outputdefinition60")==0) return Outputdefinition60Enum; @@ -996,11 +997,11 @@ else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum; else if (strcmp(name,"DoubleInput")==0) return DoubleInputEnum; else if (strcmp(name,"DoubleMatArrayParam")==0) return DoubleMatArrayParamEnum; - else if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum; else stage=9; } if(stage==9){ - if (strcmp(name,"DoubleMatParam")==0) return DoubleMatParamEnum; + if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum; + else if (strcmp(name,"DoubleMatParam")==0) return DoubleMatParamEnum; else if (strcmp(name,"DoubleParam")==0) return DoubleParamEnum; else if (strcmp(name,"DoubleVecParam")==0) return DoubleVecParamEnum; else if (strcmp(name,"Element")==0) return ElementEnum; @@ -1119,11 +1120,11 @@ else if (strcmp(name,"Materials")==0) return MaterialsEnum; else if (strcmp(name,"Matestar")==0) return MatestarEnum; else if (strcmp(name,"Matice")==0) return MaticeEnum; - else if (strcmp(name,"Matlitho")==0) return MatlithoEnum; else stage=10; } if(stage==10){ - if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum; + if (strcmp(name,"Matlitho")==0) return MatlithoEnum; + else if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum; else if (strcmp(name,"MaxAbsVx")==0) return MaxAbsVxEnum; else if (strcmp(name,"MaxAbsVy")==0) return MaxAbsVyEnum; else if (strcmp(name,"MaxAbsVz")==0) return MaxAbsVzEnum; @@ -1242,11 +1243,11 @@ else if (strcmp(name,"StringArrayParam")==0) return StringArrayParamEnum; else if (strcmp(name,"StringExternalResult")==0) return StringExternalResultEnum; else if (strcmp(name,"StringParam")==0) return StringParamEnum; - else if (strcmp(name,"SubelementFriction1")==0) return SubelementFriction1Enum; else stage=11; } if(stage==11){ - if (strcmp(name,"SubelementFriction2")==0) return SubelementFriction2Enum; + if (strcmp(name,"SubelementFriction1")==0) return SubelementFriction1Enum; + else if (strcmp(name,"SubelementFriction2")==0) return SubelementFriction2Enum; else if (strcmp(name,"SubelementMelt1")==0) return SubelementMelt1Enum; else if (strcmp(name,"SubelementMelt2")==0) return SubelementMelt2Enum; else if (strcmp(name,"SubelementMigration")==0) return SubelementMigrationEnum; Index: ../trunk-jpl/src/c/shared/Enum/Enum.vim =================================================================== --- ../trunk-jpl/src/c/shared/Enum/Enum.vim (revision 24096) +++ ../trunk-jpl/src/c/shared/Enum/Enum.vim (revision 24097) @@ -780,6 +780,7 @@ syn keyword cConstant ThicknessAcrossGradientEnum syn keyword cConstant ThicknessAlongGradientEnum syn keyword cConstant ThicknessEnum +syn keyword cConstant ThicknessOldEnum syn keyword cConstant ThicknessPositiveEnum syn keyword cConstant VelEnum syn keyword cConstant VxAverageEnum @@ -1280,6 +1281,7 @@ syn keyword cType Cfsurfacelogvel syn keyword cType Cfsurfacesquare syn keyword cType Channel +syn keyword cType classes syn keyword cType Constraint syn keyword cType Constraints syn keyword cType Contour @@ -1286,8 +1288,8 @@ syn keyword cType Contours syn keyword cType ControlInput syn keyword cType Covertree +syn keyword cType DatasetInput syn keyword cType DataSetParam -syn keyword cType DatasetInput syn keyword cType Definition syn keyword cType DependentObject syn keyword cType DoubleArrayInput @@ -1300,8 +1302,8 @@ syn keyword cType Element syn keyword cType ElementHook syn keyword cType ElementMatrix +syn keyword cType Elements syn keyword cType ElementVector -syn keyword cType Elements syn keyword cType ExponentialVariogram syn keyword cType ExternalResult syn keyword cType FemModel @@ -1308,11 +1310,12 @@ syn keyword cType FileParam syn keyword cType Friction syn keyword cType Gauss +syn keyword cType GaussianVariogram +syn keyword cType gaussobjects syn keyword cType GaussPenta syn keyword cType GaussSeg syn keyword cType GaussTetra syn keyword cType GaussTria -syn keyword cType GaussianVariogram syn keyword cType GenericExternalResult syn keyword cType GenericOption syn keyword cType GenericParam @@ -1320,6 +1323,7 @@ syn keyword cType Hook syn keyword cType Input syn keyword cType Inputs +syn keyword cType IntArrayInput syn keyword cType IntInput syn keyword cType IntMatParam syn keyword cType IntParam @@ -1327,6 +1331,7 @@ syn keyword cType IoModel syn keyword cType IssmDirectApplicInterface syn keyword cType IssmParallelDirectApplicInterface +syn keyword cType krigingobjects syn keyword cType Load syn keyword cType Loads syn keyword cType Masscon @@ -1337,6 +1342,7 @@ syn keyword cType Matestar syn keyword cType Matice syn keyword cType Matlitho +syn keyword cType matrixobjects syn keyword cType MatrixParam syn keyword cType Misfit syn keyword cType Moulin @@ -1349,8 +1355,8 @@ syn keyword cType Observation syn keyword cType Observations syn keyword cType Option +syn keyword cType Options syn keyword cType OptionUtilities -syn keyword cType Options syn keyword cType Param syn keyword cType Parameters syn keyword cType Pengrid @@ -1363,12 +1369,12 @@ syn keyword cType Quadtree syn keyword cType Regionaloutput syn keyword cType Results +syn keyword cType Riftfront syn keyword cType RiftStruct -syn keyword cType Riftfront syn keyword cType Seg syn keyword cType SegInput +syn keyword cType Segment syn keyword cType SegRef -syn keyword cType Segment syn keyword cType SpcDynamic syn keyword cType SpcStatic syn keyword cType SpcTransient @@ -1388,10 +1394,6 @@ syn keyword cType VectorParam syn keyword cType Vertex syn keyword cType Vertices -syn keyword cType classes -syn keyword cType gaussobjects -syn keyword cType krigingobjects -syn keyword cType matrixobjects syn keyword cType AdjointBalancethickness2Analysis syn keyword cType AdjointBalancethicknessAnalysis syn keyword cType AdjointHorizAnalysis @@ -1410,8 +1412,8 @@ syn keyword cType ExtrudeFromTopAnalysis syn keyword cType FreeSurfaceBaseAnalysis syn keyword cType FreeSurfaceTopAnalysis +syn keyword cType GiaIvinsAnalysis syn keyword cType GLheightadvectionAnalysis -syn keyword cType GiaIvinsAnalysis syn keyword cType HydrologyDCEfficientAnalysis syn keyword cType HydrologyDCInefficientAnalysis syn keyword cType HydrologyGlaDSAnalysis