← Back to team overview

yade-dev team mailing list archive

  1. yade-dev team
  2. Mailing list archive
  3. Message #04743

[Branch ~yade-dev/yade/trunk] Rev 2273: **** "My" files *****

  Thread Previous • Date Previous • • Thread Next 
Merge authors:
  jduriez <jduriez@c1solimara-l>
------------------------------------------------------------
revno: 2273 [merge]
committer: jduriez <jduriez@c1solimara-l>
branch nick: trunk
timestamp: Thu 2010-06-03 19:49:44 +0200
message:
  **** "My" files *****
  - Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity is now a functor as it should.
  - Related changes in few files
  - The script NormalInelasticityTest tests now the tangential component also, and it will normaly soon check the MomentLaw. It illustrates also the various ways of controlling disp$
  
  
  **** "Common files" ******
  - Change of .Length() (presented as deprecated) in .norm() in some other python scripts
  - I dared a change in the doc of "State", about the blockedDOFS (they impose, according to current version of NewtonIntegrator, zero acceleration, not zero speed, as it was alread$
  - a typo in StepDisplacer's doc
modified:
  core/State.hpp
  examples/concrete/interaction-histogram.py
  pkg/common/Engine/PartialEngine/StepDisplacer.hpp
  pkg/dem/DataClass/InteractionPhysics/NormalInelasticityPhys.hpp
  pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.cpp
  pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.hpp
  pkg/dem/PreProcessor/SimpleShear.cpp
  scripts/NormalInelasticityTest.py
  scripts/test/wm3-wrap.py


--
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 'core/State.hpp'
--- core/State.hpp	2010-04-26 13:58:23 +0000
+++ core/State.hpp	2010-06-03 17:49:44 +0000
@@ -75,7 +75,7 @@
 			((ori,se3.orientation)),
 		/* ctor */,
 		/*py*/
-		.add_property("blockedDOFs",&State::blockedDOFs_vec_get,&State::blockedDOFs_vec_set,"Degress of freedom where linear/angular velocity will be always zero, regardless of applied force/torque. List of any combination of 'x','y','z','rx','ry','rz'.")
+		.add_property("blockedDOFs",&State::blockedDOFs_vec_get,&State::blockedDOFs_vec_set,"Degress of freedom where linear/angular velocity will be always constant (equal to zero, or to an user-defined value), regardless of applied force/torque. List of any combination of 'x','y','z','rx','ry','rz'.")
 		// references must be set using wrapper funcs
 		.add_property("pos",&State::pos_get,&State::pos_set,"Current position.")
 		.add_property("ori",&State::ori_get,&State::ori_set,"Current orientation.") 

=== modified file 'examples/concrete/interaction-histogram.py'
--- examples/concrete/interaction-histogram.py	2010-05-29 20:52:10 +0000
+++ examples/concrete/interaction-histogram.py	2010-06-03 09:28:04 +0000
@@ -21,7 +21,7 @@
 	if not i.isReal: continue
 	norm=i.geom.normal
 	angle=atan(norm[ax2]/norm[ax1])
-	force=i.phys.normalForce.Length()
+	force=i.phys.normalForce.norm()
 	angles.append(angle)
 	forces.append(force)
 # easier: plain histogram

=== modified file 'pkg/common/Engine/PartialEngine/StepDisplacer.hpp'
--- pkg/common/Engine/PartialEngine/StepDisplacer.hpp	2010-05-02 16:02:29 +0000
+++ pkg/common/Engine/PartialEngine/StepDisplacer.hpp	2010-06-03 17:49:44 +0000
@@ -6,7 +6,7 @@
 class StepDisplacer: public PartialEngine {
 	public:
 		virtual void action();
-	YADE_CLASS_BASE_DOC_ATTRS(StepDisplacer,PartialEngine,"Apply generalized displacement (displacement or rotation) stewise on subscribed bodies.",
+	YADE_CLASS_BASE_DOC_ATTRS(StepDisplacer,PartialEngine,"Apply generalized displacement (displacement or rotation) stepwise on subscribed bodies.",
 		((Se3r,deltaSe3,Se3r(Vector3r::Zero(),Quaternionr::Identity()),"Difference of position/orientation that will be added. Position is added to current :yref:`State::pos` using vector addition, orientation to :yref:`State::ori` using quaternion multiplication (rotation composition)."))
 		((bool,setVelocities,false,"If true,  velocity and angularVelocity are modified in such a way that over one iteration (dt), the body will have prescribed se3 jump. In this case, se3 itself is not updated for :yref:`dynamic<Body::isDynamic>` bodies, since they would have the delta applied twice (here and in the :yref:`integrator<NewtonIntegrator>`). For non-dynamic bodies however, se3 *is* still updated."))
 	);

=== modified file 'pkg/dem/DataClass/InteractionPhysics/NormalInelasticityPhys.hpp'
--- pkg/dem/DataClass/InteractionPhysics/NormalInelasticityPhys.hpp	2010-04-25 13:18:11 +0000
+++ pkg/dem/DataClass/InteractionPhysics/NormalInelasticityPhys.hpp	2010-06-03 10:12:50 +0000
@@ -8,7 +8,8 @@
 
 #pragma once
 
-#include<yade/pkg-dem/FrictPhys.hpp>
+#include<yade/pkg-dem/FrictPhys.hpp> // ou 
+// #include <yade/pkg/dem/DataClass/InteractionPhysics/FrictPhys.hpp>
 
 /*! \brief Interaction for using Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity
 

=== modified file 'pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.cpp'
--- pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.cpp	2010-05-18 22:10:18 +0000
+++ pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.cpp	2010-06-03 17:49:44 +0000
@@ -9,8 +9,6 @@
 #include<yade/pkg-dem/NormalInelasticityLaw.hpp>
 
 #include<yade/pkg-dem/CohFrictMat.hpp>
-#include<yade/pkg-dem/ScGeom.hpp>
-#include<yade/pkg-dem/NormalInelasticityPhys.hpp>
 #include<yade/core/Omega.hpp>
 #include<yade/core/Scene.hpp>
 
@@ -18,33 +16,19 @@
 
 
 
-void Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity::action()// a remplacer par :
-// void Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity::go(shared_ptr<InteractionGeometry>& iG, shared_ptr<InteractionPhysics>& iP, Interaction* contact, Scene* scene)
+void Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity::go(shared_ptr<InteractionGeometry>& iG, shared_ptr<InteractionPhysics>& iP, Interaction* contact, Scene* scene)
 {
-// 	cout << "\n Nvlle it :"<< endl;
-	//shared_ptr<BodyContainer>& bodies = scene->bodies;				//It gave a warning. Anton Gladky.
 
 	const Real& dt = scene->dt;
 
-	InteractionContainer::iterator ii    = scene->interactions->begin();	// a supprimer pr passage au go
-	InteractionContainer::iterator iiEnd = scene->interactions->end();	// a supprimer pr passage au go
-	int nbreInteracTot=0;// a supprimer pr passage au go
-	int nbreInteracMomPlastif=0;// a supprimer pr passage au go
-	for (  ; ii!=iiEnd ; ++ii )// a supprimer pr passage au go
-	{// a supprimer pr passage au go
-        if ((*ii)->interactionGeometry && (*ii)->interactionPhysics)
-	if ((*ii)->isReal())// a supprimer pr passage au go
-		{
-		nbreInteracTot++;
-		const shared_ptr<Interaction>& contact = *ii;// supprimable
 		int id1 = contact->getId1();
 		int id2 = contact->getId2();
 // 		cout << "contact entre " << id1 << " et " << id2 << " reel ? " << contact->isReal() << endl;
 
 		State* de1 = Body::byId(id1,scene)->state.get();
 		State* de2 = Body::byId(id2,scene)->state.get();
-		ScGeom* currentContactGeometry		= YADE_CAST<ScGeom*>(contact->interactionGeometry.get());
-		NormalInelasticityPhys* currentContactPhysics = YADE_CAST<NormalInelasticityPhys*> (contact->interactionPhysics.get());
+		ScGeom* currentContactGeometry		= YADE_CAST<ScGeom*>(iG.get());
+		NormalInelasticityPhys* currentContactPhysics = YADE_CAST<NormalInelasticityPhys*> (iP.get());//remplace par :
 
 		Vector3r& shearForce 			= currentContactPhysics->shearForce;
 
@@ -99,10 +83,6 @@
 					 // probably not useful anymore
                 currentContactPhysics->normalForce = Vector3r::Zero();
                 currentContactPhysics->shearForce = Vector3r::Zero();
-
-                //return;
-                //    else
-                //    currentContactPhysics->normalForce	= -currentContactPhysics->normalAdhesion*currentContactGeometry->normal;
             	}
             else
             	{
@@ -217,7 +197,7 @@
 				if(normeMoment>normeMomentMax)
 					{
 					moment *= normeMomentMax/normeMoment;
-					nbreInteracMomPlastif++;
+// 					nbreInteracMomPlastif++;
 					}
 			}
 			scene->forces.addTorque(id1,-moment);
@@ -229,10 +209,7 @@
 //             }
 //         }
     }
-// 	cout << " En tout :"<< nbreInteracTot << "interactions (reelles)" << endl;
-    //cerr << "ncount= " << ncount << endl;//REMOVE
-// 	cout << momentAlwaysElastic << endl;
-// 	cout << "Sur " << nbreInteracTot << " interactions (reelles) " << nbreInteracMomPlastif << " se sont vues corriger leur moment" << endl;
-}
-}
-}
+}
+
+
+

=== modified file 'pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.hpp'
--- pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.hpp	2010-05-18 10:28:02 +0000
+++ pkg/dem/Engine/GlobalEngine/NormalInelasticityLaw.hpp	2010-06-03 17:49:44 +0000
@@ -8,28 +8,30 @@
 
 #pragma once
 
-#include<yade/core/GlobalEngine.hpp> // a remplacer par :
+// #include<yade/core/GlobalEngine.hpp> // a remplacer par :
 #include<yade/pkg-common/LawFunctor.hpp>
+#include<yade/pkg-dem/ScGeom.hpp>
+#include<yade/pkg-dem/NormalInelasticityPhys.hpp>
 #include <set>
 #include <boost/tuple/tuple.hpp>
 
 
 
-class Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity : public GlobalEngine // a remplacer par :
-// class Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity : public LawFunctor
+class Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity : public LawFunctor
 {
 	public :
-		void action();// a remplacer par :
-// 		virtual void go(shared_ptr<InteractionGeometry>&, shared_ptr<InteractionPhysics>&, Interaction*, Scene*);
+		virtual void go(shared_ptr<InteractionGeometry>&, shared_ptr<InteractionPhysics>&, Interaction*, Scene*);
+
+	FUNCTOR2D(ScGeom,NormalInelasticityPhys);
 
 	YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity,
-// 				LawFunctor,
-				GlobalEngine,
-				"Contact law including cohesion, moment transfer and inelastic compression behaviour\n\n This contact Law is inspired by :yref:`CohesiveFrictionalContactLaw` (inspired itselve directly from the work of Plassiard & Belheine, see the corresponding articles in (Annual Report 2006) in http://geo.hmg.inpg.fr/frederic/Discrete_Element_Group_FVD.html for example).\n\n It allows so to set moments, cohesion, tension limit and (that's the difference) inelastic unloadings in compression between bodies. All that concerned brokenBodies (this flag and the erosionactivated one) and the useless 'iter' has been suppressed.\n\n The Relationsships corresponding are Ip2_2xCohFrictMat_NormalInelasticityPhys, where the rigidities, the friction angles (with their tan()), and the orientations of the interactions are calculated. No more cohesion and tension limits are computed for all the interactions.\n\n To use it you should also use :\n- :yref:`CohFrictMat` for bodies, with *isCohesive* = 1 (A verifier ce dernier point)\n- :yref:`Ip2_2xCohFrictMat_NormalInelasticityPhys` (=> which involves interactions of :yref:`NormalInelasticityPhys` type).\n\n The effect of this law on the normal force are illustrated in scripts/NormalInelasticityTest.py",
+				LawFunctor,
+				"Contact law including cohesion, moment transfer and inelastic compression behaviour\n\n This contact Law is inspired by :yref:`CohesiveFrictionalContactLaw` (inspired itselve directly from the work of Plassiard & Belheine, see the corresponding articles in [DeghmReport2006]_ for example).\n\n It allows so to set moments, cohesion, tension limit and (that's the difference) inelastic unloadings in compression between bodies. All that concerned brokenBodies (this flag and the erosionactivated one) and the useless 'iter' has been suppressed.\n\n The Relationsships corresponding are Ip2_2xCohFrictMat_NormalInelasticityPhys, where the rigidities, the friction angles (with their tan()), and the orientations of the interactions are calculated. No more cohesion and tension limits are computed for all the interactions.\n\n To use it you should also use :\n- :yref:`CohFrictMat` for bodies, with *isCohesive* = 1 (A verifier ce dernier point)\n- :yref:`Ip2_2xCohFrictMat_NormalInelasticityPhys` (=> which involves interactions of :yref:`NormalInelasticityPhys` type).\n\n The effects of this law are illustrated in scripts/NormalInelasticityTest.py",
 				((Real,coeff_dech,1.0,"=kn(unload) / kn(load)"))
 				((bool,momentRotationLaw,true,"boolean, true=> computation of a torque (against relative rotation) exchanged between particles"))
 				((bool,momentAlwaysElastic,false,"boolean, true=> the torque (computed only if momentRotationLaw !!) is not limited by a plastic threshold"))
 				);
+	
 };
 
 REGISTER_SERIALIZABLE(Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity);

=== modified file 'pkg/dem/PreProcessor/SimpleShear.cpp'
--- pkg/dem/PreProcessor/SimpleShear.cpp	2010-05-18 10:28:02 +0000
+++ pkg/dem/PreProcessor/SimpleShear.cpp	2010-06-02 16:40:46 +0000
@@ -285,7 +285,7 @@
 	rootBody->engines.push_back(shared_ptr<Engine>(new InsertionSortCollider));
 	rootBody->engines.push_back(interactionGeometryDispatcher);
 	rootBody->engines.push_back(interactionPhysicsDispatcher);
-	rootBody->engines.push_back(shared_ptr<Engine>(new Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity));
+// 	rootBody->engines.push_back(shared_ptr<Engine>(new Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity));
 	if(gravApplied)
 		rootBody->engines.push_back(gravityCondition);
 	if(shearApplied)

=== modified file 'scripts/NormalInelasticityTest.py'
--- scripts/NormalInelasticityTest.py	2010-05-29 20:52:10 +0000
+++ scripts/NormalInelasticityTest.py	2010-06-03 17:49:44 +0000
@@ -1,10 +1,13 @@
+# -*- coding: utf-8 -*-
 
-# Script to test the constitutive law contained in NormalInelasticityLaw : consider two spheres whose penetration of the contact evolves => Monitor of the normal force
+# Script to test the constitutive law contained in NormalInelasticityLaw : consider two spheres moving one to each other (this script illustrates different ways of moving spheres)
+#- first penetration of the contact evolves => Monitor of the normal force
+#- then, test in tangential direction
 
 from yade import plot
 
 #Def of the material which will be used
-O.materials.append(CohesiveFrictionalMat(density=2600,young=4.0e9,poisson=.04,frictionAngle=.6,label='Materiau1'))
+O.materials.append(CohFrictMat(density=2600,young=4.0e9,poisson=.04,frictionAngle=.6,label='Materiau1'))
 
 #Def of the bodies of the simulations : 2 spheres, with names which will be useful after
 O.bodies.append(utils.sphere([0,0,0], 1, dynamic=False, wire=False, color=None, highlight=False)) #'Materiau1', as the latest material defined will be used
@@ -19,15 +22,17 @@
 	ForceResetter(),
 	BoundDispatcher([Bo1_Sphere_Aabb()]),
 	InsertionSortCollider(),
-	InteractionGeometryDispatcher([Ig2_Sphere_Sphere_ScGeom()]),
-	InteractionPhysicsDispatcher([Ip2_2xCohFrictMat_NormalInelasticityPhys()]),
-	PeriodicPythonRunner(iterPeriod=1,command='letMove()'),
-	Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity(coeff_dech=3)
-]
+	InteractionDispatchers(
+			      [Ig2_Sphere_Sphere_ScGeom()],
+			      [Ip2_2xCohFrictMat_NormalInelasticityPhys()],
+			      [Law2_ScGeom_NormalInelasticityPhys_NormalInelasticity(coeff_dech=3)]
+			      ),
+	PeriodicPythonRunner(iterPeriod=1,command='letMove()')
+	]
 
 
 #Def of the python command letMove() :
-# Will move "by hand" the upperSphere towards or away from the lower one. Modifying by hand only the speed of bodies is indeed not sufficient, see NewtonsIntegrator, and https://bugs.launchpad.net/yade/+bug/398089. For such purposes you could also use TranslationEngine, or you can also use StepDisplacer, which applies finite change in position/orientation in each step
+# Will move "by hand" the upperSphere towards or away from the lower one. Modifying by hand only the speed of bodies is indeed not sufficient, see NewtonsIntegrator, and https://bugs.launchpad.net/yade/+bug/398089. Alternative way is presented below
 def letMove():#Load for the first 10 iterations, unload for the 7 following iterations, then reload
 	vImposed=[0,-1,0]
 	if O.iter < 25 and O.iter>14:
@@ -41,10 +46,12 @@
 	i=O.interactions[1,0]
 	VecFn=i.phys.normalForce
 	VecDist=UpperSphere.state.pos-LowerSphere.state.pos
-	plot.addData(Normfn=VecFn.Length(),FnY=VecFn[1],step=O.iter,un=LowerSphere.shape.radius+UpperSphere.shape.radius-VecDist.Length()) # the 1e-5 because of the "shift2" which appears in the computation of the penetration depth made in Ig2_Sphere_Sphere_ScGeom. It seems that this shift2 = 1e-5
-
-
-
+	plot.addData(Normfn=VecFn.norm(),NormfnBis=VecFn.norm(),FnY=VecFn[1],step=O.iter,
+	  unPerso=LowerSphere.shape.radius+UpperSphere.shape.radius-VecDist.norm(),unVrai=i.geom.penetrationDepth,
+	  gamma=UpperSphere.state.pos[0]-LowerSphere.state.pos[0],Fx=O.forces.f(0)[0])
+
+
+# ------ Test of the law in the normal direction, using python commands to let move ------ #
 
 O.dt=1e-5
 
@@ -52,11 +59,44 @@
 O.run(2,True) #cycles "for free", so that the interaction between spheres will be defined (with his physics and so on)
 O.engines=O.engines+[PeriodicPythonRunner(iterPeriod=1,command='defData()')]
 
-#O.saveTmp('INL')
+
+
 O.run(40,True)
 
 # define of the plots to be made : un(step), and Fn(un)
-plot.plots={'step':('un',),'un':('Normfn',)}
-plot.plot()
-
-#NB : the shape of the curve Fn(un) seems to not be perfect. It is indeed not because of NormalInelasticityLaw. But of differences between the un computed here in this python script and the one which is computed in Ig2_Sphere_Sphere_ScGeom (see for example this "shift2"). Fn being linked to this last un, these slight differences explain the shape of the curves. If phys.penetrationDepth would exist in python, and thus could be directly considered, I think the curve would be perfect !
+plot.plots={'step':('unVrai',),'unPerso':('Normfn',),'unVrai':('NormfnBis',)}
+plot.plot()
+
+#NB : these different unVrai and unPerso illustrate the definition of penetrationDepth really used in the code (computed in Ig2_Sphere_Sphere_ScGeom) which is slightly different from R1 + R2 - Distance (see for example this "shift2"). According to the really used penetrationDepth, Fn evolves as it should
+
+O.saveTmp('EndComp')
+
+
+# ------ Test of the law in the tangential direction, using StepDisplacer ------ #
+
+DPos=Vector3.Zero
+Vector3.__init__(DPos,1*O.dt,0,0)
+
+O.engines=O.engines[:4]+[StepDisplacer(subscribedBodies=[1],deltaSe3=(DPos,Quaternion.Identity),setVelocities=True)]+O.engines[5:]
+O.run(1000)
+plot.plots={'step':('gamma',),'gamma':('Fx',)}
+plot.plot()
+plot.plots={'Normfn':('Fx',)}
+plot.plot()
+#Comments => 	- evolution of Fx with gamma normal (flat at the beginning because of the order of engines)
+#		- un decreases indeed during this shear, but a zoom on the curves is needed to see it.
+#		- We can observe that the force state of the sample decreases a line with a slope equal to tan(~34.5°)=tan(~0.602 rad). Why not strict equality ? Because of the measure of the slope or because something else ? To see...
+
+
+
+# ------ Test of the law for the moment, using blockedDOF_s ------ #
+O.loadTmp('EndComp')
+
+#To use blockedDOF_s, the body has to be dynamic....
+UpperSphere.isDynamic=True
+UpperSphere.state.blockedDOFs='x','rx','y','ry','z','rz'
+UpperSphere.state.angVel=Vector3(0,0,1)
+
+
+
+

=== modified file 'scripts/test/wm3-wrap.py'
--- scripts/test/wm3-wrap.py	2009-09-19 19:41:52 +0000
+++ scripts/test/wm3-wrap.py	2010-06-03 10:12:50 +0000
@@ -9,7 +9,7 @@
 x.Dot(y)==0
 x.Cross(y)==z
 # methods
-one.Length()
+one.norm()
 
 # quaternions
 # construction (implicit conversion of 3-tuple or list of length 3 to Vector3)
  Thread Previous • Date Previous • • Thread Next