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.
mark_close_nodes_process.h
Go to the documentation of this file.
1 /*
2 ==============================================================================
3 KratosULFApplication
4 A library based on:
5 Kratos
6 A General Purpose Software for Multi-Physics Finite Element Analysis
7 Version 1.0 (Released on march 05, 2007).
8 
9 Copyright 2007
10 Pooyan Dadvand, Riccardo Rossi, Pawel Ryzhakov
11 pooyan@cimne.upc.edu
12 rrossi@cimne.upc.edu
13 - CIMNE (International Center for Numerical Methods in Engineering),
14 Gran Capita' s/n, 08034 Barcelona, Spain
15 
16 
17 Permission is hereby granted, free of charge, to any person obtaining
18 a copy of this software and associated documentation files (the
19 "Software"), to deal in the Software without restriction, including
20 without limitation the rights to use, copy, modify, merge, publish,
21 distribute, sublicense and/or sell copies of the Software, and to
22 permit persons to whom the Software is furnished to do so, subject to
23 the following condition:
24 
25 Distribution of this code for any commercial purpose is permissible
26 ONLY BY DIRECT ARRANGEMENT WITH THE COPYRIGHT OWNERS.
27 
28 The above copyright notice and this permission notice shall be
29 included in all copies or substantial portions of the Software.
30 
31 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
32 EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
33 MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
34 IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
35 CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
36 TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
37 SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
38 
39 ==============================================================================
40 */
41 
42 
43 //
44 // Project Name: Kratos
45 // Last Modified by: $Author: rrossi $
46 // Date: $Date: 2007-03-06 10:30:32 $
47 // Revision: $Revision: 1.2 $
48 //
49 //
50 
51 #if !defined(KRATOS_MARK_CLOSE_NODES_PROCESS_INCLUDED )
52 #define KRATOS_MARK_CLOSE_NODES_PROCESS_INCLUDED
53 
54 
55 
56 // System includes
57 #include <string>
58 #include <iostream>
59 #include <algorithm>
60 
61 // External includes
62 
63 
64 // Project includes
65 #include "includes/define.h"
66 #include "processes/process.h"
67 #include "includes/node.h"
68 #include "includes/element.h"
69 #include "includes/model_part.h"
70 #include "custom_utilities/geometry_utilities2D.h"
72 
73 
74 namespace Kratos
75 {
76 
79 
83 
84 
88 
92 
96 
98 
104 class MarkCloseNodesProcess
105  : public Process
106 {
107 public:
110 
113 
117 
120  : mr_model_part(model_part)
121  {
122  }
123 
126  {
127  }
128 
129 
133 
134  void operator()()
135  {
136  Execute();
137  }
138 
139 
143 
144  void Execute() override
145  {
146  KRATOS_TRY
147  for(ModelPart::NodesContainerType::iterator in = mr_model_part.NodesBegin() ;
148  in != mr_model_part.NodesEnd() ; ++in)
149  {
150 
151  const double& X0 = in->X();
152  const double& Y0 = in->Y();
153  //KRATOS_WATCH("ENTERED MARKING CLOSE NODES FUCNTION!");
154 
155  for( GlobalPointersVector< Node >::iterator i = in->GetValue(NEIGHBOUR_NODES).begin();
156  i != in->GetValue(NEIGHBOUR_NODES).end(); i++)
157  {
158  const double& X1 = i->X();
159  const double& Y1 = i->Y();
160  const double& dist_sq = (X1-X0)*(X1-X0)+(Y1-Y0)*(Y1-Y0);
161  if (dist_sq<0.1*(i->GetValue(NODAL_H))*(i->GetValue(NODAL_H)) && in->GetId()>i->GetId())
162  //if (dist_sq<0.001 && in->GetId()>i->GetId())
163  {
164  i->Is(TO_ERASE)= true;
165  KRATOS_WATCH("ERASING NODE!!!!!!");
166  KRATOS_WATCH(in->GetId());
167 
168  }
169 
170  }
171 
172 
173  }
174 
175 
176  KRATOS_CATCH("")
177  }
178 
179 
183 
184 
188 
189 
193 
195  std::string Info() const override
196  {
197  return "MarkCloseNodesProcess";
198  }
199 
201  void PrintInfo(std::ostream& rOStream) const override
202  {
203  rOStream << "MarkCloseNodesProcess";
204  }
205 
207  void PrintData(std::ostream& rOStream) const override
208  {
209  }
210 
211 
215 
216 
218 
219 protected:
222 
223 
227 
228 
232 
233 
237 
238 
242 
243 
247 
248 
252 
253 
255 
256 private:
259 
260 
264  ModelPart& mr_model_part;
265 
269 
270 
274 
275 
279 
280 
284 
285 
289 
291 // MarkCloseNodesProcess& operator=(MarkCloseNodesProcess const& rOther);
292 
294 // MarkCloseNodesProcess(MarkCloseNodesProcess const& rOther);
295 
296 
298 
299 }; // Class MarkCloseNodesProcess
300 
302 
305 
306 
310 
311 
313 inline std::istream& operator >> (std::istream& rIStream,
314  MarkCloseNodesProcess& rThis);
315 
317 inline std::ostream& operator << (std::ostream& rOStream,
318  const MarkCloseNodesProcess& rThis)
319 {
320  rThis.PrintInfo(rOStream);
321  rOStream << std::endl;
322  rThis.PrintData(rOStream);
323 
324  return rOStream;
325 }
327 
328 
329 } // namespace Kratos.
330 
331 #endif // KRATOS_MARK_CLOSE_NODES_PROCESS_INCLUDED defined
332 
333 
boost::indirect_iterator< typename TContainerType::iterator > iterator
Definition: global_pointers_vector.h:79
Short class definition.
Definition: mark_close_nodes_process.h:109
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: mark_close_nodes_process.h:207
~MarkCloseNodesProcess() override
Destructor.
Definition: mark_close_nodes_process.h:125
void operator()()
Definition: mark_close_nodes_process.h:134
KRATOS_CLASS_POINTER_DEFINITION(MarkCloseNodesProcess)
Pointer definition of PushStructureProcess.
MarkCloseNodesProcess(ModelPart &model_part)
Default constructor.
Definition: mark_close_nodes_process.h:119
std::string Info() const override
Turn back information as a string.
Definition: mark_close_nodes_process.h:195
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: mark_close_nodes_process.h:201
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: mark_close_nodes_process.h:144
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
NodeIterator NodesBegin(IndexType ThisIndex=0)
Definition: model_part.h:487
NodeIterator NodesEnd(IndexType ThisIndex=0)
Definition: model_part.h:497
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_WATCH(variable)
Definition: define.h:806
#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
model_part
Definition: face_heat.py:14
X1
Definition: generate_frictional_mortar_condition.py:119
integer i
Definition: TensorModule.f:17