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.
mass_calculate_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-05-16 13:59:01 $
47 // Revision: $Revision: 1.4 $ 12 November 2007 - 3D added
48 //
49 //
50 // THIS PROCESS is INVENTED in order to CALCULATE THE MASS and
51 //STORE it node-wise.
52 //this is necessary for the projection step of ulf-frac method
53 
54 #if !defined(KRATOS_MASS_CALCULATE_PROCESS_INCLUDED )
55 #define KRATOS_MASS_CALCULATE_PROCESS_INCLUDED
56 
57 
58 
59 // System includes
60 #include <string>
61 #include <iostream>
62 #include <algorithm>
63 
64 // External includes
65 
66 
67 // Project includes
68 #include "includes/define.h"
69 #include "processes/process.h"
70 #include "includes/node.h"
71 #include "includes/element.h"
72 #include "includes/model_part.h"
73 #include "utilities/geometry_utilities.h"
74 #include "ULF_application.h"
75 
76 
77 namespace Kratos
78 {
79 
82 
86 
87 
91 
95 
99 
101 
108  : public Process
109 {
110 public:
113 
116 
120 
123  : mr_model_part(model_part)
124  {
125  }
126 
129  {
130  }
131 
132 
136 
137  void operator()()
138  {
139  Execute();
140  }
141 
142 
146 
147  void Execute() override
148  {
149  KRATOS_TRY
150 
151  ProcessInfo& proc_info = mr_model_part.GetProcessInfo();
152  double dummy;
153  //first initialize the Mass force to the old value
154 
155  //set the Mass to the old value
156  for(ModelPart::NodesContainerType::iterator in = mr_model_part.NodesBegin() ;
157  in != mr_model_part.NodesEnd() ; ++in)
158  {
159  in->FastGetSolutionStepValue(NODAL_MASS)=0.0;
160  }
161 
162  for(ModelPart::ElementsContainerType::iterator im = mr_model_part.ElementsBegin() ;
163  im != mr_model_part.ElementsEnd() ; ++im)
164  {
165  im->Calculate(NODAL_MASS,dummy,proc_info);
166  }
167 
168  KRATOS_WATCH("Execute of Mass Calulate Process");
169 
170 
171 
172  KRATOS_CATCH("")
173  }
174 
175 
179 
180 
184 
185 
189 
191  std::string Info() const override
192  {
193  return "MassCalculateProcess";
194  }
195 
197  void PrintInfo(std::ostream& rOStream) const override
198  {
199  rOStream << "MassCalculateProcess";
200  }
201 
203  void PrintData(std::ostream& rOStream) const override
204  {
205  }
206 
207 
211 
212 
214 
215 protected:
218 
219 
223 
224 
228 
229 
233 
234 
238 
239 
243 
244 
248 
249 
251 
252 private:
255 
256 
260  ModelPart& mr_model_part;
261 
262 
266 
267 
271 
272 
276 
277 
281 
282 
286 
288 // MassCalculateProcess& operator=(MassCalculateProcess const& rOther);
289 
291 // MassCalculateProcess(MassCalculateProcess const& rOther);
292 
293 
295 
296 }; // Class MassCalculateProcess
297 
299 
302 
303 
307 
308 
310 inline std::istream& operator >> (std::istream& rIStream,
311  MassCalculateProcess& rThis);
312 
314 inline std::ostream& operator << (std::ostream& rOStream,
315  const MassCalculateProcess& rThis)
316 {
317  rThis.PrintInfo(rOStream);
318  rOStream << std::endl;
319  rThis.PrintData(rOStream);
320 
321  return rOStream;
322 }
324 
325 
326 } // namespace Kratos.
327 
328 #endif // KRATOS_MASS_CALCULATE_PROCESS_INCLUDED defined
329 
330 
Short class definition.
Definition: mass_calculate_process.h:109
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: mass_calculate_process.h:147
void operator()()
Definition: mass_calculate_process.h:137
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: mass_calculate_process.h:197
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: mass_calculate_process.h:203
KRATOS_CLASS_POINTER_DEFINITION(MassCalculateProcess)
Pointer definition of MassCalculateProcess.
MassCalculateProcess(ModelPart &model_part)
Default constructor.
Definition: mass_calculate_process.h:122
~MassCalculateProcess() override
Destructor.
Definition: mass_calculate_process.h:128
std::string Info() const override
Turn back information as a string.
Definition: mass_calculate_process.h:191
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
NodeIterator NodesBegin(IndexType ThisIndex=0)
Definition: model_part.h:487
ProcessInfo & GetProcessInfo()
Definition: model_part.h:1746
ElementIterator ElementsEnd(IndexType ThisIndex=0)
Definition: model_part.h:1179
NodeIterator NodesEnd(IndexType ThisIndex=0)
Definition: model_part.h:497
The base class for all processes in Kratos.
Definition: process.h:49
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_WATCH(variable)
Definition: define.h:806
#define KRATOS_TRY
Definition: define.h:109
im
Definition: GenerateCN.py:100
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
dummy
Definition: script.py:194