[cig-commits] r11971 - in long/3D/Gale/trunk: . src/Gale/Utils/src
walter at geodynamics.org
walter at geodynamics.org
Wed May 14 23:00:09 PDT 2008
Author: walter
Date: 2008-05-14 23:00:09 -0700 (Wed, 14 May 2008)
New Revision: 11971
Modified:
long/3D/Gale/trunk/
long/3D/Gale/trunk/src/Gale/Utils/src/KineticFriction.c
long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.c
long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.h
Log:
r2145 at earth: boo | 2008-05-14 21:13:40 -0700
Icky, hard-coded, non-parallel hack to make kinetic friction work.
It does two things:
1) It does not apply the full friction force, but instead averages it
with the friction force from the previous iteration.
2) If the direction of the force is opposite from the previous
iteration, it reduces the force by a factor of 100.
These two things prevent oscillatory behavior which prevents
convergence during the non-linear iteration. Either one of these
individually does not damp down the oscillations enough. I also tried
setting the force to zero (rather than reduce by 100), but that did not
work as well.
Property changes on: long/3D/Gale/trunk
___________________________________________________________________
Name: svk:merge
- 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2144
+ 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2145
Modified: long/3D/Gale/trunk/src/Gale/Utils/src/KineticFriction.c
===================================================================
--- long/3D/Gale/trunk/src/Gale/Utils/src/KineticFriction.c 2008-05-15 06:00:06 UTC (rev 11970)
+++ long/3D/Gale/trunk/src/Gale/Utils/src/KineticFriction.c 2008-05-15 06:00:09 UTC (rev 11971)
@@ -384,6 +384,8 @@
double p;
+ double *force_vector=StaticFrictionVC_get_force_vector();
+
Stream* errorStr = Journal_Register( Error_Type, self->type );
switch(self->friction_entry.type)
{
@@ -473,6 +475,29 @@
Vec_Sub3D(v_surface_diff,v,v_perpendicular);
v_norm=Vec_Mag3D(v_surface_diff);
}
+
+ StaticFrictionVC_get_force_vector(force_vector);
+
+/* if(force_vector[elementNodes[eNode_I]]>0 */
+/* && v_norm>1.0e-11) */
+/* { */
+/* for(d=0;d<dim;++d) */
+/* { */
+/* elForceVec[ eNode_I*nodeDofCount + d]+=force_vector[elementNodes[eNode_I]]*(v_surface_diff[d]/v_norm)*area/overcount; */
+/* } */
+/* printf("kinetic slip %d %d %d %g %g %g %g\n", */
+/* lElement_I,eNode_I,elementNodes[eNode_I], */
+/* force_vector[elementNodes[eNode_I]]*area/overcount, */
+/* v_norm,v_surface_diff[0],v_surface_diff[1]); */
+/* } */
+/* else */
+/* { */
+/* printf("kinetic stick %d %d %d %g %g %g %g\n", */
+/* lElement_I,eNode_I,elementNodes[eNode_I], */
+/* force_vector[elementNodes[eNode_I]]*area/overcount, */
+/* v_norm,v_surface_diff[0],v_surface_diff[1]); */
+/* } */
+
friction_force=friction_coefficient*normal_stress*area/overcount;
/* if(!(friction_force<0 || v_norm==0)) */
@@ -482,37 +507,67 @@
{
for(d=0;d<dim;++d)
{
- elForceVec[ eNode_I*nodeDofCount + d]+=friction_force*v_surface_diff[d]/v_norm;
+ double alpha=0.5;
+ double old_force=
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d];
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]=
+ friction_force*v_surface_diff[d]/v_norm;
+ if(old_force==0)
+ alpha=1.0;
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]=
+ (friction_force*v_surface_diff[d]/v_norm)*alpha
+ + old_force*(1-alpha);
+
+ if(force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]*old_force<0)
+ {
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]/=100.0;
+ printf("kinetic reduced\n");
+ }
+/* force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]= */
+/* (friction_force*v_surface_diff[d]/v_norm-old_force) */
+/* - tangential_stress[d]*area/overcount; */
+
+ elForceVec[ eNode_I*nodeDofCount + d]+=
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d];
+
+/* elForceVec[ eNode_I*nodeDofCount + d]+=friction_force*v_surface_diff[d]/v_norm; */
/* elForceVec[ eNode_I*nodeDofCount + d]+=friction_force*v_surface_diff[d]/v_norm - tangential_stress[d]*area/overcount; */
/* node_force[ eNode_I*nodeDofCount + d]+=friction_force*v_surface_diff[d]/v_norm; */
/* node_force[ eNode_I*nodeDofCount + d]+=friction_force*v_surface_diff[d]/v_norm - tangential_stress[d]*area/overcount; */
/* for(d=0;d<dim;++d) */
- }
- printf("kinetic slip %d %d %g %g %g %g %g %g %d %g %g\n",
- lElement_I,eNode_I,friction_force,
+
+ printf("kinetic slip %d %d %d %d %g %g %g %g %g %g %g %g %d %g %g\n",
+ (lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d,
+ lElement_I,eNode_I,d,friction_force,old_force,
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d],
v_norm,v_surface_diff[0],v_surface_diff[1],
normal_stress,p,overcount,tangential_stress[0],
tangential_stress[1]);
+ }
}
else
{
- printf("kinetic stick %d %d %g %g %g %g %g %g %d %g %g\n",
- lElement_I,eNode_I,friction_force,
+ for(d=0;d<dim;++d)
+ {
+ double old_force=
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d];
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]=0;
+
+/* force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]= */
+/* - (tangential_stress[d]*area/overcount-old_force); */
+/* elForceVec[ eNode_I*nodeDofCount + d]+= */
+/* force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d]; */
+
+ printf("kinetic stick %d %d %d %d %g %g %g %g %g %g %g %g %d %g %g\n",
+ (lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d,
+ lElement_I,eNode_I,d,friction_force,old_force,
+ force_vector[(lElement_I*elementNodeCount+eNode_I)*nodeDofCount+d],
v_norm,v_surface_diff[0],v_surface_diff[1],
normal_stress,p,overcount,tangential_stress[0],
tangential_stress[1]);
- apply_force=False;
+ }
}
}
}
-
-/* if(apply_force) */
-/* { */
-/* for( eNode_I = 0 ; eNode_I < elementNodeCount; eNode_I++ ) { */
-/* for(d=0;d<dim;++d) */
-/* elForceVec[ eNode_I*nodeDofCount + d]+= */
-/* node_force[eNode_I*nodeDofCount + d]; */
-/* } */
-/* } */
}
Modified: long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.c
===================================================================
--- long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.c 2008-05-15 06:00:06 UTC (rev 11970)
+++ long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.c 2008-05-15 06:00:09 UTC (rev 11971)
@@ -39,6 +39,7 @@
extern const char* WallEnumToStr[Wall_Size];
+static double force_vector[129*65*2*4*2];
/*--------------------------------------------------------------------------------------------------------------------------
** Constructor
@@ -617,6 +618,7 @@
unsigned nDims;
Grid* grid;
FeVariable *deviatoric_stress, *pressure;
+ unsigned nnn;
nDims = Mesh_GetDimSize( self->_mesh );
grid = *(Grid**)ExtensionManager_Get( self->_mesh->info, self->_mesh,
@@ -635,6 +637,10 @@
self->deviatoric_stress_name);
+
+ for(nnn=0;nnn<129*65*2*4*2;++nnn)
+ force_vector[nnn]=0;
+
/* If this is the first time this routine is called, the
stress may not have been initialized yet. In that case,
make everything stick. */
@@ -663,6 +669,7 @@
unsigned nInc;
unsigned *inc;
unsigned nn;
+
Mesh_GetIncidence(self->_mesh, MT_VERTEX, n_i, nDims, &nInc, &inc);
p=0;
for(nn=0; nn<nInc; ++nn)
@@ -680,8 +687,19 @@
{
if(self->stick_on_first_step)
{
- IndexSet_Add(set,n_i);
- printf("static stick initial %d\n",n_i);
+ if(StaticFrictionVC_get_interface_stresses
+ (deviatoric_stress,self->_mesh,
+ grid->sizes, n_i, ijk, nDims, basis,
+ self->include_lower_boundary,
+ self->include_upper_boundary,
+ &normal_stress,tangential_stress,
+ &tangential_norm,n))
+ {
+ IndexSet_Add(set,n_i);
+ printf("static stick initial %d\n",n_i);
+ }
+ else
+ printf("static slip initial %d\n",n_i);
}
else
printf("static slip initial %d\n",n_i);
@@ -701,8 +719,10 @@
n_i,p,normal_stress,tangential_norm);
}
else
+ {
printf("static slip %d %lf %lf %lf\n",
n_i,p,normal_stress,tangential_norm);
+ }
}
}
}
@@ -1128,3 +1148,8 @@
*normal_stress=-(*normal_stress);
return True;
}
+
+double * StaticFrictionVC_get_force_vector()
+{
+ return force_vector;
+}
Modified: long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.h
===================================================================
--- long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.h 2008-05-15 06:00:06 UTC (rev 11970)
+++ long/3D/Gale/trunk/src/Gale/Utils/src/StaticFrictionVC.h 2008-05-15 06:00:09 UTC (rev 11971)
@@ -214,5 +214,7 @@
void StaticFrictionVC_get_basis_vectors(Wall wall, unsigned int sizes[],
unsigned *direction, unsigned *boundary,
unsigned basis[]);
-
+
+double * StaticFrictionVC_get_force_vector();
+
#endif /* __Discretisation_Utils_StaticFrictionVC_h__ */
More information about the cig-commits
mailing list