source: issm/oecreview/Archive/23390-24306/ISSM-24096-24097.diff@ 24307

Last change on this file since 24307 was 24307, checked in by Mathieu Morlighem, 5 years ago

NEW: adding Archive/23390-24306

File size: 14.5 KB
  • ../trunk-jpl/src/c/cores/masstransport_core.cpp

     
    6060                        depthaverage_core(femmodel);
    6161                }
    6262                femmodel->SetCurrentConfiguration(MasstransportAnalysisEnum);
     63                InputDuplicatex(femmodel,ThicknessEnum,ThicknessOldEnum);
    6364                if(stabilization==4){
    6465                        solutionsequence_fct(femmodel);
    6566                }
  • ../trunk-jpl/src/c/classes/Elements/Penta.cpp

     
    18071807
    18081808        /*Fetch dof list and allocate solution vector*/
    18091809        GetDofListLocal(&doflist,NoneApproximationEnum,GsetEnum);
    1810         IssmDouble* values    = xNew<IssmDouble>(numnodes);
     1810        IssmDouble* values = xNew<IssmDouble>(numnodes);
    18111811
    18121812        /*Use the dof list to index into the solution vector: */
    18131813        for(int i=0;i<numnodes;i++){
  • ../trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp

     
    807807        /*Only update if on base*/
    808808        if(!element->IsOnBase()) return;
    809809
    810         /*Get basal element*/
    811         int domaintype; element->FindParam(&domaintype,DomainTypeEnum);
    812    Element* basalelement=element;
    813         if(domaintype!=Domain2DhorizontalEnum) basalelement = element->SpawnBasalElement();
    814 
    815         /*Fetch number of nodes and dof for this finite element*/
    816         int numnodes    = basalelement->GetNumberOfNodes();
    817    int numvertices = basalelement->GetNumberOfVertices();
    818 
    819    /*Keep old thickness for later*/
    820         IssmDouble* oldthickness = xNew<IssmDouble>(numvertices);
    821    basalelement->GetInputListOnNodes(&oldthickness[0],ThicknessEnum);
    822 
    823810        /*Fetch dof list and allocate solution vector*/
    824811        int *doflist = NULL;
    825812        element->GetDofListLocal(&doflist,NoneApproximationEnum,GsetEnum);
     813
     814        int numnodes = element->GetNumberOfNodes();
    826815        IssmDouble* newthickness = xNew<IssmDouble>(numnodes);
    827816
    828817        /*Use the dof list to index into the solution vector: */
    829    IssmDouble minthickness = basalelement->FindParam(MasstransportMinThicknessEnum);
     818   IssmDouble minthickness = element->FindParam(MasstransportMinThicknessEnum);
    830819        for(int i=0;i<numnodes;i++){
    831820                newthickness[i]=solution[doflist[i]];
    832821                /*Check solution*/
     
    834823                if(xIsInf<IssmDouble>(newthickness[i])) _error_("Inf found in solution vector");
    835824      if(newthickness[i]<minthickness) newthickness[i]=minthickness;
    836825        }
    837 
    838826        element->AddBasalInput(ThicknessEnum,newthickness,element->GetElementType());
    839 
    840         /*Free ressources:*/
     827        xDelete<int>(doflist);
    841828        xDelete<IssmDouble>(newthickness);
    842         xDelete<int>(doflist);
    843829
     830        /*Update bed and surface accordingly*/
     831
     832        /*Get basal element*/
     833        int domaintype; element->FindParam(&domaintype,DomainTypeEnum);
     834   Element* basalelement=element;
     835        if(domaintype!=Domain2DhorizontalEnum) basalelement = element->SpawnBasalElement();
     836
     837        /*Fetch number of nodes and dof for this finite element*/
     838   int numvertices = basalelement->GetNumberOfVertices();
     839
    844840        /*Now, we need to do some "processing"*/
    845841   newthickness  = xNew<IssmDouble>(numvertices);
     842        IssmDouble* oldthickness      = xNew<IssmDouble>(numvertices);
    846843        IssmDouble* cumdeltathickness = xNew<IssmDouble>(numvertices);
    847844        IssmDouble* deltathickness    = xNew<IssmDouble>(numvertices);
    848845        IssmDouble* newbase           = xNew<IssmDouble>(numvertices);
     
    855852
    856853        /*Get previous base, thickness, surfac and current sealevel and bed:*/
    857854   basalelement->GetInputListOnVertices(&newthickness[0],ThicknessEnum);
     855        basalelement->GetInputListOnVertices(&oldthickness[0],ThicknessOldEnum);
    858856        basalelement->GetInputListOnVertices(&oldbase[0],BaseEnum);
    859857        basalelement->GetInputListOnVertices(&oldsurface[0],SurfaceEnum);
    860858        basalelement->GetInputListOnVertices(&phi[0],MaskGroundediceLevelsetEnum);
     
    868866
    869867        /*What is the delta thickness forcing the sea-level rise core: cumulated over time, hence the +=:*/
    870868        for(int i=0;i<numvertices;i++){
    871                 cumdeltathickness[i]+=newthickness[i]-oldthickness[i];
    872                 deltathickness[i]=newthickness[i]-oldthickness[i];
     869                cumdeltathickness[i] += newthickness[i]-oldthickness[i];
     870                deltathickness[i]     = newthickness[i]-oldthickness[i];
    873871        }
    874872
    875873        /*Find MasstransportHydrostaticAdjustment to figure out how to update the geometry:*/
     
    879877        IssmDouble rho_water = basalelement->FindParam(MaterialsRhoSeawaterEnum);
    880878
    881879        for(int i=0;i<numvertices;i++) {
    882                 if (phi[i]>0.){ //this is an ice sheet: just add thickness to base.
    883                         /*Update! actually, the bed has moved too!:*/
     880                if (phi[i]>0.){ //this is grounded ice: just add thickness to base.
    884881                        if(isgroundingline){
    885882                                newsurface[i] = bed[i]+newthickness[i]; //surface = bed + newthickness
    886883                                newbase[i]    = bed[i];                 //new base at new bed
     
    887884                        }
    888885                        else{
    889886                                 newsurface[i] = oldbase[i]+newthickness[i]; //surface = oldbase + newthickness
    890                                  newbase[i]     = oldbase[i];                 //same base: do nothing
     887                                 newbase[i]    = oldbase[i];                 //same base: do nothing
    891888                        }
    892889                }
    893890                else{ //this is an ice shelf: hydrostatic equilibrium*/
  • ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

     
    776776        ThicknessAcrossGradientEnum,
    777777        ThicknessAlongGradientEnum,
    778778        ThicknessEnum,
     779        ThicknessOldEnum,
    779780        ThicknessPositiveEnum,
    780781        VelEnum,
    781782        VxAverageEnum,
  • ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

     
    782782                case ThicknessAcrossGradientEnum : return "ThicknessAcrossGradient";
    783783                case ThicknessAlongGradientEnum : return "ThicknessAlongGradient";
    784784                case ThicknessEnum : return "Thickness";
     785                case ThicknessOldEnum : return "ThicknessOld";
    785786                case ThicknessPositiveEnum : return "ThicknessPositive";
    786787                case VelEnum : return "Vel";
    787788                case VxAverageEnum : return "VxAverage";
  • ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

     
    800800              else if (strcmp(name,"ThicknessAcrossGradient")==0) return ThicknessAcrossGradientEnum;
    801801              else if (strcmp(name,"ThicknessAlongGradient")==0) return ThicknessAlongGradientEnum;
    802802              else if (strcmp(name,"Thickness")==0) return ThicknessEnum;
     803              else if (strcmp(name,"ThicknessOld")==0) return ThicknessOldEnum;
    803804              else if (strcmp(name,"ThicknessPositive")==0) return ThicknessPositiveEnum;
    804805              else if (strcmp(name,"Vel")==0) return VelEnum;
    805806              else if (strcmp(name,"VxAverage")==0) return VxAverageEnum;
     
    873874              else if (strcmp(name,"Outputdefinition54")==0) return Outputdefinition54Enum;
    874875              else if (strcmp(name,"Outputdefinition55")==0) return Outputdefinition55Enum;
    875876              else if (strcmp(name,"Outputdefinition56")==0) return Outputdefinition56Enum;
    876               else if (strcmp(name,"Outputdefinition57")==0) return Outputdefinition57Enum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"Outputdefinition58")==0) return Outputdefinition58Enum;
     880              if (strcmp(name,"Outputdefinition57")==0) return Outputdefinition57Enum;
     881              else if (strcmp(name,"Outputdefinition58")==0) return Outputdefinition58Enum;
    881882              else if (strcmp(name,"Outputdefinition59")==0) return Outputdefinition59Enum;
    882883              else if (strcmp(name,"Outputdefinition5")==0) return Outputdefinition5Enum;
    883884              else if (strcmp(name,"Outputdefinition60")==0) return Outputdefinition60Enum;
     
    996997              else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
    997998              else if (strcmp(name,"DoubleInput")==0) return DoubleInputEnum;
    998999              else if (strcmp(name,"DoubleMatArrayParam")==0) return DoubleMatArrayParamEnum;
    999               else if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
    10001000         else stage=9;
    10011001   }
    10021002   if(stage==9){
    1003               if (strcmp(name,"DoubleMatParam")==0) return DoubleMatParamEnum;
     1003              if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
     1004              else if (strcmp(name,"DoubleMatParam")==0) return DoubleMatParamEnum;
    10041005              else if (strcmp(name,"DoubleParam")==0) return DoubleParamEnum;
    10051006              else if (strcmp(name,"DoubleVecParam")==0) return DoubleVecParamEnum;
    10061007              else if (strcmp(name,"Element")==0) return ElementEnum;
     
    11191120              else if (strcmp(name,"Materials")==0) return MaterialsEnum;
    11201121              else if (strcmp(name,"Matestar")==0) return MatestarEnum;
    11211122              else if (strcmp(name,"Matice")==0) return MaticeEnum;
    1122               else if (strcmp(name,"Matlitho")==0) return MatlithoEnum;
    11231123         else stage=10;
    11241124   }
    11251125   if(stage==10){
    1126               if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum;
     1126              if (strcmp(name,"Matlitho")==0) return MatlithoEnum;
     1127              else if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum;
    11271128              else if (strcmp(name,"MaxAbsVx")==0) return MaxAbsVxEnum;
    11281129              else if (strcmp(name,"MaxAbsVy")==0) return MaxAbsVyEnum;
    11291130              else if (strcmp(name,"MaxAbsVz")==0) return MaxAbsVzEnum;
     
    12421243              else if (strcmp(name,"StringArrayParam")==0) return StringArrayParamEnum;
    12431244              else if (strcmp(name,"StringExternalResult")==0) return StringExternalResultEnum;
    12441245              else if (strcmp(name,"StringParam")==0) return StringParamEnum;
    1245               else if (strcmp(name,"SubelementFriction1")==0) return SubelementFriction1Enum;
    12461246         else stage=11;
    12471247   }
    12481248   if(stage==11){
    1249               if (strcmp(name,"SubelementFriction2")==0) return SubelementFriction2Enum;
     1249              if (strcmp(name,"SubelementFriction1")==0) return SubelementFriction1Enum;
     1250              else if (strcmp(name,"SubelementFriction2")==0) return SubelementFriction2Enum;
    12501251              else if (strcmp(name,"SubelementMelt1")==0) return SubelementMelt1Enum;
    12511252              else if (strcmp(name,"SubelementMelt2")==0) return SubelementMelt2Enum;
    12521253              else if (strcmp(name,"SubelementMigration")==0) return SubelementMigrationEnum;
  • ../trunk-jpl/src/c/shared/Enum/Enum.vim

     
    780780syn keyword cConstant ThicknessAcrossGradientEnum
    781781syn keyword cConstant ThicknessAlongGradientEnum
    782782syn keyword cConstant ThicknessEnum
     783syn keyword cConstant ThicknessOldEnum
    783784syn keyword cConstant ThicknessPositiveEnum
    784785syn keyword cConstant VelEnum
    785786syn keyword cConstant VxAverageEnum
     
    12801281syn keyword cType Cfsurfacelogvel
    12811282syn keyword cType Cfsurfacesquare
    12821283syn keyword cType Channel
     1284syn keyword cType classes
    12831285syn keyword cType Constraint
    12841286syn keyword cType Constraints
    12851287syn keyword cType Contour
     
    12861288syn keyword cType Contours
    12871289syn keyword cType ControlInput
    12881290syn keyword cType Covertree
     1291syn keyword cType DatasetInput
    12891292syn keyword cType DataSetParam
    1290 syn keyword cType DatasetInput
    12911293syn keyword cType Definition
    12921294syn keyword cType DependentObject
    12931295syn keyword cType DoubleArrayInput
     
    13001302syn keyword cType Element
    13011303syn keyword cType ElementHook
    13021304syn keyword cType ElementMatrix
     1305syn keyword cType Elements
    13031306syn keyword cType ElementVector
    1304 syn keyword cType Elements
    13051307syn keyword cType ExponentialVariogram
    13061308syn keyword cType ExternalResult
    13071309syn keyword cType FemModel
     
    13081310syn keyword cType FileParam
    13091311syn keyword cType Friction
    13101312syn keyword cType Gauss
     1313syn keyword cType GaussianVariogram
     1314syn keyword cType gaussobjects
    13111315syn keyword cType GaussPenta
    13121316syn keyword cType GaussSeg
    13131317syn keyword cType GaussTetra
    13141318syn keyword cType GaussTria
    1315 syn keyword cType GaussianVariogram
    13161319syn keyword cType GenericExternalResult
    13171320syn keyword cType GenericOption
    13181321syn keyword cType GenericParam
     
    13201323syn keyword cType Hook
    13211324syn keyword cType Input
    13221325syn keyword cType Inputs
     1326syn keyword cType IntArrayInput
    13231327syn keyword cType IntInput
    13241328syn keyword cType IntMatParam
    13251329syn keyword cType IntParam
     
    13271331syn keyword cType IoModel
    13281332syn keyword cType IssmDirectApplicInterface
    13291333syn keyword cType IssmParallelDirectApplicInterface
     1334syn keyword cType krigingobjects
    13301335syn keyword cType Load
    13311336syn keyword cType Loads
    13321337syn keyword cType Masscon
     
    13371342syn keyword cType Matestar
    13381343syn keyword cType Matice
    13391344syn keyword cType Matlitho
     1345syn keyword cType matrixobjects
    13401346syn keyword cType MatrixParam
    13411347syn keyword cType Misfit
    13421348syn keyword cType Moulin
     
    13491355syn keyword cType Observation
    13501356syn keyword cType Observations
    13511357syn keyword cType Option
     1358syn keyword cType Options
    13521359syn keyword cType OptionUtilities
    1353 syn keyword cType Options
    13541360syn keyword cType Param
    13551361syn keyword cType Parameters
    13561362syn keyword cType Pengrid
     
    13631369syn keyword cType Quadtree
    13641370syn keyword cType Regionaloutput
    13651371syn keyword cType Results
     1372syn keyword cType Riftfront
    13661373syn keyword cType RiftStruct
    1367 syn keyword cType Riftfront
    13681374syn keyword cType Seg
    13691375syn keyword cType SegInput
     1376syn keyword cType Segment
    13701377syn keyword cType SegRef
    1371 syn keyword cType Segment
    13721378syn keyword cType SpcDynamic
    13731379syn keyword cType SpcStatic
    13741380syn keyword cType SpcTransient
     
    13881394syn keyword cType VectorParam
    13891395syn keyword cType Vertex
    13901396syn keyword cType Vertices
    1391 syn keyword cType classes
    1392 syn keyword cType gaussobjects
    1393 syn keyword cType krigingobjects
    1394 syn keyword cType matrixobjects
    13951397syn keyword cType AdjointBalancethickness2Analysis
    13961398syn keyword cType AdjointBalancethicknessAnalysis
    13971399syn keyword cType AdjointHorizAnalysis
     
    14101412syn keyword cType ExtrudeFromTopAnalysis
    14111413syn keyword cType FreeSurfaceBaseAnalysis
    14121414syn keyword cType FreeSurfaceTopAnalysis
     1415syn keyword cType GiaIvinsAnalysis
    14131416syn keyword cType GLheightadvectionAnalysis
    1414 syn keyword cType GiaIvinsAnalysis
    14151417syn keyword cType HydrologyDCEfficientAnalysis
    14161418syn keyword cType HydrologyDCInefficientAnalysis
    14171419syn keyword cType HydrologyGlaDSAnalysis
Note: See TracBrowser for help on using the repository browser.