Changeset 27673

Timestamp:

04/03/23 11:12:35 (2 years ago)

Author:

Mathieu Morlighem

Message:

CHG: added inversions for Coulomb2

Location:

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

Files:

• Friction.cpp (modified) (2 diffs)
• Friction.h (modified) (1 diff)

issm/trunk-jpl/src/c/classes/Loads/Friction.cpp

@@ -141 +141 @@
                                 GetAlphaSchoofComplement(palpha_complement,gauss);
                                 break;
+                        case 13:
+                                GetAlphaCoulomb2Complement(palpha_complement,gauss);
+                                break;
                         case 14:
                                 GetAlphaRegCoulombComplement(palpha_complement,gauss);
@@ -354 +357 @@
         /*Assign output pointers:*/
         *palpha_complement=alpha_complement;
+}/*}}}*/
+void Friction::GetAlphaCoulomb2Complement(IssmDouble* palpha_complement, Gauss* gauss){/*{{{*/
+
+        /* Compute the complement of Cornford's friction law for inversion
+         * d alpha2                       
+         * ------ = (C*N*v^m)/(C^(2/m)*v + (N/2)^(1/m))^m - (C^(2/m - 1)*C^2*N*v*v^m)/(C^(2/m)*v + (N/2)^(1/m))^(m + 1)
+         *  dC                           
+         */
+
+        /*diverse: */
+        IssmDouble  m, C;
+        IssmDouble  alpha_complement;
+
+        /*Recover parameters: */
+        element->GetInputValue(&C,gauss,FrictionCEnum);
+        element->GetInputValue(&m,gauss,FrictionMEnum);
+
+        /*Get effective pressure and velocity magnitude*/
+        IssmDouble N = EffectivePressure(gauss);
+        IssmDouble v = VelMag(gauss);
+
+        /*Compute alpha*/
+        if(v<1e-10){
+                alpha_complement = 0.;
+        }
+        else{
+                alpha_complement= pow(0.5*N,1./m+1)* pow(v,m-1.) * pow(v*pow(C,1./m) +pow(0.5*N,1./m) ,-m-1.);
+        }
+
+        /*Assign output pointers:*/
+        *palpha_complement=alpha_complement/2.;
 }/*}}}*/
 void Friction::GetAlpha2(IssmDouble* palpha2, Gauss* gauss){/*{{{*/

issm/trunk-jpl/src/c/classes/Loads/Friction.h

@@ -41 +41 @@
                 void  GetAlphaViscousComplement(IssmDouble* alpha_complement,Gauss* gauss);
                 void  GetAlphaSchoofComplement(IssmDouble* alpha_complement,Gauss* gauss);
+                void  GetAlphaCoulomb2Complement(IssmDouble* alpha_complement,Gauss* gauss);
                 void  GetAlphaRegCoulombComplement(IssmDouble* alpha_complement,Gauss* gauss);
                 void  GetAlphaWeertmanComplement(IssmDouble* alpha_complement,Gauss* gauss);