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.
set_h_map_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: anonymous $
46 // Date: $Date: 2009-01-15 14:50:34 $
47 // Revision: $Revision: 1.3 $
48 //
49 //
50 
51 
52 #if !defined(KRATOS_SET_H_MAP_PROCESS_INCLUDED )
53 #define KRATOS_SET_H_MAP_PROCESS_INCLUDED
54 
55 
56 
57 // System includes
58 #include <string>
59 #include <iostream>
60 #include <algorithm>
61 
62 // External includes
63 
64 
65 // Project includes
66 #include "includes/define.h"
67 #include "processes/process.h"
68 #include "includes/node.h"
69 #include "includes/element.h"
70 #include "includes/model_part.h"
71 
72 
73 namespace Kratos
74 {
75 
78 
82 
83 
87 
91 
95 
97 
105  : public Process
106 {
107 public:
110 
113 
117 
120  : mr_model_part(model_part)
121  {
122  }
123 
125  virtual ~SetHMapProcess()
126  {
127  }
128 
129 
133 
134  void operator()()
135  {
136  Execute();
137  }
138 
139 
143 
144  void CalculateOptimalH(const double h_min, const double h_max, const double mmax_dist)
145  {
146  KRATOS_TRY
147  KRATOS_WATCH("Calculate Optimal H Process is executed")
148  double max_dist=0.0;
149  double dist=0.0;
150 
151  ModelPart::NodesContainerType& rNodes = mr_model_part.Nodes();
152 
153  //we set the minimal distance, where the refinement shall start
154  double min_dist=h_min;
155  //min_dist=0.0;
156  for(ModelPart::NodesContainerType::iterator in = rNodes.begin(); in!=rNodes.end(); in++)
157  {
158 
159  dist=in->FastGetSolutionStepValue(DISTANCE);
160  if (max_dist<dist)
161  {
162  max_dist=dist;
163  }
164  }
165  KRATOS_WATCH(max_dist)
166  //if the max_dist is not specified - we shall take the largest distance that exists in the domain to be the max one
167  //if it is specified (non-zero value) - > take it
168  if (mmax_dist!=0.0)
169  max_dist=mmax_dist;
170 
171  //coeffcients that define linear distribution from minimal to maximal mesh size
172  double coef=(h_max-h_min)/(0.75*max_dist-h_min);
173  double c=h_min-coef*h_min;
174 
175  for(ModelPart::NodesContainerType::iterator in = rNodes.begin(); in!=rNodes.end(); in++)
176  {
177 
178  dist=in->FastGetSolutionStepValue(DISTANCE);
179  if (dist<0.75*max_dist && dist>=min_dist)
180  {
181  if ((coef*dist+c)>h_min && (coef*dist+c)<h_max)
182  {
183  in->FastGetSolutionStepValue(NODAL_H)=coef*dist+c;
184  }
185  }
186 
187  }
188  KRATOS_CATCH("")
189  }
190 
191 
195 
196 
200 
201 
205 
207  virtual std::string Info() const
208  {
209  return "SetHMapProcess";
210  }
211 
213  virtual void PrintInfo(std::ostream& rOStream) const
214  {
215  rOStream << "SetHMapProcess";
216  }
217 
219  virtual void PrintData(std::ostream& rOStream) const
220  {
221  }
222 
223 
227 
228 
230 
231 protected:
234 
235 
239 
240 
244 
245 
249 
250 
254 
255 
259 
260 
264 
265 
267 
268 private:
271 
272 
276  ModelPart& mr_model_part;
277 
278 
279 
283 
287 
288 
292 
293 
297 
298 
302 
304  SetHMapProcess& operator=(SetHMapProcess const& rOther) {};
305 
307  //SetHMapProcess(SetHMapProcess const& rOther);
308 
309 
311 
312 }; // Class SetHMapProcess
313 
315 
318 
319 
323 
324 
326 inline std::istream& operator >> (std::istream& rIStream,
327  SetHMapProcess& rThis);
328 
330 inline std::ostream& operator << (std::ostream& rOStream,
331  const SetHMapProcess& rThis)
332 {
333  rThis.PrintInfo(rOStream);
334  rOStream << std::endl;
335  rThis.PrintData(rOStream);
336 
337  return rOStream;
338 }
340 
341 
342 } // namespace Kratos.
343 
344 #endif // KRATOS_SET_H_MAP_PROCESS_INCLUDED defined
345 
346 
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
MeshType::NodesContainerType NodesContainerType
Nodes container. Which is a vector set of nodes with their Id's as key.
Definition: model_part.h:128
NodesContainerType & Nodes(IndexType ThisIndex=0)
Definition: model_part.h:507
The base class for all processes in Kratos.
Definition: process.h:49
virtual void Execute()
Execute method is used to execute the Process algorithms.
Definition: process.h:101
Short class definition.
Definition: set_h_map_process.h:106
virtual void PrintInfo(std::ostream &rOStream) const
Print information about this object.
Definition: set_h_map_process.h:213
virtual ~SetHMapProcess()
Destructor.
Definition: set_h_map_process.h:125
SetHMapProcess(ModelPart &model_part)
Default constructor.
Definition: set_h_map_process.h:119
KRATOS_CLASS_POINTER_DEFINITION(SetHMapProcess)
Pointer definition of SetHMapProcess.
virtual std::string Info() const
Turn back information as a string.
Definition: set_h_map_process.h:207
void operator()()
Definition: set_h_map_process.h:134
virtual void PrintData(std::ostream &rOStream) const
Print object's data.
Definition: set_h_map_process.h:219
void CalculateOptimalH(const double h_min, const double h_max, const double mmax_dist)
Definition: set_h_map_process.h:144
#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
float dist
Definition: edgebased_PureConvection.py:89
model_part
Definition: face_heat.py:14
c
Definition: generate_weakly_compressible_navier_stokes_element.py:108