9
#include<yade/core/Shape.hpp>
10
#include<yade/core/IGeom.hpp>
11
#include<yade/core/GlobalEngine.hpp>
12
#include<yade/core/Material.hpp>
13
#include<yade/pkg/common/Aabb.hpp>
14
#include<yade/pkg/common/Dispatching.hpp>
15
#include<yade/pkg/dem/FrictPhys.hpp>
16
#include<yade/pkg/common/Wall.hpp>
17
#include<yade/pkg/common/Facet.hpp>
18
#include<yade/lib/base/openmp-accu.hpp>
9
#include<core/Shape.hpp>
10
#include<core/IGeom.hpp>
11
#include<core/GlobalEngine.hpp>
12
#include<core/Material.hpp>
13
#include<pkg/common/Aabb.hpp>
14
#include<pkg/common/Dispatching.hpp>
15
#include<pkg/dem/FrictPhys.hpp>
16
#include<pkg/common/Wall.hpp>
17
#include<pkg/common/Facet.hpp>
18
#include<lib/base/openmp-accu.hpp>
20
20
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
21
21
#include <CGAL/Delaunay_triangulation_3.h>
54
54
Vector3r GetCentroid(){Initialize(); return centroid;}
55
55
Vector3r GetInertia(){Initialize(); return inertia;}
56
56
vector<int> GetSurfaceTriangulation(){Initialize(); return faceTri;}
57
vector<vector<int>> GetSurfaces() const;
58
59
bool IsInitialized(){return init;}
59
60
std::vector<Vector3r> GetOriginalVertices();
60
double GetVolume(){Initialize(); return volume;}
61
Real GetVolume(){Initialize(); return volume;}
61
62
Quaternionr GetOri(){Initialize(); return orientation;}
62
63
Polyhedron GetPolyhedron(){return P;};
63
64
void Clear(){v.clear(); P.clear(); init = 0; size = Vector3r(1.,1.,1.); faceTri.clear();};
93
94
.def("GetOri",&Polyhedra::GetOri,"return polyhedra's orientation")
94
95
.def("GetCentroid",&Polyhedra::GetCentroid,"return polyhedra's centroid")
95
96
.def("GetSurfaceTriangulation",&Polyhedra::GetSurfaceTriangulation,"triangulation of facets (for plotting)")
97
.def("GetSurfaces",&Polyhedra::GetSurfaces,"get indices of surfaces' vertices (for postprocessing)")
97
99
REGISTER_CLASS_INDEX(Polyhedra,Shape);
147
149
/*! Elastic material */
148
150
class PolyhedraMat: public Material{
150
PolyhedraMat(double N, double S, double F){Kn=N; Ks=S; frictionAngle=F;};
151
double GetStrength(){return strength;};
152
PolyhedraMat(Real N, Real S, Real F){Kn=N; Ks=S; frictionAngle=F;};
153
Real GetStrength(){return strength;};
152
154
virtual ~PolyhedraMat(){};
153
155
YADE_CLASS_BASE_DOC_ATTRS_CTOR(PolyhedraMat,Material,"Elastic material with Coulomb friction.",
154
((Real,Kn,1e8,,"Normal volumetric 'stiffness' (N/m3)."))
156
((Real,Kn,1e8,,"Normal 'stiffness' (N/m3 for Law_.._Volumetric, N/m for Law2_.._Simple)."))
155
157
((Real,Ks,1e5,,"Shear stiffness (N/m)."))
156
158
((Real,frictionAngle,.5,,"Contact friction angle (in radians)."))
157
159
((bool,IsSplitable,0,,"To be splitted ... or not"))
158
((double,strength,100,,"Stress at whis polyhedra of volume 4/3*pi [mm] breaks.")),
160
((Real,strength,100,,"Stress at whis polyhedra of volume 4/3*pi [mm] breaks.")),
159
161
/*ctor*/ createIndex();
161
163
REGISTER_CLASS_INDEX(PolyhedraMat,Material);
181
183
//***************************************************************************
182
184
#ifdef YADE_OPENGL
183
#include<yade/pkg/common/GLDrawFunctors.hpp>
184
#include<yade/lib/opengl/OpenGLWrapper.hpp>
185
#include<yade/lib/opengl/GLUtils.hpp>
185
#include<pkg/common/GLDrawFunctors.hpp>
186
#include<lib/opengl/OpenGLWrapper.hpp>
187
#include<lib/opengl/GLUtils.hpp>
186
188
#include<GL/glu.h>
187
#include<yade/pkg/dem/Shop.hpp>
189
#include<pkg/dem/Shop.hpp>
189
191
/*! Draw Polyhedra using OpenGL */
190
192
class Gl1_Polyhedra: public GlShapeFunctor{
238
240
REGISTER_SERIALIZABLE(Ip2_PolyhedraMat_PolyhedraMat_PolyhedraPhys);
240
242
//***************************************************************************
241
/*! Calculate physical response based on penetration configuration given by TTetraGeom. */
243
/*! Calculate physical response based on penetration configuration given by PolyhedraGeom. */
243
class PolyhedraVolumetricLaw: public LawFunctor{
245
class Law2_PolyhedraGeom_PolyhedraPhys_Volumetric: public LawFunctor{
244
246
OpenMPAccumulator<Real> plasticDissipation;
245
247
virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
246
248
Real elasticEnergy ();
247
249
Real getPlasticDissipation();
248
250
void initPlasticDissipation(Real initVal=0);
249
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(PolyhedraVolumetricLaw,LawFunctor,"Calculate physical response of 2 :yref:`vector<Polyhedra>` in interaction, based on penetration configuration given by :yref:`PolyhedraGeom`.",
251
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_PolyhedraGeom_PolyhedraPhys_Volumetric,LawFunctor,"Calculate physical response of 2 :yref:`vector<Polyhedra>` in interaction, based on penetration configuration given by :yref:`PolyhedraGeom`. Normal force is proportional to the volume of intersection",
252
((Real,volumePower,1.,,"Power of volume used in evaluation of normal force. Default is 1.0 - normal force is linearly proportional to volume. 1.0/3.0 would mean that normal force is proportional to the cube root of volume, approximation of penetration depth."))
250
253
((Vector3r,shearForce,Vector3r::Zero(),,"Shear force from last step"))
251
254
((bool,traceEnergy,false,,"Define the total energy dissipated in plastic slips at all contacts. This will trace only plastic energy in this law, see O.trackEnergy for a more complete energies tracing"))
252
255
((int,plastDissipIx,-1,(Attr::hidden|Attr::noSave),"Index for plastic dissipation (with O.trackEnergy)"))
253
256
((int,elastPotentialIx,-1,(Attr::hidden|Attr::noSave),"Index for elastic potential energy (with O.trackEnergy)"))
255
.def("elasticEnergy",&PolyhedraVolumetricLaw::elasticEnergy,"Compute and return the total elastic energy in all \"FrictPhys\" contacts")
256
.def("plasticDissipation",&PolyhedraVolumetricLaw::getPlasticDissipation,"Total energy dissipated in plastic slips at all FrictPhys contacts. Computed only if :yref:`Law2_ScGeom_FrictPhys_CundallStrack::traceEnergy` is true.")
257
.def("initPlasticDissipation",&PolyhedraVolumetricLaw::initPlasticDissipation,"Initialize cummulated plastic dissipation to a value (0 by default).")
258
.def("elasticEnergy",&Law2_PolyhedraGeom_PolyhedraPhys_Volumetric::elasticEnergy,"Compute and return the total elastic energy in all \"FrictPhys\" contacts")
259
.def("plasticDissipation",&Law2_PolyhedraGeom_PolyhedraPhys_Volumetric::getPlasticDissipation,"Total energy dissipated in plastic slips at all FrictPhys contacts. Computed only if :yref:`Law2_PolyhedraGeom_PolyhedraPhys_Volumetric::traceEnergy` is true.")
260
.def("initPlasticDissipation",&Law2_PolyhedraGeom_PolyhedraPhys_Volumetric::initPlasticDissipation,"Initialize cummulated plastic dissipation to a value (0 by default).")
259
262
FUNCTOR2D(PolyhedraGeom,PolyhedraPhys);
262
REGISTER_SERIALIZABLE(PolyhedraVolumetricLaw);
265
REGISTER_SERIALIZABLE(Law2_PolyhedraGeom_PolyhedraPhys_Volumetric);
265
268
//***************************************************************************
283
286
//Test if point is inside Polyhedron
284
287
bool Is_inside_Polyhedron(Polyhedron P, CGALpoint inside);
285
288
//return approximate intersection of sphere & polyhedron
286
bool Sphere_Polyhedron_intersection(Polyhedron A, double r, CGALpoint C, CGALpoint centroid, double volume, CGALvector normal, double area);
289
bool Sphere_Polyhedron_intersection(Polyhedron A, Real r, CGALpoint C, CGALpoint centroid, Real volume, CGALvector normal, Real area);
287
290
//return volume and centroid of polyhedra
288
bool P_volume_centroid(Polyhedron P, double * volume, Vector3r * centroid);
291
bool P_volume_centroid(Polyhedron P, Real * volume, Vector3r * centroid);
289
292
//CGAL - miniEigen communication
290
293
Vector3r FromCGALPoint(CGALpoint A);
291
294
Vector3r FromCGALVector(CGALvector A);
295
298
bool do_intersect(Polyhedron A, Polyhedron B);
296
299
bool do_intersect(Polyhedron A, Polyhedron B, std::vector<int> &sep_plane);
297
300
//connect triagular facets if possible
298
Polyhedron Simplify(Polyhedron P, double lim);
301
Polyhedron Simplify(Polyhedron P, Real lim);
299
302
//list of facets and edges
300
303
void PrintPolyhedron(Polyhedron P);
301
304
void PrintPolyhedron2File(Polyhedron P,FILE* X);
302
305
//normal by least square fitting of separating segments
303
306
Vector3r FindNormal(Polyhedron Int, Polyhedron PA, Polyhedron PB);
304
307
//calculate area of projection of polyhedron into the plane
305
double CalculateProjectionArea(Polyhedron Int, CGALvector CGALnormal);
308
Real CalculateProjectionArea(Polyhedron Int, CGALvector CGALnormal);
306
309
//split polyhedron
307
310
shared_ptr<Body> SplitPolyhedra(const shared_ptr<Body>& body, Vector3r direction, Vector3r point);
added
removed
pkg/dem/Polyhedra.hpp
