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.
assign_flags_to_model_part_entities_process.hpp
Go to the documentation of this file.
1 //
2 // Project Name: KratosSolidMechanicsApplication $
3 // Created by: $Author: JMCarbonell $
4 // Last modified by: $Co-Author: $
5 // Date: $Date: August 2016 $
6 // Revision: $Revision: 0.0 $
7 //
8 //
9 
10 #if !defined(KRATOS_ASSIGN_FLAGS_TO_MODEL_PART_ENTITIES_PROCESS_H_INCLUDED)
11 #define KRATOS_ASSIGN_FLAGS_TO_MODEL_PART_ENTITIES_PROCESS_H_INCLUDED
12 
13 
14 // System includes
15 
16 // External includes
17 
18 // Project includes
19 #include "includes/model_part.h"
21 #include "processes/process.h"
22 
23 namespace Kratos
24 {
25 
28 
30 
33 {
34 public:
37 
40 
44 
45 
46  AssignFlagsToModelPartEntitiesProcess(ModelPart& rModelPart, const std::string EntityType,
47  const std::vector<Flags>& rAssignFlags) : Process(Flags()) , mrModelPart(rModelPart), mEntityType(EntityType), mrTransferFlags(std::vector<Flags>()), mrAssignFlags(rAssignFlags)
48  {
50 
51  KRATOS_CATCH("");
52  }
53 
54 
55  AssignFlagsToModelPartEntitiesProcess(ModelPart& rModelPart, const std::string EntityType, const std::vector<Flags>& rAssignFlags,
56  const std::vector<Flags>& rTransferFlags) : Process(Flags()) , mrModelPart(rModelPart), mEntityType(EntityType), mrTransferFlags(rTransferFlags), mrAssignFlags(rAssignFlags)
57  {
59 
60  KRATOS_CATCH("");
61  }
62 
63 
66 
67 
71 
73  void operator()()
74  {
75  Execute();
76  }
77 
78 
82 
83 
85  void Execute() override
86  {
87  KRATOS_TRY;
88 
89  if (mEntityType == "Nodes")
90  {
91  const int nnodes = mrModelPart.Nodes().size();
92 
93  if (nnodes != 0)
94  {
95  ModelPart::NodesContainerType::iterator it_begin =
96  mrModelPart.NodesBegin();
97 
98  #pragma omp parallel for
99  for (int i = 0; i < nnodes; i++)
100  {
101  ModelPart::NodesContainerType::iterator it = it_begin + i;
102 
103  if (this->MatchTransferFlags(*(it.base())))
104  {
105  this->AssignFlags(*(it.base()));
106  }
107  }
108  }
109  }
110  else if (mEntityType == "Elements")
111  {
112  const int nelements = mrModelPart.Elements().size();
113 
114  if (nelements != 0)
115  {
116  ModelPart::ElementsContainerType::iterator it_begin =
117  mrModelPart.ElementsBegin();
118 
119  #pragma omp parallel for
120  for (int i = 0; i < nelements; i++)
121  {
122  ModelPart::ElementsContainerType::iterator it = it_begin + i;
123 
124  if (this->MatchTransferFlags(*(it.base())))
125  {
126  this->AssignFlags(*(it.base()));
127  }
128  }
129  }
130  }
131  else if (mEntityType == "Conditions")
132  {
133  const int nconditions = mrModelPart.Conditions().size();
134 
135  if (nconditions != 0)
136  {
137  ModelPart::ConditionsContainerType::iterator it_begin =
138  mrModelPart.ConditionsBegin();
139 
140  //#pragma omp parallel for
141  for (int i = 0; i < nconditions; i++)
142  {
143  ModelPart::ConditionsContainerType::iterator it = it_begin + i;
144 
145  if (this->MatchTransferFlags(*(it.base())))
146  {
147  this->AssignFlags(*(it.base()));
148  }
149  }
150  }
151  }
152 
153  KRATOS_CATCH("");
154  }
155 
158  void ExecuteInitialize() override
159  {
160  }
161 
165  {
166  }
167 
168 
171  {
172  }
173 
176  {
177  }
178 
179 
181  void ExecuteBeforeOutputStep() override
182  {
183  }
184 
185 
187  void ExecuteAfterOutputStep() override
188  {
189  }
190 
191 
194  void ExecuteFinalize() override
195  {
196  }
197 
198 
202 
203 
207 
208 
212 
214  std::string Info() const override
215  {
216  return "AssignFlagsToModelPartEntitiesProcess";
217  }
218 
220  void PrintInfo(std::ostream& rOStream) const override
221  {
222  rOStream << "AssignFlagsToModelPartEntitiesProcess";
223  }
224 
226  void PrintData(std::ostream& rOStream) const override
227  {
228  }
229 
230 
235 
236 protected:
237 
246 
249 
263 
264 private:
265 
271 
272  ModelPart& mrModelPart;
273 
274  const std::string mEntityType;
275 
276  const std::vector<Flags> mrTransferFlags;
277  const std::vector<Flags> mrAssignFlags;
278 
282 
283  bool MatchTransferFlags(const Node::Pointer& pNode)
284  {
285 
286  for(unsigned int i = 0; i<mrTransferFlags.size(); i++)
287  {
288  if( pNode->IsNot(mrTransferFlags[i]) )
289  return false;
290  }
291 
292  return true;
293 
294  }
295 
296  void AssignFlags(const Node::Pointer& pNode)
297  {
298 
299  for(unsigned int i = 0; i<mrAssignFlags.size(); i++)
300  pNode->Set(mrAssignFlags[i]);
301 
302  }
303 
304 
305  bool MatchTransferFlags(const Element::Pointer& pElement)
306  {
307 
308  for(unsigned int i = 0; i<mrTransferFlags.size(); i++)
309  {
310  if( pElement->IsNot(mrTransferFlags[i]) )
311  return false;
312  }
313 
314  return true;
315 
316  }
317 
318  void AssignFlags(const Element::Pointer& pElement)
319  {
320 
321  for(unsigned int i = 0; i<mrAssignFlags.size(); i++)
322  pElement->Set(mrAssignFlags[i]);
323 
324  }
325 
326 
327  bool MatchTransferFlags(const Condition::Pointer& pCondition)
328  {
329 
330  for(unsigned int i = 0; i<mrTransferFlags.size(); i++)
331  {
332  if( pCondition->IsNot(mrTransferFlags[i]) )
333  return false;
334  }
335 
336  return true;
337 
338  }
339 
340  void AssignFlags(const Condition::Pointer& pCondition)
341  {
342 
343  for(unsigned int i = 0; i<mrAssignFlags.size(); i++)
344  pCondition->Set(mrAssignFlags[i]);
345 
346  }
347 
354 
357 
358 
368 
369 }; // Class AssignFlagsToModelPartEntitiesProcess
370 
371 
373 
376 
377 
381 
382 
384 inline std::istream& operator >> (std::istream& rIStream,
386 
388 inline std::ostream& operator << (std::ostream& rOStream,
390 {
391  rThis.PrintInfo(rOStream);
392  rOStream << std::endl;
393  rThis.PrintData(rOStream);
394 
395  return rOStream;
396 }
398 
399 
400 } // namespace Kratos.
401 
402 #endif // KRATOS_ASSIGN_FLAGS_TO_MODEL_PART_ENTITIES_PROCESS_H_INCLUDED defined
The base class for assigning a value to scalar variables or array_1d components processes in Kratos.
Definition: assign_flags_to_model_part_entities_process.hpp:33
void ExecuteAfterOutputStep() override
this function will be executed at every time step AFTER writing the output
Definition: assign_flags_to_model_part_entities_process.hpp:187
void ExecuteInitializeSolutionStep() override
this function will be executed at every time step BEFORE performing the solve phase
Definition: assign_flags_to_model_part_entities_process.hpp:170
AssignFlagsToModelPartEntitiesProcess(ModelPart &rModelPart, const std::string EntityType, const std::vector< Flags > &rAssignFlags, const std::vector< Flags > &rTransferFlags)
Definition: assign_flags_to_model_part_entities_process.hpp:55
void ExecuteFinalizeSolutionStep() override
this function will be executed at every time step AFTER performing the solve phase
Definition: assign_flags_to_model_part_entities_process.hpp:175
void ExecuteBeforeSolutionLoop() override
Definition: assign_flags_to_model_part_entities_process.hpp:164
void ExecuteFinalize() override
Definition: assign_flags_to_model_part_entities_process.hpp:194
AssignFlagsToModelPartEntitiesProcess(AssignFlagsToModelPartEntitiesProcess const &rOther)
Copy constructor.
void ExecuteBeforeOutputStep() override
this function will be executed at every time step BEFORE writing the output
Definition: assign_flags_to_model_part_entities_process.hpp:181
void Execute() override
Execute method is used to execute the AssignFlagsToModelPartEntitiesProcess algorithms.
Definition: assign_flags_to_model_part_entities_process.hpp:85
~AssignFlagsToModelPartEntitiesProcess() override
Destructor.
Definition: assign_flags_to_model_part_entities_process.hpp:65
std::string Info() const override
Turn back information as a string.
Definition: assign_flags_to_model_part_entities_process.hpp:214
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: assign_flags_to_model_part_entities_process.hpp:226
void operator()()
This operator is provided to call the process as a function and simply calls the Execute method.
Definition: assign_flags_to_model_part_entities_process.hpp:73
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: assign_flags_to_model_part_entities_process.hpp:220
AssignFlagsToModelPartEntitiesProcess(ModelPart &rModelPart, const std::string EntityType, const std::vector< Flags > &rAssignFlags)
Definition: assign_flags_to_model_part_entities_process.hpp:46
void ExecuteInitialize() override
Definition: assign_flags_to_model_part_entities_process.hpp:158
KRATOS_CLASS_POINTER_DEFINITION(AssignFlagsToModelPartEntitiesProcess)
Pointer definition of AssignFlagsToModelPartEntitiesProcess.
Definition: flags.h:58
void Set(const Flags ThisFlag)
Definition: flags.cpp:33
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
ElementIterator ElementsBegin(IndexType ThisIndex=0)
Definition: model_part.h:1169
ConditionIterator ConditionsBegin(IndexType ThisIndex=0)
Definition: model_part.h:1361
NodeIterator NodesBegin(IndexType ThisIndex=0)
Definition: model_part.h:487
ConditionsContainerType & Conditions(IndexType ThisIndex=0)
Definition: model_part.h:1381
ElementsContainerType & Elements(IndexType ThisIndex=0)
Definition: model_part.h:1189
NodesContainerType & Nodes(IndexType ThisIndex=0)
Definition: model_part.h:507
The base class for all processes in Kratos.
Definition: process.h:49
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
int nnodes
Definition: sensitivityMatrix.py:24
namespace
Definition: array_1d.h:793
integer i
Definition: TensorModule.f:17