Changeset 7408

Timestamp:

02/09/11 16:54:26 (14 years ago)

Author:

seroussi

Message:

added shelf_dampening in penta

Location:

issm/trunk/src/c

Files:

• modules/ModelProcessorx/CreateParameters.cpp (modified) (1 diff)
• objects/Elements/Penta.cpp (modified) (4 diffs)
• objects/IoModel.cpp (modified) (2 diffs)
• objects/IoModel.h (modified) (1 diff)

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

@@ -72 +72 @@
         parameters->AddObject(new IntParam(StabilizeConstraintsEnum,iomodel->stabilize_constraints));
         parameters->AddObject(new DoubleParam(StokesReconditioningEnum,iomodel->stokesreconditioning));
+        parameters->AddObject(new IntParam(ShelfDampeningEnum,iomodel->shelf_dampening));
         parameters->AddObject(new DoubleParam(ViscosityOvershootEnum,iomodel->viscosity_overshoot));
         parameters->AddObject(new BoolParam(WaitOnLockEnum,iomodel->waitonlock));

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

@@ -2625 +2625 @@
         /*Intermediaries*/
         int         i,j,ig;
-        int         approximation;
+        int         approximation,shelf_dampening;
         double      gravity,rho_water,bed,water_pressure;
+        double      damper,normal_vel,vx,vy,vz,dt;
         double          xyz_list_tria[NUMVERTICES2D][3];
         double      xyz_list[NUMVERTICES][3];
         double          bed_normal[3];
+        double      dz[3];
         double      l1l6[6]; //for the six nodes of the penta
         double      Jdet2d;
@@ -2637 +2639 @@
         if(!IsOnBed() || !IsOnShelf()) return NULL;
         inputs->GetParameterValue(&approximation,ApproximationEnum);
+        this->parameters->FindParam(&shelf_dampening,ShelfDampeningEnum);
         if(approximation!=StokesApproximationEnum && approximation!=MacAyealStokesApproximationEnum && approximation!=PattynStokesApproximationEnum) return NULL;
         ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters,StokesApproximationEnum);
@@ -2645 +2648 @@
         GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
         Input* bed_input=inputs->GetInput(BedEnum); _assert_(bed_input);
+        Input* vx_input=inputs->GetInput(VxEnum);   _assert_(vx_input);
+        Input* vy_input=inputs->GetInput(VyEnum);   _assert_(vy_input);
+        Input* vz_input=inputs->GetInput(VzEnum);   _assert_(vz_input);
 
         for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) xyz_list_tria[i][j]=xyz_list[i][j];
@@ -2657 +2663 @@
                 GetNodalFunctionsP1(l1l6, gauss);
 
+                BedNormal(&bed_normal[0],xyz_list_tria);
                 bed_input->GetParameterValue(&bed, gauss);
-                BedNormal(&bed_normal[0],xyz_list_tria);
+                if(shelf_dampening){ //add dampening to avoid too high vertical velocities when not in hydrostatic equilibrium
+                        bed_input->GetParameterDerivativeValue(&dz[0],&xyz_list[0][0],gauss);
+                        vx_input->GetParameterValue(&vx, gauss);
+                        vy_input->GetParameterValue(&vy, gauss);
+                        vz_input->GetParameterValue(&vz, gauss);
+                        dt=0;
+                        normal_vel=bed_normal[0]*vx+bed_normal[1]*vy+bed_normal[2]*vz;
+                        damper=gravity*rho_water*pow(1+pow(dz[0],2)+pow(dz[1],2),0.5)*normal_vel*dt;
+                }
+                else damper=0;
+                vz_input->GetParameterValue(&vz, gauss);
                 water_pressure=gravity*rho_water*bed;
 
-                for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) pe->values[i*NDOF4+j]+=water_pressure*gauss->weight*Jdet2d*l1l6[i]*bed_normal[j];
+                for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) pe->values[i*NDOF4+j]+=(water_pressure+damper)*gauss->weight*Jdet2d*l1l6[i]*bed_normal[j];
         }
 

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

@@ -191 +191 @@
         IoModelFetchData(&this->prognostic_DG,iomodel_handle,"prognostic_DG");
         IoModelFetchData(&this->stokesreconditioning,iomodel_handle,"stokesreconditioning");
+        IoModelFetchData(&this->shelf_dampening,iomodel_handle,"shelf_dampening");
         IoModelFetchData(&this->waitonlock,iomodel_handle,"waitonlock");
         IoModelFetchData(&this->kff,iomodel_handle,"kff");
@@ -353 +354 @@
         this->prognostic_DG=0;
         this->stokesreconditioning=0;
+        this->shelf_dampening=0;
         this->waitonlock=0;
 

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

@@ -119 +119 @@
                 double  viscosity_overshoot;
                 double  stokesreconditioning;
+                int     shelf_dampening;
                 double  cm_noisedmp;
                 double* cm_min;