yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #12613
[Branch ~yade-pkg/yade/git-trunk] Rev 3839: Polyhedra modification
------------------------------------------------------------
revno: 3839
committer: Jan Stransky <jan.stransky@xxxxxxxxxxx>
timestamp: Wed 2016-04-13 18:42:10 +0200
message:
Polyhedra modification
added Bo1_Polyhedra_Aabb.aabbEnlargeFactor
added Ig2_Polyhedra_Polyhedra_PolyhedraGeom.interactionDetectionFactor
added Ig2_Polyhedra_Polyhedra_ScGeom class
added Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom class to handle both ScGeom and PolyhedraGeom between polyhedral particles
added example/test scripts
added:
examples/polyhedra/tests/
examples/polyhedra/tests/interactinDetectionFactor.py
examples/polyhedra/tests/scGeom.py
modified:
pkg/dem/Polyhedra.cpp
pkg/dem/Polyhedra.hpp
pkg/dem/Polyhedra_Ig2.cpp
pkg/dem/Polyhedra_Ig2.hpp
--
lp:yade
https://code.launchpad.net/~yade-pkg/yade/git-trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-pkg/yade/git-trunk/+edit-subscription
=== added directory 'examples/polyhedra/tests'
=== added file 'examples/polyhedra/tests/interactinDetectionFactor.py'
--- examples/polyhedra/tests/interactinDetectionFactor.py 1970-01-01 00:00:00 +0000
+++ examples/polyhedra/tests/interactinDetectionFactor.py 2016-04-13 16:42:10 +0000
@@ -0,0 +1,27 @@
+from yade import plot, polyhedra_utils
+from yade import qt
+
+m = PolyhedraMat()
+O.materials.append(m)
+
+t1,t2 = [polyhedra_utils.polyhedra(m,size=(1,1,1),seed=i) for i in (5,6)]
+O.bodies.append((t1,t2))
+t2.state.pos = (3,0,0)
+t2.state.vel = (-1,0,0)
+
+factor = 2
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Polyhedra_Aabb(aabbEnlargeFactor=factor)]),
+ InteractionLoop(
+ [Ig2_Polyhedra_Polyhedra_PolyhedraGeom(interactionDetectionFactor=factor)],
+ [Ip2_PolyhedraMat_PolyhedraMat_PolyhedraPhys()],
+ [Law2_PolyhedraGeom_PolyhedraPhys_Volumetric()]
+ ),
+ NewtonIntegrator(),
+]
+
+O.dt=0.00001
+
+try: qt.View()
+except: pass
=== added file 'examples/polyhedra/tests/scGeom.py'
--- examples/polyhedra/tests/scGeom.py 1970-01-01 00:00:00 +0000
+++ examples/polyhedra/tests/scGeom.py 2016-04-13 16:42:10 +0000
@@ -0,0 +1,28 @@
+from yade import plot, polyhedra_utils
+from yade import qt
+
+scGeom = 1
+
+m = FrictMat() if scGeom else PolyhedraMat()
+O.materials.append(m)
+
+t1,t2 = [polyhedron(((-1,-1,-1),(+1,-1,-1),(-1,+1,-1),(+1,+1,-1),(-1,-1,+1),(+1,-1,+1),(-1,+1,+1),(+1,+1,+1))) for i in (0,1)]
+O.bodies.append((t1,t2))
+t2.state.pos = (3,0,0)
+t2.state.vel = (-1,0,0)
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Polyhedra_Aabb()]),
+ InteractionLoop(
+ [Ig2_Polyhedra_Polyhedra_ScGeom() if scGeom else Ig2_Polyhedra_Polyhedra_PolyhedraGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys() if scGeom else Ip2_PolyhedraMat_PolyhedraMat_PolyhedraPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack() if scGeom else Law2_PolyhedraGeom_PolyhedraPhys_Volumetric()]
+ ),
+ NewtonIntegrator(),
+]
+
+O.dt=0.00001
+
+try: qt.View()
+except: pass
=== modified file 'pkg/dem/Polyhedra.cpp'
--- pkg/dem/Polyhedra.cpp 2016-03-24 22:08:26 +0000
+++ pkg/dem/Polyhedra.cpp 2016-04-13 16:42:10 +0000
@@ -250,6 +250,10 @@
mincoords = Vector3r(min(mincoords[0],v_g[0]),min(mincoords[1],v_g[1]),min(mincoords[2],v_g[2]));
maxcoords = Vector3r(max(maxcoords[0],v_g[0]),max(maxcoords[1],v_g[1]),max(maxcoords[2],v_g[2]));
}
+ if (aabbEnlargeFactor>0) {
+ mincoords *= aabbEnlargeFactor;
+ maxcoords *= aabbEnlargeFactor;
+ }
aabb->min=se3.position+mincoords;
aabb->max=se3.position+maxcoords;
}
=== modified file 'pkg/dem/Polyhedra.hpp'
--- pkg/dem/Polyhedra.hpp 2016-03-24 22:08:26 +0000
+++ pkg/dem/Polyhedra.hpp 2016-04-13 16:42:10 +0000
@@ -153,7 +153,9 @@
public:
void go(const shared_ptr<Shape>& ig, shared_ptr<Bound>& bv, const Se3r& se3, const Body*);
FUNCTOR1D(Polyhedra);
- YADE_CLASS_BASE_DOC(Bo1_Polyhedra_Aabb,BoundFunctor,"Create/update :yref:`Aabb` of a :yref:`Polyhedra`");
+ YADE_CLASS_BASE_DOC_ATTRS(Bo1_Polyhedra_Aabb,BoundFunctor,"Create/update :yref:`Aabb` of a :yref:`Polyhedra`",
+ ((Real,aabbEnlargeFactor,((void)"deactivated",-1),,"see :yref:`Bo1_Sphere_Aabb.aabbEnlargeFactor`"))
+ );
};
REGISTER_SERIALIZABLE(Bo1_Polyhedra_Aabb);
=== modified file 'pkg/dem/Polyhedra_Ig2.cpp'
--- pkg/dem/Polyhedra_Ig2.cpp 2016-04-12 04:37:00 +0000
+++ pkg/dem/Polyhedra_Ig2.cpp 2016-04-13 16:42:10 +0000
@@ -5,7 +5,9 @@
#include "Polyhedra_Ig2.hpp"
-YADE_PLUGIN(/* self-contained in hpp: */ (Ig2_Polyhedra_Polyhedra_PolyhedraGeom) (Ig2_Wall_Polyhedra_PolyhedraGeom) (Ig2_Facet_Polyhedra_PolyhedraGeom) (Ig2_Sphere_Polyhedra_ScGeom) );
+YADE_PLUGIN(/* self-contained in hpp: */ (Ig2_Polyhedra_Polyhedra_PolyhedraGeom) (Ig2_Wall_Polyhedra_PolyhedraGeom) (Ig2_Facet_Polyhedra_PolyhedraGeom) (Ig2_Sphere_Polyhedra_ScGeom)
+ (Ig2_Polyhedra_Polyhedra_ScGeom) (Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom)
+);
//**********************************************************************************
/*! Create Polyhedra (collision geometry) from colliding Polyhedras. */
@@ -29,11 +31,13 @@
//move and rotate 1st the CGAL structure Polyhedron
Matrix3r rot_mat = (se31.orientation).toRotationMatrix();
Vector3r trans_vec = se31.position;
+ const Real& s = interactionDetectionFactor;
Transformation t_rot_trans(
- rot_mat(0,0), rot_mat(0,1), rot_mat(0,2), trans_vec[0],
- rot_mat(1,0), rot_mat(1,1), rot_mat(1,2), trans_vec[1],
- rot_mat(2,0), rot_mat(2,1), rot_mat(2,2), trans_vec[2], 1.
+ s*rot_mat(0,0), s*rot_mat(0,1), s*rot_mat(0,2), trans_vec[0],
+ s*rot_mat(1,0), s*rot_mat(1,1), s*rot_mat(1,2), trans_vec[1],
+ s*rot_mat(2,0), s*rot_mat(2,1), s*rot_mat(2,2), trans_vec[2],
+ 1.
);
Polyhedron PA = A->GetPolyhedron();
@@ -44,7 +48,12 @@
//move and rotate 2nd the CGAL structure Polyhedron
rot_mat = (se32.orientation).toRotationMatrix();
trans_vec = se32.position + shift2;
- t_rot_trans = Transformation(rot_mat(0,0),rot_mat(0,1),rot_mat(0,2), trans_vec[0],rot_mat(1,0),rot_mat(1,1),rot_mat(1,2),trans_vec[1],rot_mat(2,0),rot_mat(2,1),rot_mat(2,2),trans_vec[2],1.);
+ t_rot_trans = Transformation(
+ s*rot_mat(0,0), s*rot_mat(0,1), s*rot_mat(0,2), trans_vec[0],
+ s*rot_mat(1,0), s*rot_mat(1,1), s*rot_mat(1,2), trans_vec[1],
+ s*rot_mat(2,0), s*rot_mat(2,1), s*rot_mat(2,2), trans_vec[2],
+ 1.
+ );
Polyhedron PB = B->GetPolyhedron();
@@ -59,7 +68,6 @@
bang->sep_plane.assign(3,0);
bang->contactPoint = Vector3r(0,0,0);
bang->isShearNew = true;
- interaction->geom = bang;
}else{
// use data from old interaction
bang=YADE_PTR_CAST<PolyhedraGeom>(interaction->geom);
@@ -79,14 +87,16 @@
bang->equivalentPenetrationDepth=0;
bang->penetrationVolume=min(A->GetVolume(),B->GetVolume());
bang->normal = (A->GetVolume()>B->GetVolume() ? 1 : -1)*(se32.position+shift2-se31.position);
- return true;
+ return !isNew;
}
if ((!Is_inside_Polyhedron(PA, ToCGALPoint(centroid))) or
(!Is_inside_Polyhedron(PB, ToCGALPoint(centroid)))) {
bang->equivalentPenetrationDepth=0;
- return true;
+ return !isNew;
}
+
+ if (isNew) interaction->geom = bang;
//find normal direction
Vector3r normal = FindNormal(Int, PA, PB);
@@ -99,6 +109,7 @@
bang->contactPoint=centroid;
bang->penetrationVolume=volume;
bang->equivalentPenetrationDepth=volume/area;
+ scene = Omega::instance().getScene().get();
bang->precompute(state1,state2,scene,interaction,normal,bang->isShearNew,shift2);
bang->normal=normal;
@@ -428,4 +439,76 @@
return true;
}
+
+
+//**********************************************************************************
+/*! Plyhedra -> ScGeom. */
+bool Ig2_Polyhedra_Polyhedra_ScGeom::go(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& interaction) {
+ const Se3r& se31=state1.se3;
+ const Se3r& se32=state2.se3;
+ shared_ptr<ScGeom> geom;
+ bool isNew = !interaction->geom;
+ if (isNew) {
+ Ig2_Polyhedra_Polyhedra_PolyhedraGeom ppGeom = Ig2_Polyhedra_Polyhedra_PolyhedraGeom();
+ ppGeom.interactionDetectionFactor = interactionDetectionFactor;
+ bool pp = ppGeom.go(shape1,shape2,state2,state1,shift2,force,interaction);
+ if (!pp) {
+ return false;
+ }
+ shared_ptr<PolyhedraGeom> pGeom = YADE_PTR_CAST<PolyhedraGeom>(interaction->geom);
+ geom = shared_ptr<ScGeom>(new ScGeom());
+ geom->radius1 = (pGeom->contactPoint-se31.position).norm();
+ geom->radius2 = (pGeom->contactPoint-se32.position+shift2).norm();
+ interaction->geom=geom;
+ } else {
+ geom = YADE_PTR_CAST<ScGeom>(interaction->geom);
+ }
+ const Real& radius1 = geom->radius1;
+ const Real& radius2 = geom->radius2;
+ Vector3r normal=(se32.position+shift2)-se31.position;
+ Real norm=normal.norm(); normal/=norm; // normal is unit vector now
+ Real penetrationDepth=radius1+radius2-norm;
+ geom->contactPoint=se31.position+(radius1-0.5*penetrationDepth)*normal;//0.5*(pt1+pt2);
+ geom->penetrationDepth=penetrationDepth;
+ scene = Omega::instance().getScene().get();
+ geom->precompute(state1,state2,scene,interaction,normal,isNew,shift2,false);
+ return true;
+}
+
+bool Ig2_Polyhedra_Polyhedra_ScGeom::goReverse(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& interaction) {
+ return go(shape1,shape2,state2,state1,-shift2,force,interaction);
+}
+
+
+
+
+bool Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom::go(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& interaction) {
+ bool isNew = !interaction->geom;
+ if (isNew) {
+ if (createScGeom) {
+ return ig2scGeom->go(shape1,shape2,state1,state2,shift2,force,interaction);
+ } else {
+ return ig2polyhedraGeom->go(shape1,shape2,state1,state2,shift2,force,interaction);
+ }
+ }
+ //
+ ScGeom* scGeom = dynamic_cast<ScGeom*>(interaction->geom.get());
+ if (scGeom) {
+ return ig2scGeom->go(shape1,shape2,state1,state2,shift2,force,interaction);
+ }
+ //
+ PolyhedraGeom* pGeom = dynamic_cast<PolyhedraGeom*>(interaction->geom.get());
+ if (pGeom) {
+ return ig2polyhedraGeom->go(shape1,shape2,state1,state2,shift2,force,interaction);
+ }
+ //
+ LOG_ERROR("TODO, should not happen");
+ return false;
+}
+
+bool Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom::goReverse(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1, const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& interaction) {
+ return go(shape1,shape2,state2,state1,-shift2,force,interaction);
+}
+
+
#endif // YADE_CGAL
=== modified file 'pkg/dem/Polyhedra_Ig2.hpp'
--- pkg/dem/Polyhedra_Ig2.hpp 2016-03-24 22:08:26 +0000
+++ pkg/dem/Polyhedra_Ig2.hpp 2016-04-13 16:42:10 +0000
@@ -16,7 +16,9 @@
const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
FUNCTOR2D(Polyhedra,Polyhedra);
DEFINE_FUNCTOR_ORDER_2D(Polyhedra,Polyhedra);
- YADE_CLASS_BASE_DOC(Ig2_Polyhedra_Polyhedra_PolyhedraGeom,IGeomFunctor,"Create/update geometry of collision between 2 Polyhedras");
+ YADE_CLASS_BASE_DOC_ATTRS(Ig2_Polyhedra_Polyhedra_PolyhedraGeom,IGeomFunctor,"Create/update geometry of collision between 2 Polyhedras",
+ ((Real,interactionDetectionFactor,1,,"see :yref:`Ig2_Sphere_Sphere_ScGeom.interactionDetectionFactor`"))
+ );
DECLARE_LOGGER;
};
REGISTER_SERIALIZABLE(Ig2_Polyhedra_Polyhedra_PolyhedraGeom);
@@ -70,5 +72,44 @@
};
REGISTER_SERIALIZABLE(Ig2_Sphere_Polyhedra_ScGeom);
+
+//***************************************************************************
+/*! Plyhedra -> ScGeom. */
+class Ig2_Polyhedra_Polyhedra_ScGeom: public IGeomFunctor
+{
+ public:
+ virtual ~Ig2_Polyhedra_Polyhedra_ScGeom(){};
+ virtual bool go(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1,
+ const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
+ virtual bool goReverse(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1,
+ const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
+ FUNCTOR2D(Polyhedra,Polyhedra);
+ DEFINE_FUNCTOR_ORDER_2D(Polyhedra,Polyhedra);
+ YADE_CLASS_BASE_DOC_ATTRS(Ig2_Polyhedra_Polyhedra_ScGeom,IGeomFunctor,"EXPERIMENTAL. Ig2 functor creating ScGeom from two Polyhedra shapes. The radii are computed as a distance of contact point (computed using Ig2_Polyhedra_Polyhedra_PolyhedraGeom) and center of particle. Tested only for face-face contacts (like brick wall).",
+ ((Real,interactionDetectionFactor,1,,"see Ig2_Sphere_Sphere_ScGeom.interactionDetectionFactor"))
+ );
+ DECLARE_LOGGER;
+};
+REGISTER_SERIALIZABLE(Ig2_Polyhedra_Polyhedra_ScGeom);
+
+class Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom: public IGeomFunctor
+{
+ public:
+ virtual ~Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom(){};
+ virtual bool go(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1,
+ const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
+ virtual bool goReverse(const shared_ptr<Shape>& shape1, const shared_ptr<Shape>& shape2, const State& state1,
+ const State& state2, const Vector3r& shift2, const bool& force, const shared_ptr<Interaction>& c);
+ FUNCTOR2D(Polyhedra,Polyhedra);
+ DEFINE_FUNCTOR_ORDER_2D(Polyhedra,Polyhedra);
+ YADE_CLASS_BASE_DOC_ATTRS(Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom,IGeomFunctor,"EXPERIMENTAL. A hacky helper Ig2 functor combining two Polyhedra shapes and creating, according to the settings, either ScGeom or PolyhedraGeom.",
+ ((bool,createScGeom,true,,"If true, creates ScGeom on new contacts. Creates PolyhedraGeom otherwise. On existing contacts Ig2_Polyhedra_Polyhedra_PolyhedraGeom or Ig2_Polyhedra_Polyhedra_ScGeom is used according to present IGeom isntance."))
+ ((shared_ptr<Ig2_Polyhedra_Polyhedra_PolyhedraGeom>,ig2polyhedraGeom,new Ig2_Polyhedra_Polyhedra_PolyhedraGeom,,"Helper Ig2 functor for PolyhedraGeom (to be able to modify its settings)"))
+ ((shared_ptr<Ig2_Polyhedra_Polyhedra_ScGeom>,ig2scGeom,new Ig2_Polyhedra_Polyhedra_ScGeom,,"Helper Ig2 functor for ScGeom (to be able to modify its settings)"))
+ );
+ DECLARE_LOGGER;
+};
+REGISTER_SERIALIZABLE(Ig2_Polyhedra_Polyhedra_PolyhedraGeomOrScGeom);
+
//***************************************************************************
#endif // YADE_CGAL
