Changeset 6946

Timestamp:

01/05/11 11:28:39 (14 years ago)

Author:

seroussi

Message:

started to add hydrostatic adjustment

Location:

issm/trunk/src/c

Files:

• EnumDefinitions/EnumDefinitions.h (modified) (2 diffs)
• EnumDefinitions/EnumToString.cpp (modified) (2 diffs)
• EnumDefinitions/StringToEnum.cpp (modified) (2 diffs)
• modules/ModelProcessorx/CreateParameters.cpp (modified) (1 diff)
• objects/Elements/Penta.cpp (modified) (3 diffs)
• objects/Elements/Tria.cpp (modified) (3 diffs)
• objects/IoModel.cpp (modified) (2 diffs)
• objects/IoModel.h (modified) (1 diff)

issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

@@ -199 +199 @@
         FrictionEnum,
         GeothermalFluxEnum,
+        HydrostaticAdjustmentEnum,
         InternalEnum,
         KflagEnum,
@@ -310 +311 @@
         ResidualEnum,
         AbsoluteEnum,
+        IncrementalEnum,
         /*}}}*/
         /*Material{{{1*/

issm/trunk/src/c/EnumDefinitions/EnumToString.cpp

@@ -171 +171 @@
                 case FrictionEnum : return "Friction";
                 case GeothermalFluxEnum : return "GeothermalFlux";
+                case HydrostaticAdjustmentEnum : return "HydrostaticAdjustment";
                 case InternalEnum : return "Internal";
                 case KflagEnum : return "Kflag";
@@ -272 +273 @@
                 case ResidualEnum : return "Residual";
                 case AbsoluteEnum : return "Absolute";
+                case IncrementalEnum : return "Incremental";
                 case RhoIceEnum : return "RhoIce";
                 case RhoWaterEnum : return "RhoWater";

issm/trunk/src/c/EnumDefinitions/StringToEnum.cpp

@@ -169 +169 @@
         else if (strcmp(name,"Friction")==0) return FrictionEnum;
         else if (strcmp(name,"GeothermalFlux")==0) return GeothermalFluxEnum;
+        else if (strcmp(name,"HydrostaticAdjustment")==0) return HydrostaticAdjustmentEnum;
         else if (strcmp(name,"Internal")==0) return InternalEnum;
         else if (strcmp(name,"Kflag")==0) return KflagEnum;
@@ -270 +271 @@
         else if (strcmp(name,"Residual")==0) return ResidualEnum;
         else if (strcmp(name,"Absolute")==0) return AbsoluteEnum;
+        else if (strcmp(name,"Incremental")==0) return IncrementalEnum;
         else if (strcmp(name,"RhoIce")==0) return RhoIceEnum;
         else if (strcmp(name,"RhoWater")==0) return RhoWaterEnum;

issm/trunk/src/c/modules/ModelProcessorx/CreateParameters.cpp

@@ -55 +55 @@
         parameters->AddObject(new BoolParam(TimeAdaptEnum,iomodel->time_adapt)); 
         parameters->AddObject(new DoubleParam(CflCoefficientEnum,iomodel->cfl_coefficient)); 
+        parameters->AddObject(new IntParam(HydrostaticAdjustmentEnum,iomodel->hydrostatic_adjustment)); 
         parameters->AddObject(new DoubleParam(PenaltyOffsetEnum,iomodel->penalty_offset));
         parameters->AddObject(new DoubleParam(SparsityEnum,iomodel->sparsity));

issm/trunk/src/c/objects/Elements/Penta.cpp

@@ -5805 +5805 @@
 
         /*Intermediaries*/
-        int    i;
+        int hydroadjustment;
         double rho_ice,rho_water;
 
@@ -5813 +5813 @@
         /*If shelf: hydrostatic equilibrium*/
         if (this->IsOnShelf()){
+
+                /*Fing HydrostaticAdjustment to figure out how to update the geometry:*/
+                this->parameters->FindParam(&hydroadjustment,HydrostaticAdjustmentEnum);
 
                 /*recover material parameters: */
@@ -5818 +5821 @@
                 rho_water=matpar->GetRhoWater();
 
-                /*Create New Surface: s = (1-rho_ice/rho_water) h*/
-                this->inputs->DuplicateInput(ThicknessEnum,SurfaceEnum);     //1: copy thickness into surface
-                InputScale(SurfaceEnum,(1-rho_ice/rho_water)); //2: surface = surface * (1-di)
-
-                /*Create New Bed b = -rho_ice/rho_water h*/
-                this->inputs->DuplicateInput(ThicknessEnum,BedEnum);         //1: copy thickness into bed
-                InputScale(BedEnum, -rho_ice/rho_water);       //2: bed = bed * (-di)
+                if(hydroadjustment==AbsoluteEnum){
+                        /*Create New Surface: s = (1-rho_ice/rho_water) h*/
+                        this->inputs->DuplicateInput(ThicknessEnum,SurfaceEnum);     //1: copy thickness into surface
+                        InputScale(SurfaceEnum,(1-rho_ice/rho_water)); //2: surface = surface * (1-di)
+
+                        /*Create New Bed b = -rho_ice/rho_water h*/
+                        this->inputs->DuplicateInput(ThicknessEnum,BedEnum);         //1: copy thickness into bed
+                        InputScale(BedEnum, -rho_ice/rho_water);       //2: bed = bed * (-di)
+                }
+                else if(hydroadjustment==IncrementalEnum){
+                        /*The bed is update with di*Thickness*/
+                        this->inputs->AXPY(SurfaceEnum,1.0-rho_ice/rho_water,ThicknessEnum);     //surface = surface + (1-di) * thickness
+
+                        /*The surface is update with (1-di)*Thickness*/
+                        this->inputs->AXPY(BedEnum,rho_ice/rho_water,ThicknessEnum);     //bed = bde + di * thickness
+                }
+                else _error_("Hydrostatic adjustment %i (%s) not supported yet",hydroadjustment,EnumToString(hydroadjustment));
         }
 

issm/trunk/src/c/objects/Elements/Tria.cpp

@@ -5129 +5129 @@
 
         /*Intermediaries*/
+        int hydroadjustment;
         double rho_ice,rho_water;
 
@@ -5136 +5137 @@
         /*If shelf: hydrostatic equilibrium*/
         if (this->IsOnShelf()){
+
+                /*Fing HydrostaticAdjustment to figure out how to update the geometry:*/
+                this->parameters->FindParam(&hydroadjustment,HydrostaticAdjustmentEnum);
 
                 /*recover material parameters: */
@@ -5141 +5145 @@
                 rho_water=matpar->GetRhoWater();
 
-                /*Create New Surface: s = (1-rho_ice/rho_water) h*/
-                this->inputs->DuplicateInput(ThicknessEnum,SurfaceEnum);     //1: copy thickness into surface
-                InputScale(SurfaceEnum,(1.0-rho_ice/rho_water)); //2: surface = surface * (1-di)
-
-                /*Create New Bed b = -rho_ice/rho_water h*/
-                this->inputs->DuplicateInput(ThicknessEnum,BedEnum);         //1: copy thickness into bed
-                InputScale(BedEnum, -rho_ice/rho_water);       //2: bed = bed * (-di)
+                if(hydroadjustment==AbsoluteEnum){
+                        /*Create New Surface: s = (1-rho_ice/rho_water) h*/
+                        this->inputs->DuplicateInput(ThicknessEnum,SurfaceEnum);     //1: copy thickness into surface
+                        InputScale(SurfaceEnum,(1.0-rho_ice/rho_water)); //2: surface = surface * (1-di)
+
+                        /*Create New Bed b = -rho_ice/rho_water h*/
+                        this->inputs->DuplicateInput(ThicknessEnum,BedEnum);         //1: copy thickness into bed
+                        InputScale(BedEnum, -rho_ice/rho_water);       //2: bed = bed * (-di)
+                }
+                else if(hydroadjustment==IncrementalEnum){
+                        /*The bed is update with di*Thickness*/
+                        this->inputs->AXPY(SurfaceEnum,1.0-rho_ice/rho_water,ThicknessEnum);     //surface = surface + (1-di) * thickness
+
+                        /*The surface is update with (1-di)*Thickness*/
+                        this->inputs->AXPY(BedEnum,rho_ice/rho_water,ThicknessEnum);     //bed = bde + di * thickness
+                }
+                else _error_("Hydrostatic adjustment %i (%s) not supported yet",hydroadjustment,EnumToString(hydroadjustment));
         }
 

issm/trunk/src/c/objects/IoModel.cpp

@@ -178 +178 @@
         IoModelFetchData(&this->time_adapt,iomodel_handle,"time_adapt");
         IoModelFetchData(&this->cfl_coefficient,iomodel_handle,"cfl_coefficient");
+        IoModelFetchData(&this->hydrostatic_adjustment,iomodel_handle,"hydrostatic_adjustment");
         IoModelFetchData(&this->penalty_offset,iomodel_handle,"penalty_offset");
         IoModelFetchData(&this->penalty_melting,iomodel_handle,"penalty_melting");
@@ -326 +327 @@
         this->time_adapt=0;
         this->cfl_coefficient=0;
+        this->hydrostatic_adjustment=0;
         this->penalty_offset=0;
         this->penalty_melting=0;

issm/trunk/src/c/objects/IoModel.h

@@ -146 +146 @@
                 int     time_adapt;
                 double  cfl_coefficient;
+                int     hydrostatic_adjustment;
                 double  penalty_offset;
                 double  penalty_melting;