[cig-commits] r12994 - in long/3D/Gale/trunk: . src/StgFEM/plugins/StandardConditionFunctions

[walter at geodynamics.org]

Author: walter
Date: 2008-10-04 11:44:12 -0700 (Sat, 04 Oct 2008)
New Revision: 12994

Modified:
   long/3D/Gale/trunk/
   long/3D/Gale/trunk/src/StgFEM/plugins/StandardConditionFunctions/StandardConditionFunctions.c
Log:
 r2352 at earth:  boo | 2008-09-24 05:07:34 -0700
 Fix a sign error and a factor of two bug in frcition.  Does not seem to work.



Property changes on: long/3D/Gale/trunk
___________________________________________________________________
Name: svk:merge
   - 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2351
   + 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2352

Modified: long/3D/Gale/trunk/src/StgFEM/plugins/StandardConditionFunctions/StandardConditionFunctions.c
===================================================================
--- long/3D/Gale/trunk/src/StgFEM/plugins/StandardConditionFunctions/StandardConditionFunctions.c	2008-10-04 18:44:09 UTC (rev 12993)
+++ long/3D/Gale/trunk/src/StgFEM/plugins/StandardConditionFunctions/StandardConditionFunctions.c	2008-10-04 18:44:12 UTC (rev 12994)
@@ -1580,15 +1580,17 @@
   FeVariable_GetValueAtNode(velocityField,node_lI,v);
   FeVariable_GetValueAtNode(viscosityField,node_lI,&viscosity);
 
+  node_gI = Mesh_DomainToGlobal( feMesh, MT_VERTEX, node_lI );
+  RegularMeshUtils_Node_1DTo3D( feMesh, node_gI, inds );
+  
   /* Assume 2D, boundaries along x or y axis */
   normal_stress=stress[normal_dim];
+  if(inds[normal_dim]!=0)
+    normal_stress=-normal_stress;
+
   tangent_stress=stress[2];
-
   yield_stress=(normal_stress+p)*coefficient;
 
-  node_gI = Mesh_DomainToGlobal( feMesh, MT_VERTEX, node_lI );
-  RegularMeshUtils_Node_1DTo3D( feMesh, node_gI, inds );
-  
   grid = *(Grid**)ExtensionManager_Get( feMesh->info, feMesh, 
                                         ExtensionManager_GetHandle( feMesh->info, "vertexGrid" ) );
   
@@ -1616,9 +1618,11 @@
   else if(*(context->nonLinearIteration_I)==0
           || (fabs(tangent_stress)<yield_stress && v[tangent_dim]==boundary_v))
     {
-      *result=v[tangent_dim];
+      *result=boundary_v;
+/*       *result=v[tangent_dim]; */
 
-      printf("friction %d %d old\n",tangent_dim,inds[tangent_dim]);
+      printf("friction %d %d %g %g old\n",tangent_dim,inds[tangent_dim],
+             *result,tangent_stress);
     }      
   else
     {
@@ -1641,9 +1645,9 @@
       /* Make sure that we always make the velocity closer to the
          boundary velocity */
       if(v_off[tangent_dim]>boundary_v)
-        *result=v_off[tangent_dim]-delta*yield_stress/(2*viscosity);
+        *result=v_off[tangent_dim]-delta*yield_stress/viscosity;
       else
-        *result=v_off[tangent_dim]+delta*yield_stress/(2*viscosity);
+        *result=v_off[tangent_dim]+delta*yield_stress/viscosity;
       
       /* Average the results so that we don't oscillate a lot, looking
          for a solution. */
@@ -1662,12 +1666,18 @@
 
       /* If we overshot, make the node stick */
       if((*result-boundary_v)*(v_off[tangent_dim]-boundary_v)<0)
-        *result=boundary_v;
+        {
+          printf("friction %d %d %g %g %g overshoot\n",
+                 tangent_dim,inds[tangent_dim],
+                 *result, boundary_v, v_off[tangent_dim]);
+          *result=boundary_v;
+        }
 
-      printf("friction %d %d %g %g %g %g %g %g slide\n",
+
+      printf("friction %d %d %g %g %g %g %g %g %g slide\n",
              tangent_dim, inds[tangent_dim],
              *result,delta,
-             yield_stress,normal_stress,p,v_off[tangent_dim]);
+             yield_stress,normal_stress,p,v_off[tangent_dim],tangent_stress);
 
 
     }