Changeset 27271

Timestamp:

09/08/22 13:36:36 (3 years ago)

Author:

jwilner

Message:

CHG: use CalvingFromRate to determine the direction of calving so that Levermann is consistent with other calving laws

Location:

issm/trunk-jpl/src/c/classes/Elements

Files:

• Penta.cpp (modified) (5 diffs)
• Tria.cpp (modified) (5 diffs)

issm/trunk-jpl/src/c/classes/Elements/Penta.cpp

@@ -334 +334 @@
 void       Penta::CalvingRateLevermann(){/*{{{*/
 
-        IssmDouble  vx,vy,vel;
         IssmDouble  strainparallel;
         IssmDouble  propcoeff;
         IssmDouble  strainperpendicular;
-        IssmDouble  calvingratex[NUMVERTICES];
-        IssmDouble  calvingratey[NUMVERTICES];
         IssmDouble  calvingrate[NUMVERTICES];
 
@@ -355 +352 @@
 
                 /* Get the value we need*/
-                vx_input->GetInputValue(&vx,&gauss);
-                vy_input->GetInputValue(&vy,&gauss);
-                vel=vx*vx+vy*vy;
                 strainparallel_input->GetInputValue(&strainparallel,&gauss);
                 strainperpendicular_input->GetInputValue(&strainperpendicular,&gauss);
@@ -367 +361 @@
                         calvingrate[iv]=0;
                 }
-                calvingratex[iv]=calvingrate[iv]*vx/(sqrt(vel)+1.e-14);
-                calvingratey[iv]=calvingrate[iv]*vy/(sqrt(vel)+1.e-14);
         }
 
         /*Add input*/
-        this->AddBasalInput(CalvingratexEnum,&calvingratex[0],P1DGEnum);
-        this->AddBasalInput(CalvingrateyEnum,&calvingratey[0],P1DGEnum);
         this->AddBasalInput(CalvingCalvingrateEnum,&calvingrate[0],P1DGEnum);
+        this->CalvingFromRate();
 }/*}}}*/
 void       Penta::CalvingFluxLevelset(){/*{{{*/
@@ -2932 +2923 @@
                 /*Get calving speed*/
                 switch(calvinglaw){
+                        /*"Contiuous" calving*/
                         case DefaultCalvingEnum:
                         case CalvingVonmisesEnum:
@@ -2943 +2935 @@
                                 break;
 
+                        /*"Discrete" calving*/
                         case CalvingMinthicknessEnum:
                         case CalvingHabEnum:

issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

@@ -526 +526 @@
 void       Tria::CalvingRateLevermann(){/*{{{*/
 
-        IssmDouble  vx,vy,vel;
         IssmDouble  strainparallel;
         IssmDouble  propcoeff,bed;
         IssmDouble  strainperpendicular;
-        IssmDouble  calvingratex[NUMVERTICES];
-        IssmDouble  calvingratey[NUMVERTICES];
         IssmDouble  calvingrate[NUMVERTICES];
 
         /*Retrieve all inputs and parameters we will need*/
-        Input *vx_input                  = this->GetInput(VxEnum);                      _assert_(vx_input);
-        Input *vy_input                  = this->GetInput(VyEnum);                      _assert_(vy_input);
         Input *bs_input                  = this->GetInput(BaseEnum);                    _assert_(bs_input);
         Input *strainparallel_input      = this->GetInput(StrainRateparallelEnum);      _assert_(strainparallel_input);
@@ -548 +543 @@
 
                 /* Get the value we need*/
-                vx_input->GetInputValue(&vx,&gauss);
-                vy_input->GetInputValue(&vy,&gauss);
-                vel=vx*vx+vy*vy;
                 strainparallel_input->GetInputValue(&strainparallel,&gauss);
                 strainperpendicular_input->GetInputValue(&strainperpendicular,&gauss);
@@ -562 +554 @@
                 else
                         calvingrate[iv]=0.;
-
-                calvingratex[iv]=calvingrate[iv]*vx/(sqrt(vel)+1.e-14);
-                calvingratey[iv]=calvingrate[iv]*vy/(sqrt(vel)+1.e-14);
         }
 
         /*Add input*/
-        this->AddInput(CalvingratexEnum,&calvingratex[0],P1DGEnum);
-        this->AddInput(CalvingrateyEnum,&calvingratey[0],P1DGEnum);
         this->AddInput(CalvingCalvingrateEnum,&calvingrate[0],P1DGEnum);
+        this->CalvingFromRate();
 }/*}}}*/
 void       Tria::CalvingFluxLevelset(){/*{{{*/
@@ -4358 +4346 @@
                 /*Get calving speed*/
                 switch(calvinglaw){
+
+                        /*RATE calving laws*/
                         case DefaultCalvingEnum:
                         case CalvingVonmisesEnum:
                         case CalvingTestEnum:
                         case CalvingParameterizationEnum:
+                        case CalvingLevermannEnum:
                                 calvingratex_input->GetInputValue(&c[0],&gauss);
                                 calvingratey_input->GetInputValue(&c[1],&gauss);
@@ -4376 +4367 @@
                                 break;
 
-                        case CalvingLevermannEnum:
-                                calvingratex_input->GetInputValue(&c[0],&gauss);
-                                if(dim==2) calvingratey_input->GetInputValue(&c[1],&gauss);
-                                meltingrate_input->GetInputValue(&meltingrate,&gauss);
-                                norm_calving=0.;
-                                for(i=0;i<dim;i++) norm_calving+=pow(c[i],2);
-                                norm_calving=sqrt(norm_calving)+1.e-14;
-                                for(i=0;i<dim;i++) m[i]=meltingrate*c[i]/norm_calving;
-                                break;
-
+                        /*Discrete calving laws*/
                         case CalvingMinthicknessEnum:
                         case CalvingHabEnum: