Changeset 26423

Timestamp:

09/01/21 18:34:26 (4 years ago)

Author:

Mathieu Morlighem

Message:

BUG: minor

Location:

issm/trunk-jpl/src/c/analyses

Files:

• LevelsetAnalysis.cpp (modified) (2 diffs)
• StressbalanceAnalysis.cpp (modified) (2 diffs)

issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp

@@ -221 +221 @@
         basalelement->FindParam(&migrationmax,MigrationMaxEnum);
 
-        h = element->CharacteristicLength();
+        h = basalelement->CharacteristicLength();
 
         Input* mf_vx_input        = NULL;
@@ -388 +388 @@
         Input* mf_vy_input    = basalelement->GetInput(MovingFrontalVyEnum); _assert_(mf_vy_input);
 
-        IssmDouble h=element->CharacteristicLength();
+        IssmDouble h=basalelement->CharacteristicLength();
 
         /* Start  looping on the number of gaussian points: */

issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

@@ -1878 +1878 @@
         /*Start looping on Gaussian points*/
         Gauss* gauss=element->NewGauss(xyz_list,xyz_list_front,3);
+    
+    double time;
+    element->FindParam(&time,TimeEnum);
+    
         while(gauss->next()){
                 thickness_input->GetInputValue(&thickness,gauss);
@@ -1890 +1894 @@
                 ice_pressure   = 1.0/2.0*gravity*rho_ice*thickness*thickness;
                 pressure = ice_pressure + water_pressure;
+        
+        /*double alpha = 0.1/2*(1+sin(2*pi*time/(365*24*3600)))'*/
+        /*double alpha = 0.01/2*(1+sin(2*pi*time/(365*24*3600)))'*/
+        /*double alpha = 0.05/2*(1+sin(2*pi*time/(365*24*3600)))'*/
+        
+        if(base<sealevel){
+        double alpha = 0*0.02/2*(1+sin(2*PI*time/(365*24*3600)));
+                  alpha = 0.1;
+        pressure = ice_pressure -( 0.5*rho_ice*gravity*pow(alpha*thickness,2) + 0.5*rho_water*gravity*base*base-0.5*rho_ice*rho_ice/rho_water*gravity*pow(alpha*thickness,2));
+                  //pressure = ice_pressure + 0.5*rho_ice*gravity*pow(alpha*thickness,2) + 0.5*rho_water*gravity*(surface_under_water*surface_under_water - base_under_water*base_under_water)-0.5*rho_ice*rho_ice/rho_water*gravity*pow(alpha*thickness,2);
+                  }
 
                 for(int i=0;i<numnodes;i++){