yade-dev team mailing list archive

[noreply@xxxxxxxxxxxxx]

------------------------------------------------------------
revno: 2584
committer: help
branch nick: yade
timestamp: Thu 2010-12-02 19:27:56 +0000
message:
  Modify CohFrict* classes as discussed with Bruno on the list.
modified:
  pkg/dem/CohFrictMat.hpp
  pkg/dem/CohFrictPhys.hpp
  pkg/dem/CohesiveFrictionalContactLaw.cpp
  pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.cpp
  pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.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 'pkg/dem/CohFrictMat.hpp'
--- pkg/dem/CohFrictMat.hpp	2010-11-07 11:46:20 +0000
+++ pkg/dem/CohFrictMat.hpp	2010-12-02 19:27:56 +0000
@@ -20,7 +20,6 @@
 
 /// Serialization
 	YADE_CLASS_BASE_DOC_ATTRS_CTOR(CohFrictMat,FrictMat,"",
-		((bool,isBroken,true,,""))
 		((bool,isCohesive,true,,""))
 		((Real,normalCohesion,10000000,,""))
 		((Real,shearCohesion,10000000,,""))

=== modified file 'pkg/dem/CohFrictPhys.hpp'
--- pkg/dem/CohFrictPhys.hpp	2010-11-26 15:59:45 +0000
+++ pkg/dem/CohFrictPhys.hpp	2010-12-02 19:27:56 +0000
@@ -22,6 +22,8 @@
 		((bool,cohesionBroken,true,,"is cohesion active? will be set false when a fragile contact is broken"))
 		((bool,fragile,true,,"do cohesion disapear when contact strength is exceeded?"))
 		((Real,kr,0,,"rotational stiffness [N.m/rad]"))
+		((Real,ktw,0,,"twist stiffness [N.m/rad]"))
+		((Real,maxRollPl,0.0,,"Coefficient to determine the maximum plastic rolling moment."))
 		((Real,normalAdhesion,0,,"tensile strength"))
 		((Real,shearAdhesion,0,,"cohesive part of the shear strength (a frictional term might be added depending on :yref:`Law2_ScGeom6D_CohFrictPhys_CohesionMoment::always_use_moment_law`)"))
 		((Real,unp,0,,"plastic normal displacement, only used for tensile behaviour and if :yref:`CohFrictPhys::fragile`=false."))

=== modified file 'pkg/dem/CohesiveFrictionalContactLaw.cpp'
--- pkg/dem/CohesiveFrictionalContactLaw.cpp	2010-11-26 15:59:45 +0000
+++ pkg/dem/CohesiveFrictionalContactLaw.cpp	2010-12-02 19:27:56 +0000
@@ -93,11 +93,20 @@
 				Quaternionr q_twist_delta(q_twist_creeped * q_twist.conjugate());
 				currentContactGeometry->twistCreep = currentContactGeometry->twistCreep * q_twist_delta;
 			}
-			currentContactPhysics->moment_twist = (currentContactGeometry->getTwist()*currentContactPhysics->kr)*currentContactGeometry->normal;
+			currentContactPhysics->moment_twist = (currentContactGeometry->getTwist()*currentContactPhysics->ktw)*currentContactGeometry->normal;
 			currentContactPhysics->moment_bending = currentContactGeometry->getBending() * currentContactPhysics->kr;
+			
+			// limit rolling moment to the plastic value, if required
+			Real RollMax = currentContactPhysics->maxRollPl*currentContactPhysics->normalForce.norm();
+			if (RollMax>0.){ // apply plasticity
+				Real scalarRoll = currentContactPhysics->moment_bending.norm();		
+				Real ratio = RollMax/scalarRoll;
+				currentContactPhysics->moment_bending *= ratio;		
+			}
+			
 			Vector3r moment = currentContactPhysics->moment_twist + currentContactPhysics->moment_bending;
 			scene->forces.addTorque(id1,-moment);
-			scene->forces.addTorque(id2, moment);
+			scene->forces.addTorque(id2, moment);			
 		}
 		/// Moment law END       ///
 	}

=== modified file 'pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.cpp'
--- pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.cpp	2010-11-07 11:46:20 +0000
+++ pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.cpp	2010-12-02 19:27:56 +0000
@@ -39,6 +39,7 @@
 			Real Vb 	= sdec2->poisson;
 			Real Da 	= geom->radius1;
 			Real Db 	= geom->radius2;
+			Real Dmean  = (Da+Db)/2.; // mean radius
 			Real fa 	= sdec1->frictionAngle;
 			Real fb 	= sdec2->frictionAngle;
 			Real Kn = 2.0*Ea*Da*Eb*Db/(Ea*Da+Eb*Db);//harmonic average of two stiffnesses
@@ -50,7 +51,9 @@
 			// Jean-Patrick Plassiard, Noura Belhaine, Frederic
 			// Victor Donze, "Calibration procedure for spherical
 			// discrete elements using a local moemnt law".
-			Real Kr = Da*Db*Ks*2.0; // just like "2.0" above - it's an arbitrary parameter
+			Real Kr = pow(Dmean,2)*Ks*alphaKr;
+			Real Ktw = pow(Dmean,2)*Ks*alphaKtw;
+			
 			contactPhysics->frictionAngle			= std::min(fa,fb);
 			contactPhysics->tangensOfFrictionAngle		= std::tan(contactPhysics->frictionAngle);
 
@@ -63,7 +66,11 @@
 			}
 			contactPhysics->kn = Kn;
 			contactPhysics->ks = Ks;
+			
 			contactPhysics->kr = Kr;
+			contactPhysics->ktw = Ktw;
+			contactPhysics->maxRollPl = etaRoll*Dmean;
+			
 			contactPhysics->momentRotationLaw=(sdec1->momentRotationLaw && sdec2->momentRotationLaw);
 			//contactPhysics->elasticRollingLimit = elasticRollingLimit;
 

=== modified file 'pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.hpp'
--- pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.hpp	2010-11-12 18:44:15 +0000
+++ pkg/dem/Ip2_2xCohFrictMat_CohFrictPhys.hpp	2010-12-02 19:27:56 +0000
@@ -22,6 +22,9 @@
 		"Generates cohesive-frictional interactions with moments. Used in the contact law :yref:`Law2_ScGeom6D_CohFrictPhys_CohesionMoment`.",
 		((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."))
+		((Real,alphaKr,0.0,,"Dimensionless coefficient used for the rolling stiffness."))
+		((Real,alphaKtw,0.0,,"Dimensionless coefficient used for the twist stiffness."))
+		((Real,etaRoll,-1.,,"Dimensionless coefficient used to calculate the plastic rolling moment (if negative, plasticity will not be applied)."))
 		,
 		cohesionDefinitionIteration = -1; 
 		);