KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
create_and_destroy.h
Go to the documentation of this file.
1 //
2 // Authors:
3 // Miguel Angel Celigueta maceli@cimne.upc.edu
4 //
5 //README::::look to the key word "VERSION" if you want to find all the points where you have to change something so that you can pass from a kdtree to a bin data search structure;
6 
7 #ifndef CREATE_AND_DESTROY_H
8 #define CREATE_AND_DESTROY_H
9 
10 
11 // System includes
12 #include <string>
13 #include <iostream>
14 
15 // Project includes
16 #include "includes/define.h"
18 #include "includes/model_part.h"
19 #include "includes/kratos_flags.h"
20 #include "utilities/timer.h"
21 #include "utilities/openmp_utils.h"
22 #include "utilities/quaternion.h"
29 
30 
31 namespace Kratos {
32 
33 class KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor {
35 
36 public:
37 
38  static const std::size_t space_dim = 3;
46  unsigned int mMaxNodeId;
47 
49 
52 
54 
55  ParticleCreatorDestructor(AnalyticWatcher::Pointer p_watcher);
56 
57  ParticleCreatorDestructor(AnalyticWatcher::Pointer p_watcher, Parameters settings);
58 
60  virtual ~ParticleCreatorDestructor();
61 
62  int FindMaxNodeIdInModelPart(ModelPart& r_modelpart);
63  void FindAndSaveMaxNodeIdInModelPart(ModelPart& r_modelpart);
64  int FindMaxElementIdInModelPart(ModelPart& r_modelpart);
65  int FindMaxConditionIdInModelPart(ModelPart& r_modelpart);
66  void RenumberElementIdsFromGivenValue(ModelPart& r_modelpart, const int initial_id);
67  void DestroyMarkedParticles(ModelPart& r_model_part);
68  virtual double SelectRadius(bool initial,
69  ModelPart& r_sub_model_part_with_parameters,
70  std::map<std::string, std::unique_ptr<RandomVariable>>& r_random_variables_map);
71 
72  void NodeCreatorWithPhysicalParameters(ModelPart& r_modelpart,
73  Node ::Pointer& pnew_node,
74  int aId,
75  Node ::Pointer & reference_node,
76  double radius,
77  Properties& params,
78  ModelPart& r_sub_model_part_with_parameters,
79  bool has_sphericity,
80  bool has_rotation,
81  bool initial);
82 
83  void NodeForClustersCreatorWithPhysicalParameters(ModelPart& r_modelpart,
84  Node ::Pointer& pnew_node,
85  int aId,
86  Node ::Pointer& reference_node,
87  Properties& params,
88  ModelPart& r_sub_model_part_with_parameters,
89  bool has_sphericity,
90  bool has_rotation,
91  bool initial);
92 
93  SphericParticle* ElementCreatorWithPhysicalParameters(ModelPart& r_modelpart,
94  int r_Elem_Id,
95  Node ::Pointer reference_node,
96  Element::Pointer injector_element,
97  Properties::Pointer r_params,
98  ModelPart& r_sub_model_part_with_parameters,
99  std::map<std::string, std::unique_ptr<RandomVariable>>& r_random_variables_map,
100  const Element& r_reference_element,
101  PropertiesProxy* p_fast_properties,
102  bool has_sphericity,
103  bool has_rotation,
104  bool initial,
105  ElementsContainerType& array_of_injector_elements);
106 
107  SphericParticle* AddInitialDataToNewlyCreatedElementAndNode(ModelPart& r_modelpart,
108  Properties::Pointer r_params,
109  const double radius, Node::Pointer& pnew_node,
110  Element::Pointer& p_particle);
111 
112 
113  SphericParticle* CreateSphericParticleRaw(ModelPart& r_modelpart,
114  int r_Elem_Id,
115  const array_1d<double, 3 >& coordinates,
116  Properties::Pointer r_params,
117  const double radius,
118  const Element& r_reference_element);
119 
120  SphericParticle* CreateSphericParticleRaw(ModelPart& r_modelpart,
121  int r_Elem_Id,
122  Node ::Pointer reference_node,
123  Properties::Pointer r_params,
124  const double radius,
125  const Element& r_reference_element);
126 
127  SphericParticle* CreateSphericParticleRaw(ModelPart& r_modelpart,
128  int r_Elem_Id,
129  Node ::Pointer reference_node,
130  Properties::Pointer r_params,
131  const double radius,
132  const std::string& element_type);
133 
134  SphericParticle* CreateSphericParticleRaw(ModelPart& r_modelpart,
135  Node ::Pointer reference_node,
136  Properties::Pointer r_params,
137  const double radius,
138  const std::string& element_type);
139 
140  SphericParticle* CreateSphericParticleRaw(ModelPart& r_modelpart,
141  int r_Elem_Id,
142  const array_1d<double, 3 >& coordinates,
143  Properties::Pointer r_params,
144  const double radius,
145  const std::string& element_type);
146 
147  SphericParticle* CreateSphericParticleRaw(ModelPart& r_modelpart,
148  const array_1d<double, 3 >& coordinates,
149  Properties::Pointer r_params,
150  const double radius,
151  const std::string& element_type);
152 
153  Element::Pointer CreateSphericParticle(ModelPart& r_modelpart,
154  int r_Elem_Id,
155  const array_1d<double, 3 >& coordinates,
156  Properties::Pointer r_params,
157  const double radius,
158  const Element& r_reference_element);
159 
160  Element::Pointer CreateSphericParticle(ModelPart& r_modelpart,
161  int r_Elem_Id,
162  Node ::Pointer reference_node,
163  Properties::Pointer r_params,
164  const double radius,
165  const Element& r_reference_element);
166 
167  Element::Pointer CreateSphericParticle(ModelPart& r_modelpart,
168  int r_Elem_Id,
169  Node ::Pointer reference_node,
170  Properties::Pointer r_params,
171  const double radius,
172  const std::string& element_type);
173 
174  Element::Pointer CreateSphericParticle(ModelPart& r_modelpart,
175  Node ::Pointer reference_node,
176  Properties::Pointer r_params,
177  const double radius,
178  const std::string& element_type);
179 
180  Element::Pointer CreateSphericParticle(ModelPart& r_modelpart,
181  int r_Elem_Id,
182  const array_1d<double, 3 >& coordinates,
183  Properties::Pointer r_params,
184  const double radius,
185  const std::string& element_type);
186 
187  Element::Pointer CreateSphericParticle(ModelPart& r_modelpart,
188  const array_1d<double, 3 >& coordinates,
189  Properties::Pointer r_params,
190  const double radius,
191  const std::string& element_type);
192 
193 
194  Cluster3D* ClusterCreatorWithPhysicalParameters(ModelPart& r_modelpart,
195  ModelPart& r_clusters_modelpart,
196  int r_Elem_Id,
197  Node ::Pointer reference_node,
198  Element::Pointer injector_element,
199  Properties::Pointer r_params,
200  ModelPart& r_sub_model_part_with_parameters,
201  const Element& r_reference_element,
202  PropertiesProxy* p_fast_properties,
203  bool has_sphericity,
204  bool has_rotation,
205  ElementsContainerType& array_of_injector_elements,
206  int& number_of_added_spheres,
207  const bool mStrategyForContinuum,
208  std::vector<SphericParticle*>& new_component_spheres);
209 
210 
211  void NodeCreatorForClusters(ModelPart& r_modelpart,
212  Node ::Pointer& pnew_node,
213  int aId,
214  array_1d<double, 3 >& reference_coordinates,
215  double radius,
216  Properties& params);
217 
218  void CentroidCreatorForRigidBodyElements(ModelPart& r_modelpart,
219  Node::Pointer& pnew_node,
220  int aId,
221  array_1d<double, 3>& reference_coordinates);
222 
223  SphericParticle* SphereCreatorForClusters(ModelPart& r_modelpart,
224  Node ::Pointer& pnew_node,
225  int r_Elem_Id,
226  double radius,
227  array_1d<double, 3 >& reference_coordinates,
228  double cluster_mass,
229  Properties::Pointer r_params,
230  const Element& r_reference_element,
231  const int cluster_id,
232  PropertiesProxy* p_fast_properties);
233 
234  SphericParticle* SphereCreatorForBreakableClusters(ModelPart& r_modelpart,
235  Node ::Pointer& pnew_node,
236  int r_Elem_Id,
237  double radius,
238  array_1d<double, 3>& reference_coordinates,
239  Properties::Pointer r_params,
240  const Element& r_reference_element,
241  const int cluster_id,
242  PropertiesProxy* p_fast_properties);
243 
244  void CalculateSurroundingBoundingBox(ModelPart& r_balls_model_part,
245  ModelPart& r_clusters_model_part,
246  ModelPart& r_rigid_faces_model_part,
247  ModelPart& r_dem_inlet_model_part,
248  double scale_factor,
249  bool automatic);
250 
251  void UpdateSurroundingBoundingBox(ModelPart& spheres_model_part);
252 
253  template<class TParticleType>
254  bool CheckParticlePreservationCriteria(const Element::Pointer p_element, const double current_time);
255 
256  template<class TParticleType>
257  void DestroyParticles(ModelPart& r_model_part);
258 
259  void DestroyParticleElements(ModelPart& r_model_part, Flags flag_for_destruction);
260 
261  template<class TParticleType>
262  void DestroyParticles(ModelPart::MeshType& rMesh, const double current_time);
263  void DestroyContactElements(ModelPart& r_model_part);
264  void MarkIsolatedParticlesForErasing(ModelPart& r_model_part);
266  void RemoveUnusedNodesOfTheClustersModelPart(ModelPart& r_clusters_modelpart);
267 
268  template<class TParticleType>
269  void MarkDistantParticlesForErasing(ModelPart& r_model_part);
270  void MarkParticlesForErasingGivenScalarVariableValue(ModelPart& r_model_part, const Variable<double>& rVariable, double value, double tol);
271  void MarkParticlesForErasingGivenVectorVariableModulus(ModelPart& r_model_part, const Variable<array_1d<double, 3> >& rVariable, double value, double tol);
272 
273  template<class TParticleType>
274  void MarkParticlesForErasingGivenBoundingBox(ModelPart& r_model_part, array_1d<double, 3> low_point, array_1d<double, 3> high_point);
275  void MarkParticlesForErasingGivenCylinder(ModelPart& r_model_part, array_1d<double, 3 > center, array_1d<double, 3 > axis_vector, const double radius);
276  void MarkContactElementsForErasing(ModelPart& r_model_part, ModelPart& mcontacts_model_part);
277 
278  template<class TParticleType>
279  void DestroyParticlesOutsideBoundingBox(ModelPart& r_model_part);
280  void MoveParticlesOutsideBoundingBoxBackInside(ModelPart& r_model_part);
281  void DestroyContactElementsOutsideBoundingBox(ModelPart& r_model_part, ModelPart& mcontacts_model_part);
282  Element::Pointer GetAnalyticReplacement(const Element& sample_element, Geometry<Node >::PointsArrayType nodelist, Element::Pointer p_elem_to_be_replaced, ModelPart& spheres_model_part);
283  static double rand_normal(const double mean, const double stddev, const double max_radius, const double min_radius);
284  static double rand_lognormal(const double mean, const double stddev, const double max_radius, const double min_radius);
285 
286  array_1d<double, 3> GetHighNode();
287  array_1d<double, 3> GetLowNode();
288  array_1d<double, 3> GetStrictHighNode();
289  array_1d<double, 3> GetStrictLowNode();
290 
291  double GetDiameter();
292  double GetStrictDiameter();
293  void SetHighNode(array_1d<double, 3> node);
294  void SetLowNode(array_1d<double, 3> node);
295  unsigned int GetCurrentMaxNodeId();
296  unsigned int* pGetCurrentMaxNodeId();
297  void SetMaxNodeId(unsigned int id);
298 
300  virtual std::string Info() const;
301 
303  virtual void PrintInfo(std::ostream& rOStream) const;
304 
306  virtual void PrintData(std::ostream& rOStream) const;
307 
308 protected:
309 
310 
311 private:
312 
313  void SetDoSearchNeighbourElements(bool true_or_false){mDoSearchNeighbourElements = true_or_false;}
314 
315  array_1d<double, 3 > mHighPoint;
316  array_1d<double, 3 > mLowPoint;
317  array_1d<double, 3 > mStrictHighPoint;
318  array_1d<double, 3 > mStrictLowPoint;
319  double mDiameter;
320  double mStrictDiameter;
321  double mScaleFactor=1.0;
322  bool mDoSearchNeighbourElements;
323  AnalyticWatcher::Pointer mpAnalyticWatcher;
324  Parameters mSettings;
325  void Clear(ModelPart::NodesContainerType::iterator node_it, int step_data_size);
326  inline void ClearVariables(ModelPart::NodesContainerType::iterator node_it, Variable<array_1d<double, 3 > >& rVariable);
327  inline void ClearVariables(ParticleIterator particle_it, Variable<double>& rVariable);
328 
331 
332 }; // Class ParticleCreatorDestructor
333 
334 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::DestroyParticles<SphericParticle>(ModelPart&);
335 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::DestroyParticles<Cluster3D>(ModelPart&);
336 extern template bool KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::CheckParticlePreservationCriteria<SphericParticle>(const Element::Pointer, const double);
337 extern template bool KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::CheckParticlePreservationCriteria<Cluster3D>(const Element::Pointer, const double);
338 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::DestroyParticles<SphericParticle>(ModelPart::MeshType&, const double);
339 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::DestroyParticles<Cluster3D>(ModelPart::MeshType&, const double);
340 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::MarkDistantParticlesForErasing<SphericParticle>(ModelPart&);
341 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::MarkDistantParticlesForErasing<Cluster3D>(ModelPart&);
342 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::MarkParticlesForErasingGivenBoundingBox<SphericParticle>(ModelPart&, array_1d<double, 3 >, array_1d<double, 3 >);
343 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::MarkParticlesForErasingGivenBoundingBox<Cluster3D>(ModelPart&, array_1d<double, 3 >, array_1d<double, 3 >);
344 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::DestroyParticlesOutsideBoundingBox<SphericParticle>(ModelPart&);
345 extern template void KRATOS_API(DEM_APPLICATION) ParticleCreatorDestructor::DestroyParticlesOutsideBoundingBox<Cluster3D>(ModelPart&);
346 
347 } // namespace Kratos
348 
349 #endif // CREATE_AND_DESTROY_H defined
350 
351 
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
std::string Info() const override
Turn back information as a string.
Definition: periodic_interface_process.hpp:93
Definition: cluster3D.h:29
Definition: discrete_particle_configure.h:46
ContainerType::iterator IteratorType
Definition: discrete_particle_configure.h:68
SearchType::ElementsContainerType::ContainerType ContainerType
Definition: discrete_particle_configure.h:63
Base class for all Elements.
Definition: element.h:60
Definition: explicit_solver_strategy.h:70
Definition: flags.h:58
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
MeshType::ElementsContainerType ElementsContainerType
Element container. A vector set of Elements with their Id's as key.
Definition: model_part.h:168
Mesh< NodeType, PropertiesType, ElementType, ConditionType > MeshType
Definition: model_part.h:123
MeshType::NodesContainerType NodesContainerType
Nodes container. Which is a vector set of nodes with their Id's as key.
Definition: model_part.h:128
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
Definition: create_and_destroy.h:33
Configure::IteratorType ParticleIterator
Definition: create_and_destroy.h:42
unsigned int mMaxNodeId
Definition: create_and_destroy.h:46
PointerVectorSet< Element, IndexedObject > ElementsContainerType
Definition: create_and_destroy.h:43
DiscreteParticleConfigure< space_dim > Configure
WARNING: generalize to 2d.
Definition: create_and_destroy.h:39
ModelPart::NodesContainerType NodesArrayType
Definition: create_and_destroy.h:45
ModelPart::ElementsContainerType ElementsArrayType
Definition: create_and_destroy.h:44
void MarkInitialNeighboursThatAreBeingRemoved(ModelPart &r_model_part)
Configure::ContainerType ParticlePointerVector
Definition: create_and_destroy.h:40
ParticlePointerVector::iterator ParticlePointerIterator
Definition: create_and_destroy.h:41
KRATOS_CLASS_POINTER_DEFINITION(ParticleCreatorDestructor)
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
boost::indirect_iterator< typename TContainerType::iterator > iterator
Definition: pointer_vector.h:89
A sorted associative container similar to an STL set, but uses a vector to store pointers to its data...
Definition: pointer_vector_set.h:72
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
Definition: properties_proxies.h:18
Definition: spheric_particle.h:31
#define KRATOS_API(...)
Definition: kratos_export_api.h:40
Kratos::ModelPart ModelPart
Definition: kratos_wrapper.h:31
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
int tol
Definition: hinsberg_optimization.py:138
element_type
Choosing element type for lagrangian_model_part.
Definition: lagrangian_droplet_test.py:56
float radius
Definition: mesh_to_mdpa_converter.py:18
subroutine initial(STRESS, T, DSTRAN, DEPS, NTENS, NDI, NSHR)
Definition: TensorModule.f:377
Definition: mesh_converter.cpp:38