yade-dev team mailing list archive

[noreply@xxxxxxxxxxxxx]

------------------------------------------------------------
revno: 2512
committer: bchareyre <bchareyre@dt-rv020>
branch nick: yade
timestamp: Tue 2010-10-26 16:09:37 +0200
message:
  - Fix cylinder reg. test with respect to r2509(2510).
  - Add output for nan quaternion products #ifdef Q_DEBUG.
  - Document some flags in Ip2_2xCohFrictMat_CohFrictPhys.hpp. 
modified:
  pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.hpp
  pkg/dem/ScGeom.cpp
  py/tests/cohesive-chain.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_2xCohFrictMat_CohFrictPhys.hpp'
--- pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.hpp	2010-10-25 14:52:21 +0000
+++ pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.hpp	2010-10-26 14:09:37 +0000
@@ -20,8 +20,8 @@
 		
 		YADE_CLASS_BASE_DOC_ATTRS_CTOR(Ip2_2xCohFrictMat_CohFrictPhys,IPhysFunctor,
 		"Generates cohesive-frictional interactions with moments. Used in the contact law :yref:`Law2_ScGeom_CohFrictPhys_CohesionMoment`.",
-		((bool,setCohesionNow,false,,""))
-		((bool,setCohesionOnNewContacts,false,,""))
+		((bool,setCohesionNow,false,,"If true, assign cohesion to all existing contacts in current time-step. The flag is turned false automatically, so that assignment is done in the current timestep only."))
+		((bool,setCohesionOnNewContacts,false,,"If true, assign cohesion at all new contacts. If false, only existing contacts can be cohesive (also see :yref:`Ip2_2xCohFrictMat_CohFrictPhys::setCohesionNow`), and new contacts are only frictional."))
 		,
 		cohesionDefinitionIteration = -1; 
 		);

=== modified file 'pkg/dem/ScGeom.cpp'
--- pkg/dem/ScGeom.cpp	2010-10-26 13:41:30 +0000
+++ pkg/dem/ScGeom.cpp	2010-10-26 14:09:37 +0000
@@ -94,7 +94,16 @@
 		if (creep) delta = delta * twistCreep;
 		AngleAxisr aa(delta); // axis of rotation - this is the Moment direction UNIT vector; // angle represents the power of resistant ELASTIC moment
 		//Eigen::AngleAxisr(q) returns nan's when q close to identity, next tline fixes the pb.
+// #define Q_DEBUG
+#ifdef Q_DEBUG
+		if (isnan(aa.angle())) {
+			cerr<<"NaN found after quaternion product"<<endl;
+			cerr<<"rbp1.ori * (initialOrientation1.conjugate())) * (initialOrientation2 * (rbp2.ori.conjugate()) is:"<<endl;
+			cerr<<rbp1.ori<<" * "<<initialOrientation1.conjugate()<<" * "<<initialOrientation2<<" * "<<rbp2.ori.conjugate()<<endl<<"with sub-products :"<<endl<<rbp1.ori * (initialOrientation1.conjugate())<<" * "<<initialOrientation2 * (rbp2.ori.conjugate())<<endl;
+		}
+#else
 		if (isnan(aa.angle())) aa.angle()=0;
+#endif
 		if (aa.angle() > Mathr::PI) aa.angle() -= Mathr::TWO_PI;   // angle is between 0 and 2*pi, but should be between -pi and pi
 		twist = (aa.angle() * aa.axis().dot(normal));
 		bending = Vector3r(aa.angle()*aa.axis() - twist*normal);

=== modified file 'py/tests/cohesive-chain.py'
--- py/tests/cohesive-chain.py	2010-10-20 15:38:06 +0000
+++ py/tests/cohesive-chain.py	2010-10-26 14:09:37 +0000
@@ -22,7 +22,7 @@
 		poisson=5
 		density=2.60e3 
 		frictionAngle=radians(30)
-		O.materials.append(CohFrictMat(young=young,poisson=poisson,density=density,frictionAngle=frictionAngle))
+		O.materials.append(CohFrictMat(young=young,poisson=poisson,density=density,frictionAngle=frictionAngle,normalCohesion=1e13,shearCohesion=1e13,momentRotationLaw=True))
 		O.dt=1e-3
 		O.engines=[
 			ForceResetter(),
@@ -31,8 +31,8 @@
 			Bo1_Sphere_Aabb()]),
 		InteractionLoop(
 			[Ig2_ChainedCylinder_ChainedCylinder_ScGeom(),Ig2_Sphere_ChainedCylinder_CylScGeom()],
-			[Ip2_2xCohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=True,normalCohesion=1e13,shearCohesion=1e13)],
-			[Law2_ScGeom_CohFrictPhys_CohesionMoment(momentRotationLaw=True)]),
+			[Ip2_2xCohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=True)],
+			[Law2_ScGeom_CohFrictPhys_CohesionMoment()]),
 		## Apply gravity
 		GravityEngine(gravity=[0,-9.81,0]),
 		## Motion equation