yade-dev team mailing list archive

[noreply@xxxxxxxxxxxxx]

------------------------------------------------------------
revno: 2571
committer: Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
branch nick: yade
timestamp: Thu 2010-11-25 13:16:32 +0100
message:
  - FlowEngine : rewrite updateVolume so that cell volumes are not computed twice.
  - Remove irrelevant remark in SCPE
  - Some missing "6D" here and ther in decorations (could affect only devirtualized build) 
modified:
  pkg/common/Cylinder.cpp
  pkg/common/Cylinder.hpp
  pkg/dem/FlowEngine.cpp
  pkg/dem/Ig2_Sphere_Sphere_ScGeom.cpp
  pkg/dem/SampleCapillaryPressureEngine.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/common/Cylinder.cpp'
--- pkg/common/Cylinder.cpp	2010-11-12 18:44:15 +0000
+++ pkg/common/Cylinder.cpp	2010-11-25 12:16:32 +0000
@@ -9,10 +9,11 @@
 ChainedCylinder::~ChainedCylinder(){}
 ChainedState::~ChainedState(){}
 CylScGeom::~CylScGeom(){}
-
-
-YADE_PLUGIN(
-	(Cylinder)(ChainedCylinder)(ChainedState)(CylScGeom)(Ig2_Sphere_ChainedCylinder_CylScGeom)(Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D)
+// Ig2_Sphere_ChainedCylinder_CylScGeom::~Ig2_Sphere_ChainedCylinder_CylScGeom() {}
+// Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D::~Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D() {}
+
+
+YADE_PLUGIN((Cylinder)(ChainedCylinder)(ChainedState)(CylScGeom)(Ig2_Sphere_ChainedCylinder_CylScGeom)(Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D)
 	#ifdef YADE_OPENGL
 		(Gl1_Cylinder)(Gl1_ChainedCylinder)
 	#endif
@@ -64,7 +65,7 @@
 	Real dist = direction.dot(branch);
 // 	if ((dist<-interactionDetectionFactor*cylinder->radius) && !c->isReal()) return (false);//position _before_ start of cylinder
 	if ((dist<-interactionDetectionFactor*cylinder->radius) && isNew) return (false);//position _before_ start of cylinder
-		
+
 	//Check sphere-cylinder distance
 	Vector3r projectedP = cylinderSt->pos+shift2 + direction*dist;
 	branch = projectedP-sphereSt->pos;
@@ -74,7 +75,7 @@
 	shared_ptr<CylScGeom> scm;
 	if(!isNew) scm=YADE_PTR_CAST<CylScGeom>(c->geom);
 	else { scm=shared_ptr<CylScGeom>(new CylScGeom()); c->geom=scm; }
-	
+
 	scm->radius1=sphere->radius;
 	scm->radius2=cylinder->radius;
 	scm->id3=cylinderSt->chains[cylinderSt->chainNumber][cylinderSt->rank+1];
@@ -160,9 +161,9 @@
 	const ChainedState& bchain2 = revert? *pChain1 : *pChain2;
 	if (bchain2.rank-bchain1.rank != 1) {/*cerr<<"Mutual contacts in same chain between not adjacent elements, not handled*/ return false;}
 	if (pChain2->chainNumber!=pChain1->chainNumber) {cerr<<"PROBLEM0124"<<endl; return false;}
-	
+
 	ChainedCylinder *bs1=static_cast<ChainedCylinder*>(revert? cm2.get():cm1.get());
-	
+
 	shared_ptr<ScGeom6D> scm;
 	bool isNew = !c->geom;
 	if(!isNew) scm=YADE_PTR_CAST<ScGeom6D>(c->geom);
@@ -176,8 +177,8 @@
 	scm->contactPoint=bchain2.pos;
 	//bs1->segment used for fast BBs and projections + display
 	bs1->segment= bchain2.pos-bchain1.pos;
-#ifdef YADE_OPENGL 
-	//bs1->length and s1->chainedOrientation used for display only, 
+#ifdef YADE_OPENGL
+	//bs1->length and s1->chainedOrientation used for display only,
 	bs1->length=length;
 	bs1->chainedOrientation.setFromTwoVectors(Vector3r::UnitZ(),bchain1.ori.conjugate()*segt);
 #endif
@@ -301,10 +302,10 @@
 //    gluDeleteQuadric(quadObj);
 //    glPopMatrix();
 // //    GLERROR;
-// 	
+//
 // // 	glEndList();}
 // // 	glCallList(glCylinderList);
-// 
+//
 // }
 
 //!##################	BOUNDS FUNCTOR   #####################

=== modified file 'pkg/common/Cylinder.hpp'
--- pkg/common/Cylinder.hpp	2010-11-12 18:44:15 +0000
+++ pkg/common/Cylinder.hpp	2010-11-25 12:16:32 +0000
@@ -1,9 +1,9 @@
 #pragma once
+#include<yade/pkg/common/Dispatching.hpp>
 #include<yade/core/Shape.hpp>
 #include<yade/core/State.hpp>
 #include<yade/core/Body.hpp>
 #include<yade/pkg/dem/ScGeom.hpp>
-#include<yade/pkg/common/Dispatching.hpp>
 #include<yade/pkg/common/Sphere.hpp>
 #include<yade/core/Scene.hpp>
 #ifdef YADE_OPENGL
@@ -13,28 +13,31 @@
 
 class Cylinder: public Sphere{
 	public:
-		Cylinder(Real _radius, Real _length): length(_length) { /*segment=Vector3r(0,0,1)*_length;*/ radius=_radius; createIndex(); }
+// 		Cylinder(Real _radius, Real _length): length(_length) { /*segment=Vector3r(0,0,1)*_length;*/ radius=_radius; createIndex(); }
 		virtual ~Cylinder ();
 	YADE_CLASS_BASE_DOC_ATTRS_CTOR(Cylinder,Sphere,"Geometry of a cylinder, as Minkowski sum of line and sphere.",
 // 		((Real,radius,NaN,,"Radius [m]"))
 		((Real,length,NaN,,"Length [m]"))
 		((Vector3r,segment,Vector3r::Zero(),,"Length vector")),
 		createIndex();
+
 		/*ctor*/
+		segment=Vector3r(0,0,1)*length;
 	);
 	REGISTER_CLASS_INDEX(Cylinder,Sphere);
 };
+REGISTER_SERIALIZABLE(Cylinder);
 
 class ChainedCylinder: public Cylinder{
 	public:
-		ChainedCylinder(Real _radius, Real _length): Cylinder(_radius,_length){}
+// 		ChainedCylinder(Real _radius, Real _length);/*: Cylinder(_radius,_length){}*/
 		virtual ~ChainedCylinder ();
-		 
+
 		//Keep pointers or copies of connected states?
 // 		ChainedState st1, st2;
-		
-
-	YADE_CLASS_BASE_DOC_ATTRS_CTOR(ChainedCylinder,Cylinder,"Geometry of a deformable chained cylinder, using geometry :yref:`MinkCylinder`.",
+
+
+	YADE_CLASS_BASE_DOC_ATTRS_CTOR(ChainedCylinder,Cylinder,"Geometry of a deformable chained cylinder, using geometry :yref:`Cylinder`.",
   		((Real,initLength,0,,"tensile-free length, used as reference for tensile strain"))
   		((Quaternionr,chainedOrientation,Quaternionr::Identity(),,"Deviation of node1 orientation from node-to-node vector"))
 		,createIndex();/*ctor*/
@@ -43,6 +46,7 @@
 	);
 	REGISTER_CLASS_INDEX(ChainedCylinder,Cylinder);
 };
+REGISTER_SERIALIZABLE(ChainedCylinder);
 
 class CylScGeom: public ScGeom{
 	public:
@@ -57,6 +61,7 @@
 	);
 	REGISTER_CLASS_INDEX(CylScGeom,ScGeom);
 };
+REGISTER_SERIALIZABLE(CylScGeom);
 
 
 class ChainedState: public State{
@@ -87,12 +92,14 @@
 		.def_readwrite("currentChain",&ChainedState::currentChain,"Current active chain (where newly created chained bodies will be appended).")
 		.def("addToChain",&ChainedState::addToChain,(python::arg("bodyId")),"Add body to current active chain")
 	);
+	REGISTER_CLASS_INDEX(ChainedState,State);
 };
-
+REGISTER_SERIALIZABLE(ChainedState);
 
 
 class Ig2_Sphere_ChainedCylinder_CylScGeom: public IGeomFunctor{
 	public:
+// 		virtual ~Ig2_Sphere_ChainedCylinder_CylScGeom ();
 		virtual bool go(const shared_ptr<Shape>& cm1, const shared_ptr<Shape>& cm2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
 		virtual bool goReverse(	const shared_ptr<Shape>& cm1, const shared_ptr<Shape>& cm2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
 	#ifdef YADE_DEVIRT_FUNCTORS
@@ -105,9 +112,11 @@
 	FUNCTOR2D(Sphere,ChainedCylinder);
 	DEFINE_FUNCTOR_ORDER_2D(Sphere,ChainedCylinder);
 };
+REGISTER_SERIALIZABLE(Ig2_Sphere_ChainedCylinder_CylScGeom);
 
 class Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D: public IGeomFunctor{
 	public:
+// 		virtual ~Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D ()  {};
 		virtual bool go(const shared_ptr<Shape>& cm1, const shared_ptr<Shape>& cm2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
 		virtual bool goReverse(	const shared_ptr<Shape>& cm1, const shared_ptr<Shape>& cm2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
 	#ifdef YADE_DEVIRT_FUNCTORS
@@ -121,7 +130,7 @@
 	// needed for the dispatcher, even if it is symmetric
 	DEFINE_FUNCTOR_ORDER_2D(ChainedCylinder,ChainedCylinder);
 };
-
+REGISTER_SERIALIZABLE(Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D);
 
 
 #ifdef YADE_OPENGL
@@ -184,6 +193,7 @@
 		((Real,aabbEnlargeFactor,((void)"deactivated",-1),,"Relative enlargement of the bounding box; deactivated if negative.\n\n.. note::\n\tThis attribute is used to create distant interaction, but is only meaningful with an :yref:`IGeomFunctor` which will not simply discard such interactions: :yref:`Ig2_Cylinder_Cylinder_Dem3DofGeom::distFactor` / :yref:`Ig2_Cylinder_Cylinder_ScGeom::interactionDetectionFactor` should have the same value as :yref:`aabbEnlargeFactor<Bo1_Cylinder_Aabb::aabbEnlargeFactor>`."))
 	);
 };
+REGISTER_SERIALIZABLE(Bo1_Cylinder_Aabb);
 
 class Bo1_ChainedCylinder_Aabb : public BoundFunctor
 {
@@ -194,6 +204,7 @@
 		((Real,aabbEnlargeFactor,((void)"deactivated",-1),,"Relative enlargement of the bounding box; deactivated if negative.\n\n.. note::\n\tThis attribute is used to create distant interaction, but is only meaningful with an :yref:`IGeomFunctor` which will not simply discard such interactions: :yref:`Ig2_Cylinder_Cylinder_Dem3DofGeom::distFactor` / :yref:`Ig2_Cylinder_Cylinder_ScGeom::interactionDetectionFactor` should have the same value as :yref:`aabbEnlargeFactor<Bo1_Cylinder_Aabb::aabbEnlargeFactor>`."))
 	);
 };
+REGISTER_SERIALIZABLE(Bo1_ChainedCylinder_Aabb);
 
 // Keep this : Cylinders and ChainedCylinders will have different centers maybe.
 // class Bo1_ChainedCylinder_Aabb : public Bo1_Cylinder_Aabb
@@ -206,14 +217,6 @@
 // 	);
 // };
 
-REGISTER_SERIALIZABLE(Bo1_Cylinder_Aabb);
-REGISTER_SERIALIZABLE(Bo1_ChainedCylinder_Aabb);
-REGISTER_SERIALIZABLE(Cylinder);
-REGISTER_SERIALIZABLE(ChainedCylinder);
-REGISTER_SERIALIZABLE(ChainedState);
-REGISTER_SERIALIZABLE(CylScGeom);
-REGISTER_SERIALIZABLE(Ig2_Sphere_ChainedCylinder_CylScGeom);
-REGISTER_SERIALIZABLE(Ig2_ChainedCylinder_ChainedCylinder_ScGeom6D);
 #ifdef YADE_OPENGL
 REGISTER_SERIALIZABLE(Gl1_Cylinder);
 REGISTER_SERIALIZABLE(Gl1_ChainedCylinder);

=== modified file 'pkg/dem/FlowEngine.cpp'
--- pkg/dem/FlowEngine.cpp	2010-11-19 17:18:56 +0000
+++ pkg/dem/FlowEngine.cpp	2010-11-25 12:16:32 +0000
@@ -25,11 +25,13 @@
 {
 }
 
+
 void FlowEngine::action()
 {
 	if (!flow) {
 		flow = shared_ptr<FlowSolver> (new FlowSolver);
 		first=true;
+		cerr <<"first = true"<<endl;
 		Update_Triangulation=false;/*IS=0.f;*/
 		eps_vol_max=0.f;
 		Eps_Vol_Cumulative=0.f;
@@ -64,7 +66,7 @@
 			eps_vol_max=0.f;
 			UpdateVolumes ( );
 			Eps_Vol_Cumulative += eps_vol_max;
-			if (Eps_Vol_Cumulative > EpsVolPercent_RTRG || retriangulationLastIter>100) {
+			if ((Eps_Vol_Cumulative > EpsVolPercent_RTRG || retriangulationLastIter>1000) && retriangulationLastIter>10) {
 				Update_Triangulation = true;
 				Eps_Vol_Cumulative=0;
 				retriangulationLastIter=0;
@@ -355,35 +357,28 @@
 	Real deltaT = scene->dt;
 	Real invDeltaT = 1/scene->dt;
 	CGT::Finite_cells_iterator cell_end = flow->T[flow->currentTes].Triangulation().finite_cells_end();
-	for ( CGT::Finite_cells_iterator cell = flow->T[flow->currentTes].Triangulation().finite_cells_begin(); cell != cell_end; cell++ )
-	{
-		switch ( cell->info().fictious() )
-		{
+	double newVol=0; double dVol;
+	eps_vol_max=0;
+	for ( CGT::Finite_cells_iterator cell = flow->T[flow->currentTes].Triangulation().finite_cells_begin(); cell != cell_end; cell++ ) {
+		switch ( cell->info().fictious()) {
 			case ( 3 ):
-			{
-				cell->info().dv() = ( Volume_cell_triple_fictious ( cell ) - cell->info().volume() ) *invDeltaT;
-				eps_vol_max = max(eps_vol_max, (abs(cell->info().dv()*deltaT))/cell->info().volume());
-				cell->info().volume() = Volume_cell_triple_fictious ( cell );
-			}break;
+				newVol= Volume_cell_triple_fictious ( cell );
+				break;
 			case ( 2 ) :
-			{
-				cell->info().dv() = ( Volume_cell_double_fictious ( cell )-cell->info().volume() ) *invDeltaT;
-				eps_vol_max = max(eps_vol_max, (abs(cell->info().dv()*deltaT))/cell->info().volume());
-				cell->info().volume() = Volume_cell_double_fictious ( cell );
-			}break;
+				newVol = Volume_cell_double_fictious ( cell );
+				break;
 			case ( 1 ) :
-			{
-				cell->info().dv() = ( Volume_cell_single_fictious ( cell )-cell->info().volume() ) *invDeltaT;
-				eps_vol_max = max(eps_vol_max, (abs(cell->info().dv()*deltaT))/cell->info().volume());
-				cell->info().volume() = Volume_cell_single_fictious ( cell );
-			}break;
+				newVol = Volume_cell_single_fictious ( cell );
+				break;
 			case ( 0 ) :
-			{
-				cell->info().dv() = ( Volume_cell ( cell )-cell->info().volume() ) *invDeltaT;
-				eps_vol_max = max(eps_vol_max, (abs(cell->info().dv()*deltaT))/cell->info().volume());
-				cell->info().volume() = Volume_cell ( cell );
-			}break;
+				newVol = Volume_cell ( cell );
+				break;
 		}
+		dVol=newVol - cell->info().volume();
+		eps_vol_max = max(eps_vol_max, abs(dVol/newVol));
+		cell->info().dv() = dVol*invDeltaT;
+		cell->info().volume() = newVol;
+// 		if (Debug) cerr<<"v/dv : "<<cell->info().volume()<<" "<<cell->info().dv()<<" ("<<cell->info().fictious()<<")"<<endl;
 	}
 }
 

=== modified file 'pkg/dem/Ig2_Sphere_Sphere_ScGeom.cpp'
--- pkg/dem/Ig2_Sphere_Sphere_ScGeom.cpp	2010-11-07 11:46:20 +0000
+++ pkg/dem/Ig2_Sphere_Sphere_ScGeom.cpp	2010-11-25 12:16:32 +0000
@@ -60,7 +60,7 @@
 #ifdef YADE_DEVIRT_FUNCTORS
 bool Ig2_Sphere_Sphere_ScGeom6D::go(const shared_ptr<Shape>& cm1, const shared_ptr<Shape>& cm2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c){ throw runtime_error("Do not call Ig2_Sphere_Sphere_ScGeom6D::go, use getStaticFunctorPtr and call that function instead."); }
 bool Ig2_Sphere_Sphere_ScGeom6D::goStatic(IGeomFunctor* _self, const shared_ptr<Shape>& cm1, const shared_ptr<Shape>& cm2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c){
-	const Ig2_Sphere_Sphere_ScGeom6D* self=static_cast<Ig2_Sphere_Sphere_ScGeom*>(_self);
+	const Ig2_Sphere_Sphere_ScGeom6D* self=static_cast<Ig2_Sphere_Sphere_ScGeom6D*>(_self);
 	const Real& interactionDetectionFactor=self->interactionDetectionFactor;
 #else
 bool Ig2_Sphere_Sphere_ScGeom6D::go(	const shared_ptr<Shape>& cm1,

=== modified file 'pkg/dem/SampleCapillaryPressureEngine.hpp'
--- pkg/dem/SampleCapillaryPressureEngine.hpp	2010-11-07 11:46:20 +0000
+++ pkg/dem/SampleCapillaryPressureEngine.hpp	2010-11-25 12:16:32 +0000
@@ -26,15 +26,15 @@
 		std::string Phase1End;
 		//! is this the beginning of the simulation, after reading the scene?
 		bool firstRun;
-		
+
 		shared_ptr<Law2_ScGeom_CapillaryPhys_Capillarity>  capillaryCohesiveLaw;
 		//Law2_ScGeom_CapillaryPhys_Capillarity* capillaryCohesiveLaw; // which one is right?
-		
+
 		virtual ~SampleCapillaryPressureEngine();
 		void updateParameters();
 		virtual void action();
-		
-	YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(SampleCapillaryPressureEngine,TriaxialStressController,"Rk: this engine has to be tested withthe new formalism. It produces the isotropic compaction of an assembly and allows to controlled the capillary pressure inside (uses Law2_ScGeom_CapillaryPhys_Capillarity).",
+
+	YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(SampleCapillaryPressureEngine,TriaxialStressController,"It produces the isotropic compaction of an assembly and allows to controlled the capillary pressure inside (uses Law2_ScGeom_CapillaryPhys_Capillarity).",
 		((Real,Pressure,0,,"Value of the capillary pressure Uc=Ugas-Uliquid (see Law2_ScGeom_CapillaryPhys_Capillarity). [Pa]"))
 		((bool,pressureVariationActivated,1,,"Is the capillary pressure varying?"))
 		((bool,fusionDetection,1,,"Is the detection of menisci overlapping activated?"))