← Back to team overview

yade-dev team mailing list archive

[Branch ~yade-dev/yade/trunk] Rev 2536: Fixes https://bugs.launchpad.net/yade/+bug/672473 (an error in the script finally...)

 

------------------------------------------------------------
revno: 2536
committer: jduriez <jduriez@c1solimara-l>
branch nick: yade
timestamp: Mon 2010-11-08 16:26:46 +0100
message:
  Fixes https://bugs.launchpad.net/yade/+bug/672473 (an error in the script finally...)
  
  By the way, more cleaning in NormalInelasticity... files about relative rotations.
modified:
  pkg/dem/Ip2_2xNormalInelasticMat_NormalInelasticityPhys.cpp
  pkg/dem/NormalInelasticityLaw.cpp
  pkg/dem/NormalInelasticityLaw.hpp
  pkg/dem/NormalInelasticityPhys.hpp
  scripts/normalInelasticity-test.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 'pkg/dem/Ip2_2xNormalInelasticMat_NormalInelasticityPhys.cpp'
--- pkg/dem/Ip2_2xNormalInelasticMat_NormalInelasticityPhys.cpp	2010-10-13 16:23:08 +0000
+++ pkg/dem/Ip2_2xNormalInelasticMat_NormalInelasticityPhys.cpp	2010-11-08 15:26:46 +0000
@@ -67,18 +67,9 @@
 			contactPhysics->prevNormal 			= geom->normal;
 			
 			contactPhysics->knLower = Kn;
-			contactPhysics->kn = Kn;
-			
+			contactPhysics->kn = Kn;			
 			contactPhysics->ks = Ks;
-			contactPhysics->initialOrientation1	= Body::byId(interaction->getId1())->state->ori;
-			contactPhysics->initialOrientation2	= Body::byId(interaction->getId2())->state->ori;
-			contactPhysics->initialPosition1    = Body::byId(interaction->getId1())->state->pos;
-			contactPhysics->initialPosition2    = Body::byId(interaction->getId2())->state->pos;
 			contactPhysics->kr = Kr;
-			contactPhysics->initialContactOrientation.setFromTwoVectors(Vector3r(1.0,0.0,0.0),geom->normal);
-			contactPhysics->currentContactOrientation = contactPhysics->initialContactOrientation;
-			contactPhysics->orientationToContact1   = contactPhysics->initialOrientation1.conjugate() * contactPhysics->initialContactOrientation;
-			contactPhysics->orientationToContact2	= contactPhysics->initialOrientation2.conjugate() * contactPhysics->initialContactOrientation;
 		}
 		
 	}

=== modified file 'pkg/dem/NormalInelasticityLaw.cpp'
--- pkg/dem/NormalInelasticityLaw.cpp	2010-11-08 10:10:00 +0000
+++ pkg/dem/NormalInelasticityLaw.cpp	2010-11-08 15:26:46 +0000
@@ -165,10 +165,10 @@
 
 		if(momentRotationLaw)
 		{
-			Vector3r moment_twist( (currentContactGeometry->getTwist()*currentContactPhysics->kr)*currentContactGeometry->normal );
-			Vector3r moment_bending( currentContactGeometry->getBending() * currentContactPhysics->kr );
+			currentContactPhysics->moment_twist = (currentContactGeometry->getTwist()*currentContactPhysics->kr)*currentContactGeometry->normal ;
+			currentContactPhysics->moment_bending = currentContactGeometry->getBending() * currentContactPhysics->kr;
 
-			Vector3r moment = moment_twist + moment_bending;
+			moment = currentContactPhysics->moment_twist + currentContactPhysics->moment_bending;
 
 // 	Limitation by plastic threshold
 			if (!momentAlwaysElastic)

=== modified file 'pkg/dem/NormalInelasticityLaw.hpp'
--- pkg/dem/NormalInelasticityLaw.hpp	2010-11-03 15:50:02 +0000
+++ pkg/dem/NormalInelasticityLaw.hpp	2010-11-08 15:26:46 +0000
@@ -18,16 +18,20 @@
 
 class Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity : public LawFunctor
 {
+	private :
+		Vector3r moment; // contact torque of the interaction
 	public :
 		virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
 
 	FUNCTOR2D(ScGeom,NormalInelasticityPhys);
 
-	YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity,
+	YADE_CLASS_BASE_DOC_ATTRS_CTOR(Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity,
 				LawFunctor,
-				"Contact law used to simulate granulate filler in rock joints [Duriez2009a]_, [Duriez2010]_. It includes possibility of cohesion, moment transfer and inelastic compression behaviour (to reproduce the normal inelasticity observed for rock joints, for the latter).\n\n The moment transfer relation corresponds to the adaptation of the work of Plassiard & Belheine (see in [DeghmReport2006]_ for example), which was realized by J. Kozicki, and is now coded in ScGeom6D.\n\n As others :yref:`LawFunctor`, it uses pre-computed data of the interactions (rigidities, friction angles -with their tan()-, orientations of the interactions); this work is done here in :yref:`Ip2_2xNormalInelasticMat_NormalInelasticityPhys`.\n\n To use this you should also use :yref:`NormalInelasticMat` as material type of the bodies.\n\n The effects of this law are illustrated in scripts/normalInelasticityTest.py",
+				"Contact law used to simulate granulate filler in rock joints [Duriez2009a]_, [Duriez2010]_. It includes possibility of cohesion, moment transfer and inelastic compression behaviour (to reproduce the normal inelasticity observed for rock joints, for the latter).\n\n The moment transfer relation corresponds to the adaptation of the work of Plassiard & Belheine (see in [DeghmReport2006]_ for example), which was realized by J. Kozicki, and is now coded in :yref:`ScGeom6D`.\n\n As others :yref:`LawFunctor`, it uses pre-computed data of the interactions (rigidities, friction angles -with their tan()-, orientations of the interactions); this work is done here in :yref:`Ip2_2xNormalInelasticMat_NormalInelasticityPhys`.\n\n To use this you should also use :yref:`NormalInelasticMat` as material type of the bodies.\n\n The effects of this law are illustrated in scripts/normalInelasticityTest.py",
 				((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"))
+				,
+				moment=Vector3r::Zero();
 				);
 	
 };

=== modified file 'pkg/dem/NormalInelasticityPhys.hpp'
--- pkg/dem/NormalInelasticityPhys.hpp	2010-11-03 15:50:02 +0000
+++ pkg/dem/NormalInelasticityPhys.hpp	2010-11-08 15:26:46 +0000
@@ -17,23 +17,19 @@
 		virtual ~NormalInelasticityPhys();
 
 	YADE_CLASS_BASE_DOC_ATTRS_CTOR(NormalInelasticityPhys,FrictPhys,
-				 "Physics (of interaction) for using :yref:`Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity` : with inelastic unloadings",
-				 ((Real,unMax,0.0,,"the maximum value of penetration depth of the history of this interaction"))
-				 ((Real,previousun,0.0,,"the value of this un at the last time step"))
-				 ((Real,previousFn,0.0,,"the value of the normal force at the last time step"))
-				 ((Quaternionr,initialOrientation1,Quaternionr::Identity(),,""))
-				 ((Quaternionr,initialOrientation2,Quaternionr::Identity(),,""))
-				 ((Quaternionr,orientationToContact1,Quaternionr::Identity(),,""))
-				 ((Quaternionr,orientationToContact2,Quaternionr::Identity(),,""))
-				 ((Quaternionr,currentContactOrientation,Quaternionr::Identity(),,""))
-				 ((Quaternionr,initialContactOrientation,Quaternionr::Identity(),,""))
-				 ((Vector3r,initialPosition1,Vector3r::Zero(),,""))
-				 ((Vector3r,initialPosition2,Vector3r::Zero(),,""))
-				 ((Real,forMaxMoment,1.0,,"parameter stored for each interaction, and allowing to compute the maximum value of the exchanged torque : TorqueMax= forMaxMoment * NormalForce"))
-				 ((Real,kr,0.0,,"the rolling stiffness of the interaction"))
-				 ((Real,knLower,0.0,,"the stifness corresponding to a virgin load for example")),
-				 createIndex();
-				 );
+				"Physics (of interaction) for using :yref:`Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity` : with inelastic unloadings",
+				((Real,unMax,0.0,,"the maximum value of penetration depth of the history of this interaction"))
+				((Real,previousun,0.0,,"the value of this un at the last time step"))
+				((Real,previousFn,0.0,,"the value of the normal force at the last time step"))
+				((Real,forMaxMoment,1.0,,"parameter stored for each interaction, and allowing to compute the maximum value of the exchanged torque : TorqueMax= forMaxMoment * NormalForce"))
+				((Real,kr,0.0,,"the rolling stiffness of the interaction"))
+				((Real,knLower,0.0,,"the stifness corresponding to a virgin load for example"))
+				// internal attributes
+				((Vector3r,moment_twist,Vector3r(0,0,0),(Attr::noSave | Attr::readonly),"Twist moment. Defined here, being initialized as it should be, to be used in :yref:`Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity`"))
+				((Vector3r,moment_bending,Vector3r(0,0,0),(Attr::noSave | Attr::readonly),"Bending moment. Defined here, being initialized as it should be, to be used in :yref:`Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity`"))
+				,
+				createIndex();
+				);
 	REGISTER_CLASS_INDEX(NormalInelasticityPhys,FrictPhys);
 };
 

=== modified file 'scripts/normalInelasticity-test.py'
--- scripts/normalInelasticity-test.py	2010-11-08 10:10:00 +0000
+++ scripts/normalInelasticity-test.py	2010-11-08 15:26:46 +0000
@@ -28,7 +28,7 @@
 	ForceResetter(),
 	InsertionSortCollider([Bo1_Sphere_Aabb()]),
 	InteractionLoop(
-			      [Ig2_Sphere_Sphere_ScGeom()],
+			      [Ig2_Sphere_Sphere_ScGeom6D()],
 			      [Ip2_2xNormalInelasticMat_NormalInelasticityPhys(betaR=0.24)],
 			      [Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity()]
 			      ),
@@ -70,9 +70,10 @@
 
 
 # define of the plots to be made : un(step), and Fn(un)
-plot.plots={'step':('unTrue',),'unPerso':('normFn',),'unTrue':('normFnBis',)}
+plot.plots={'step':('unTrue','torque',),'unPerso':('normFn',),'unTrue':('normFnBis',)}
 plot.plot()
 raw_input()
+print 'Type Return to go ahead'
 print ''
 #NB : these different unTrue 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
 
@@ -89,12 +90,14 @@
 O.engines=O.engines[:3]+[StepDisplacer(ids=[1],mov=dpos,setVelocities=True)]+O.engines[4:]
 O.run(1000)
 print 'End of tangential loading'
-plot.plots={'step':('gamma',),'gamma':('fx',)}
+plot.plots={'step':('gamma','torque',),'gamma':('fx',)}
 plot.plot()
 raw_input()
+print 'Type Return to go ahead'
 plot.plots={'normFn':('fx',)}
 plot.plot()
 raw_input()
+print 'Type Return to go ahead'
 #pylab.show() #to pause on the plot window. Effective only first time
 
 print ''
@@ -116,7 +119,8 @@
 upperSphere.state.angVel=Vector3(0,0,1)
 upperSphere.state.vel=Vector3(0,0,0)
 i=O.interactions[1,0]
-O.engines=O.engines[:4]+[NewtonIntegrator()]+O.engines[5:]#+[PyRunner(iterPeriod=1,command='printInfo()')]
+
+O.engines=O.engines[:3]+[NewtonIntegrator()]+O.engines[4:]#+[PyRunner(iterPeriod=1,command='printInfo()')]
 
 
 def printInfo():
@@ -129,5 +133,6 @@
 
 plot.plots={'step':('torque',)}
 plot.plot()
+print 'Type Return to go ahead'
 
 #-- Comments : TO DO