← Back to team overview

yade-dev team mailing list archive

[Branch ~yade-dev/yade/trunk] Rev 2972: - attempt to fix https://bugs.launchpad.net/yade/+bug/897237

 

------------------------------------------------------------
revno: 2972
committer: Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
branch nick: trunk
timestamp: Wed 2011-11-30 14:30:02 +0100
message:
  - attempt to fix https://bugs.launchpad.net/yade/+bug/897237
  - small fix in the theoretical background
modified:
  doc/sphinx/formulation.rst
  pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.cpp
  pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.hpp


--
lp:yade
https://code.launchpad.net/~yade-dev/yade/trunk

Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-dev/yade/trunk/+edit-subscription
=== modified file 'doc/sphinx/formulation.rst'
--- doc/sphinx/formulation.rst	2011-06-08 16:21:31 +0000
+++ doc/sphinx/formulation.rst	2011-11-30 13:30:02 +0000
@@ -694,7 +694,7 @@
 General mass-spring system
 """"""""""""""""""""""""""
 				
-In a general mass-spring system, the highest frequency occurs if two connected masses $m_i$, $m_j$ are in opposite motion; let us suppose they have equal velocities (which is conservative) and they are connected by a spring with stiffness $K_{i}$: displacement $\Delta x_i$ of $m_i$ will be accompained by $\Delta x_j=-\Delta x_i$ of $m_j$, giving $\Delta F_i=-K_{i}(\Delta x_i-(-\Delta x_i))=-2K_{i}\Delta x_i$. That results in apparent stiffness $K_{i}^{(2)}=2K_{i}$, giving maximum angular frequency of the whole system
+In a general mass-spring system, the highest frequency occurs if two connected masses $m_i$, $m_j$ are in opposite motion; let us suppose they have equal velocities (which is conservative) and they are connected by a spring with stiffness $K_{i}$: displacement $\Delta x_i$ of $m_i$ will be accompained by $\Delta x_j=-\Delta x_i$ of $m_j$, giving $\Delta F_i=-K_{i}(\Delta x_i-(-\Delta x_i))=-2K_{i}\Delta x_i$. That results in apparent stiffness $K_{i}^{(2)}=2K_{i}$, giving maximum eigenfrequency of the whole system
 
 .. math:: \omega_{\rm max}=\max_i\sqrt{K_i^{(2)}/m_i}.
 			
@@ -705,7 +705,7 @@
 
 	\Dtcr=\frac{2}{\omega_{\rm max}}=\min_i\, 2\sqrt{\frac{m_i}{K_i^{(2)}}}=\min_i\, 2\sqrt{\frac{m_i}{2K_i}}=\min_i \sqrt{2}\sqrt{\frac{m_i}{K_i}}.
 
-This equation can be used for all 6 degrees of freedom (DOF) in translation and rotation, by considering generalized mass and stiffness matrices $M$ and $K$, and replacing fractions $\frac{m_i}{K_i}$ by eigen values of $K.M^{-1}$. The critical timestep is then associated to the eigen mode with highest frequency : 
+This equation can be used for all 6 degrees of freedom (DOF) in translation and rotation, by considering generalized mass and stiffness matrices $M$ and $K$, and replacing fractions $\frac{m_i}{K_i}$ by eigen values of $M.K^{-1}$. The critical timestep is then associated to the eigen mode with highest frequency :
 
 .. math::
 	:label: eq-dtcr-axes

=== modified file 'pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.cpp'
--- pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.cpp	2011-11-15 08:36:34 +0000
+++ pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.cpp	2011-11-30 13:30:02 +0000
@@ -73,14 +73,19 @@
                 if (!bodiesMenisciiList.initialized) bodiesMenisciiList.prepare(scene);
                 //bodiesMenisciiList.display();
         }
-
-        InteractionContainer::iterator ii    = scene->interactions->begin();
+ 
+ 	InteractionContainer::iterator ii    = scene->interactions->begin();
         InteractionContainer::iterator iiEnd = scene->interactions->end();
-
+	bool hertzInitialized = false;
         for(  ; ii!=iiEnd ; ++ii ) {
-		
                 if ((*ii)->isReal()) {
                         const shared_ptr<Interaction>& interaction = *ii;
+			if (!hertzInitialized) {//NOTE: We are assuming that only one type is used in one simulation here
+				if (CapillaryPhys::getClassIndexStatic()==interaction->phys->getClassIndex()) hertzOn=false;
+				else if (MindlinCapillaryPhys::getClassIndexStatic()==interaction->phys->getClassIndex()) hertzOn=true;
+				else LOG_ERROR("The capillary law is not implemented for interactions using"<<interaction->phys->getClassName());}
+			hertzInitialized = true;
+
                         unsigned int id1 = interaction->getId1();
                         unsigned int id2 = interaction->getId2();
                         
@@ -120,7 +125,7 @@
                         /// intergranular distance
                         Real D = alpha*((b2->state->pos-b1->state->pos).norm()-(currentContactGeometry->radius1+ currentContactGeometry->radius2)); // scGeom->penetrationDepth could probably be used here?
 
-                        if ((currentContactGeometry->penetrationDepth>=0)||(D<=0)) { //||(scene->iter < 1) ) // a simplified way to define meniscii everywhere
+                        if ((currentContactGeometry->penetrationDepth>=0)|| D<=0 || createDistantMeniscii) { //||(scene->iter < 1) ) // a simplified way to define meniscii everywhere
                                 D=0; // defines Fcap when spheres interpenetrate //FIXME : D<0 leads to wrong interpolation has D<0 has no solution in the interpolation : this is not physically interpretable!! even if, interpenetration << grain radius.
                                 if (!hertzOn){
 					if (fusionDetection && !cundallContactPhysics->meniscus) bodiesMenisciiList.insert((*ii));
@@ -162,7 +167,7 @@
 					}
 					if (!Vinterpol){
 						if (fusionDetection) bodiesMenisciiList.remove((*ii));
-						scene->interactions->requestErase(id1,id2);
+						if (D>0) scene->interactions->requestErase(id1,id2);
 					}
 		
 					/// wetting angles

=== modified file 'pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.hpp'
--- pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.hpp	2011-11-15 08:36:34 +0000
+++ pkg/dem/Law2_ScGeom_CapillaryPhys_Capillarity.hpp	2011-11-30 13:30:02 +0000
@@ -93,11 +93,12 @@
 		void action();
 		void postLoad(Law2_ScGeom_CapillaryPhys_Capillarity&);
 		
-	YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom_CapillaryPhys_Capillarity,GlobalEngine,"This law allows one to take into account capillary forces/effects between spheres coming from the presence of interparticular liquid bridges (menisci).\n\nThe control parameter is the capillary pressure (or suction) Uc = ugas - Uliquid. Liquid bridges properties (volume V, extent over interacting grains delta1 and delta2) are computed as a result of the defined capillary pressure and of the interacting geometry (spheres radii and interparticular distance).\n\nReferences: in english [Scholtes2009b]_; more detailed, but in french [Scholtes2009d]_.\n\nThe law needs ascii files M(r=i) with i=R1/R2 to work (see https://yade-dem.org/index.php/CapillaryTriaxialTest). These ASCII files contain a set of results from the resolution of the Laplace-Young equation for different configurations of the interacting geometry.\n\nIn order to allow capillary forces between distant spheres, it is necessary to enlarge the bounding boxes using the Ig2's :yref:`interactionDetectionFactor<Ig2_Sphere_Sphere_ScGeom::interactionDetectionFactor>`, and also to disable interactions removal by the constitutive law (:yref:`Law2<Law2_ScGeom_FrictPhys_CundallStrack::neverErase>=True`). The only combinations of laws supported are currently capillary law + :yref:`Law2_ScGeom_FrictPhys_CundallStrack` and capillary law + :yref:`Law2_ScGeom_MindlinPhys_Mindlin` (and the other variants of Hertz-Mindlin).\n\nSee CapillaryPhys-example.py for an example script.",
+	YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom_CapillaryPhys_Capillarity,GlobalEngine,"This law allows one to take into account capillary forces/effects between spheres coming from the presence of interparticular liquid bridges (menisci).\n\nThe control parameter is the capillary pressure (or suction) Uc = ugas - Uliquid. Liquid bridges properties (volume V, extent over interacting grains delta1 and delta2) are computed as a result of the defined capillary pressure and of the interacting geometry (spheres radii and interparticular distance).\n\nReferences: in english [Scholtes2009b]_; more detailed, but in french [Scholtes2009d]_.\n\nThe law needs ascii files M(r=i) with i=R1/R2 to work (see https://yade-dem.org/index.php/CapillaryTriaxialTest). These ASCII files contain a set of results from the resolution of the Laplace-Young equation for different configurations of the interacting geometry.\n\nIn order to allow capillary forces between distant spheres, it is necessary to enlarge the bounding boxes using :yref:`Bo1_Sphere_Aabb::aabbEnlargeFactor` and make the Ig2 define define distant interactions via :yref:`interactionDetectionFactor<Ig2_Sphere_Sphere_ScGeom::interactionDetectionFactor>`. It is also necessary to disable interactions removal by the constitutive law (:yref:`Law2<Law2_ScGeom_FrictPhys_CundallStrack::neverErase>=True`). The only combinations of laws supported are currently capillary law + :yref:`Law2_ScGeom_FrictPhys_CundallStrack` and capillary law + :yref:`Law2_ScGeom_MindlinPhys_Mindlin` (and the other variants of Hertz-Mindlin).\n\nSee CapillaryPhys-example.py for an example script.",
 	((Real,CapillaryPressure,0.,,"Value of the capillary pressure Uc defines as Uc=Ugas-Uliquid"))
 	((bool,fusionDetection,false,,"If true potential menisci overlaps are checked"))
 	((bool,binaryFusion,true,,"If true, capillary forces are set to zero as soon as, at least, 1 overlap (menisci fusion) is detected"))
-	((bool,hertzOn,false,,"Has to be true, if hertz model is set by user (Ip2_FrictMat_FrictMat_MindlinCapillaryPhys)"))
+	((bool,hertzOn,false,,"|yupdate| true if hertz model is used"))
+	((bool,createDistantMeniscii,false,,"Generate meniscii between distant spheres? Else only maintain the existing one. For modeling a wetting path this flag should always be false. For a drying path it should be true for one step (initialization) then false, as in the logic of [Scholtes2009c]_"))
 	 );
 };