yade-dev team mailing list archive

[gladky_anton@xxxxxxxxxxxxxxx]

Author: gladky_anton
Date: 2009-07-09 17:24:23 +0200 (Thu, 09 Jul 2009)
New Revision: 1850

Modified:
   trunk/pkg/dem/RockPM.cpp
   trunk/pkg/dem/RockPM.hpp
Log:
1. Test
2. RockPM.* updated



Modified: trunk/pkg/dem/RockPM.cpp
===================================================================
--- trunk/pkg/dem/RockPM.cpp	2009-07-08 16:03:12 UTC (rev 1849)
+++ trunk/pkg/dem/RockPM.cpp	2009-07-09 15:24:23 UTC (rev 1850)
@@ -16,12 +16,11 @@
 **************************************************************************/
 
 #include"RockPM.hpp"
-#include"yade/pkg-dem/ElasticContactLaw.hpp"
+#include<yade/core/MetaBody.hpp>
 #include<yade/pkg-dem/BodyMacroParameters.hpp>
-#include<yade/pkg-dem/ElasticContactInteraction.hpp>
+#include<yade/pkg-common/Sphere.hpp>
 #include<yade/pkg-dem/DemXDofGeom.hpp>
-#include<yade/core/Omega.hpp>
-#include<yade/core/MetaBody.hpp>
+#include<yade/extra/Shop.hpp>
 
 
 YADE_PLUGIN("Law2_Dem3DofGeom_RockPMPhys_Rpm", "RpmMat","Ip2_RpmMat_RpmMat_RpmPhys","RpmPhys");
@@ -37,44 +36,27 @@
 	
 	Real displN=geom->displacementN();
 	
-	if(displN>0){rootBody->interactions->requestErase(contact->getId1(),contact->getId2()); return; }
+	if((displN<=0)){
+/*Normal Interaction*/
+		phys->normalForce=phys->kn*displN*geom->normal;
+		Real maxFsSq=phys->normalForce.SquaredLength()*pow(phys->tanFrictionAngle,2);
+		Vector3r trialFs=phys->ks*geom->displacementT();
+		if(trialFs.SquaredLength()>maxFsSq){ geom->slipToDisplacementTMax(sqrt(maxFsSq)); trialFs*=sqrt(maxFsSq/(trialFs.SquaredLength()));} 
+		applyForceAtContactPoint(phys->normalForce+trialFs,geom->contactPoint,contact->getId1(),geom->se31.position,contact->getId2(),geom->se32.position,rootBody);
+/*Normal Interaction_____*/
+		return;
+	} else if ((displN>0)&&(phys->isCohesive)){
+		if (displN<(phys->lengthMaxTension)) {
+			rootBody->interactions->requestErase(contact->getId1(),contact->getId2()); return;
+		}
+		rootBody->interactions->requestErase(contact->getId1(),contact->getId2()); return;
+	} else {
+		rootBody->interactions->requestErase(contact->getId1(),contact->getId2()); return;
+	}
 	
-	phys->normalForce=phys->kn*displN*geom->normal;
-	
-	Real maxFsSq=phys->normalForce.SquaredLength()*pow(phys->tangensOfFrictionAngle,2);
-	
-	Vector3r trialFs=phys->ks*geom->displacementT();
-	
-	if(trialFs.SquaredLength()>maxFsSq){ geom->slipToDisplacementTMax(sqrt(maxFsSq)); trialFs*=sqrt(maxFsSq/(trialFs.SquaredLength()));} 
-	
-	applyForceAtContactPoint(phys->normalForce+trialFs,geom->contactPoint,contact->getId1(),geom->se31.position,contact->getId2(),geom->se32.position,rootBody);
-	
 }
 
-/*
-void Law2_Dem3DofGeom_SimplePhys_Simple::go(shared_ptr<InteractionGeometry>& ig, shared_ptr<InteractionPhysics>& ip, Interaction* contact, MetaBody* rootBody){
-	Dem3DofGeom* geom=static_cast<Dem3DofGeom*>(ig.get());
-	
-	RpmPhys* phys=static_cast<RpmPhys*>(ip.get());
-	
-	Real displN=geom->displacementN();
-	
-	if(displN>0){rootBody->interactions->requestErase(contact->getId1(),contact->getId2()); return; }
-	
-	phys->normalForce=phys->kn*displN*geom->normal;
-	
-	Real maxFsSq=phys->normalForce.SquaredLength()*pow(phys->tangensOfFrictionAngle,2);
-	
-	Vector3r trialFs=phys->ks*geom->displacementT();
-	
-	if(trialFs.SquaredLength()>maxFsSq){ geom->slipToDisplacementTMax(sqrt(maxFsSq)); trialFs*=sqrt(maxFsSq/(trialFs.SquaredLength()));} 
-	
-	applyForceAtContactPoint(phys->normalForce+trialFs,geom->contactPoint,contact->getId1(),geom->se31.position,contact->getId2(),geom->se32.position,rootBody);
-	
-}
-*/
 
-
 CREATE_LOGGER(Ip2_RpmMat_RpmMat_RpmPhys);
 
 void Ip2_RpmMat_RpmMat_RpmPhys::go(const shared_ptr<PhysicalParameters>& pp1, const shared_ptr<PhysicalParameters>& pp2, const shared_ptr<Interaction>& interaction){
@@ -82,7 +64,10 @@
 
 	Dem3DofGeom* contGeom=YADE_CAST<Dem3DofGeom*>(interaction->interactionGeometry.get());
 	assert(contGeom);
-
+	
+	const shared_ptr<RpmMat>& rpm1=YADE_PTR_CAST<RpmMat>(pp1);
+	const shared_ptr<RpmMat>& rpm2=YADE_PTR_CAST<RpmMat>(pp2);
+	
 	const shared_ptr<BodyMacroParameters>& elast1=static_pointer_cast<BodyMacroParameters>(pp1);
 	const shared_ptr<BodyMacroParameters>& elast2=static_pointer_cast<BodyMacroParameters>(pp2);
 
@@ -96,6 +81,15 @@
 	contPhys->crossSection=S12;
 	contPhys->kn=contPhys->E*contPhys->crossSection;
 	contPhys->ks=contPhys->G*contPhys->crossSection;
+	
+	contPhys->lengthMaxCompression=S12*(rpm2->stressCompressMax)/(contPhys->kn);
+	contPhys->lengthMaxTension=S12*(rpm2->stressTensionMax)/(contPhys->kn);
+	
+	//If 2 bodies were not damaged previously, and they belong to 1 example (stone), they can be cohesive
+	if ((rpm1->exampleNumber==rpm2->exampleNumber)&&(!(rpm1->isDamaged||rpm2->isDamaged))) {
+		contPhys->isCohesive=true;
+	}
+	
 	interaction->interactionPhysics=contPhys;
 }
 

Modified: trunk/pkg/dem/RockPM.hpp
===================================================================
--- trunk/pkg/dem/RockPM.hpp	2009-07-08 16:03:12 UTC (rev 1849)
+++ trunk/pkg/dem/RockPM.hpp	2009-07-09 15:24:23 UTC (rev 1850)
@@ -22,6 +22,8 @@
 mechanical behavior of mining rocks.
 */
 
+#pragma once
+
 #include<yade/pkg-common/RigidBodyParameters.hpp>
 #include<yade/pkg-dem/BodyMacroParameters.hpp>
 #include<yade/pkg-common/InteractionPhysicsEngineUnit.hpp>
@@ -43,9 +45,12 @@
 /* This class holds information associated with each body */
 class RpmMat: public BodyMacroParameters {
 	public:
-		Real normYungConn, normInitLengthConn;
-		RpmMat(): normYungConn(0.) {createIndex();};
-		REGISTER_ATTRIBUTES(BodyMacroParameters, (normYungConn));
+		int exampleNumber; //Number of "stone"
+		bool isDamaged;
+		Real stressCompressMax, stressTensionMax; //Parameters for damaging
+
+		RpmMat(): exampleNumber(0.), isDamaged(false), stressCompressMax(0), stressTensionMax(0) {createIndex();};
+		REGISTER_ATTRIBUTES(BodyMacroParameters, (exampleNumber) (isDamaged) (stressCompressMax) (stressTensionMax));
 		REGISTER_CLASS_AND_BASE(RpmMat,BodyMacroParameters);
 };
 REGISTER_SERIALIZABLE(RpmMat);
@@ -76,16 +81,19 @@
 class RpmPhys: public NormalShearInteraction {
 	private:
 	public:
-		Real crossSection, E, G,tangensOfFrictionAngle, tanFrictionAngle;
+		Real crossSection,E,G,tanFrictionAngle,lengthMaxCompression,lengthMaxTension;
+		bool isCohesive;
 
-		RpmPhys(): NormalShearInteraction(),crossSection(0),E(0),G(0),tangensOfFrictionAngle(0), tanFrictionAngle(0) { createIndex();}
+		RpmPhys(): NormalShearInteraction(),crossSection(0),E(0),G(0), tanFrictionAngle(0), lengthMaxCompression(0), lengthMaxTension(0), isCohesive(false) { createIndex();}
 		virtual ~RpmPhys();
 
 		REGISTER_ATTRIBUTES(NormalShearInteraction,
 			(E)
 			(G)
-			(tangensOfFrictionAngle)
 			(tanFrictionAngle)
+			(isCohesive)
+			(lengthMaxCompression)
+			(lengthMaxTension)
 		);
 	REGISTER_CLASS_AND_BASE(RpmPhys,NormalShearInteraction);
 };