Changeset 26133

Timestamp:

03/23/21 12:18:31 (4 years ago)

Author:

tsantos

Message:

CHG: working on mono layer HO

File:

• issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp (modified) (3 diffs)

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

@@ -2886 +2886 @@
 
         /*Intermediaries */
-        IssmDouble  thickness,Jdet,slope[2];
+        IssmDouble  thickness,Jdet,slope[2],n;
         IssmDouble* xyz_list = NULL;
 
@@ -2900 +2900 @@
         Input*     thickness_input=element->GetInput(ThicknessEnum); _assert_(thickness_input);
         Input*     surface_input  =element->GetInput(SurfaceEnum);   _assert_(surface_input);
+        Input*     n_input        =element->GetInput(MaterialsRheologyNEnum); _assert_(n_input);
         IssmDouble rhog = element->FindParam(MaterialsRhoIceEnum)*element->FindParam(ConstantsGEnum);
 
@@ -2910 +2911 @@
                 thickness_input->GetInputValue(&thickness,gauss);
                 surface_input->GetInputDerivativeValue(&slope[0],xyz_list,gauss);
+                n_input->GetInputValue(&n,gauss);
 
                 for(int i=0;i<numnodes;i++){
-                        pe->values[i*2+0]+=-rhog*thickness*slope[0]*Jdet*gauss->weight*basis[i];
-                        pe->values[i*2+1]+=-rhog*thickness*slope[1]*Jdet*gauss->weight*basis[i];
+                        pe->values[i+0]+=-rhog*thickness*slope[0]*Jdet*gauss->weight*basis[i]; //F1
+                        pe->values[i+3]+=-rhog*thickness*slope[0]*Jdet*gauss->weight*basis[i]*(n+1)/(n+2); //F2
+                        pe->values[i+6]+=-rhog*thickness*slope[1]*Jdet*gauss->weight*basis[i]; //F3
+                        pe->values[i+9]+=-rhog*thickness*slope[1]*Jdet*gauss->weight*basis[i]*(n+1)/(n+2); //F4
                 }
         }