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.
dam_fix_temperature_condition_process.hpp
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main authors: Lorenzo Gracia
11 //
12 //
13 
14 #if !defined(KRATOS_DAM_FIX_TEMPERATURE_CONDITION_PROCESS)
15 #define KRATOS_DAM_FIX_TEMPERATURE_CONDITION_PROCESS
16 
17 #include <cmath>
18 
19 // Project includes
20 #include "includes/kratos_flags.h"
22 #include "processes/process.h"
23 
24 // Application includes
26 
27 namespace Kratos
28 {
29 
31 {
32 
33  public:
35 
37 
38  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
39 
42  Parameters &rParameters) : Process(Flags()), mrModelPart(rModelPart)
43  {
45 
46  //only include validation with c++11 since raw_literals do not exist in c++03
47  Parameters default_parameters(R"(
48  {
49  "model_part_name" : "PLEASE_CHOOSE_MODEL_PART_NAME",
50  "variable_name" : "PLEASE_PRESCRIBE_VARIABLE_NAME",
51  "is_fixed" : false,
52  "value" : 0.0,
53  "table" : 0,
54  "interval":[
55  0.0,
56  0.0
57  ]
58  } )");
59 
60  // Some values need to be mandatorily prescribed since no meaningful default value exist. For this reason try accessing to them
61  // So that an error is thrown if they don't exist
62  rParameters["variable_name"];
63  rParameters["model_part_name"];
64 
65  // Now validate agains defaults -- this also ensures no type mismatch
66  rParameters.ValidateAndAssignDefaults(default_parameters);
67 
68  mVariableName = rParameters["variable_name"].GetString();
69  mIsFixed = rParameters["is_fixed"].GetBool();
70  mTemperature = rParameters["value"].GetDouble();
71 
72  mTimeUnitConverter = mrModelPart.GetProcessInfo()[TIME_UNIT_CONVERTER];
73  mTableId = rParameters["table"].GetInt();
74 
75  if (mTableId != 0)
77 
78  KRATOS_CATCH("");
79  }
80 
82 
85 
86  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
87  void Execute() override
88  {
89  }
90 
91  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
92 
93  void ExecuteInitialize() override
94  {
95 
96  KRATOS_TRY;
97 
99  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
100 
101  if (nnodes != 0)
102  {
103  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
104 
105  #pragma omp parallel for
106  for (int i = 0; i < nnodes; i++)
107  {
108  ModelPart::NodesContainerType::iterator it = it_begin + i;
109 
110  if (mIsFixed)
111  {
112  it->Fix(var);
113  }
114 
115  it->FastGetSolutionStepValue(var) = mTemperature;
116  }
117  }
118 
119  KRATOS_CATCH("");
120  }
121 
122  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
123 
125  {
126 
127  KRATOS_TRY;
128 
130 
131  // Getting the values of table in case that it exist
132  if (mTableId != 0)
133  {
134  double time = mrModelPart.GetProcessInfo()[TIME];
136  mTemperature = mpTable->GetValue(time);
137  }
138 
139  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
140 
141  if (nnodes != 0)
142  {
143  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
144 
145  #pragma omp parallel for
146  for (int i = 0; i < nnodes; i++)
147  {
148  ModelPart::NodesContainerType::iterator it = it_begin + i;
149 
150  if (mIsFixed)
151  {
152  it->Fix(var);
153  }
154 
155  it->FastGetSolutionStepValue(var) = mTemperature;
156  }
157  }
158 
159  KRATOS_CATCH("");
160  }
161 
162  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
163 
165  {
166 
167  KRATOS_TRY;
168 
170 
171  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
172 
173  if (nnodes != 0)
174  {
175 
176  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
177 
178  #pragma omp parallel for
179  for (int i = 0; i < nnodes; i++)
180  {
181  ModelPart::NodesContainerType::iterator it = it_begin + i;
182  it->Free(var);
183  }
184  }
185 
186  KRATOS_CATCH("");
187  }
188 
190  std::string Info() const override
191  {
192  return "FixTemperatureConditionProcess";
193  }
194 
196  void PrintInfo(std::ostream &rOStream) const override
197  {
198  rOStream << "FixTemperatureConditionProcess";
199  }
200 
202  void PrintData(std::ostream &rOStream) const override
203  {
204  }
205 
207 
208  protected:
210 
212  std::string mVariableName;
213  std::string mGravityDirection;
214  bool mIsFixed;
215  double mTemperature;
217  TableType::Pointer mpTable;
218  int mTableId;
219 
220  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
221 
222  private:
225 
226 }; //Class
227 
229 inline std::istream &operator>>(std::istream &rIStream,
231 
233 inline std::ostream &operator<<(std::ostream &rOStream,
235 {
236  rThis.PrintInfo(rOStream);
237  rOStream << std::endl;
238  rThis.PrintData(rOStream);
239 
240  return rOStream;
241 }
242 
243 } /* namespace Kratos.*/
244 
245 #endif /* KRATOS_DAM_FIX_TEMPERATURE_CONDITION_PROCESS defined */
Definition: dam_fix_temperature_condition_process.hpp:31
TableType::Pointer mpTable
Definition: dam_fix_temperature_condition_process.hpp:217
std::string mGravityDirection
Definition: dam_fix_temperature_condition_process.hpp:213
virtual ~DamFixTemperatureConditionProcess()
Destructor.
Definition: dam_fix_temperature_condition_process.hpp:84
std::string mVariableName
Definition: dam_fix_temperature_condition_process.hpp:212
void ExecuteInitialize() override
This function is designed for being called at the beginning of the computations right after reading t...
Definition: dam_fix_temperature_condition_process.hpp:93
int mTableId
Definition: dam_fix_temperature_condition_process.hpp:218
bool mIsFixed
Definition: dam_fix_temperature_condition_process.hpp:214
void ExecuteInitializeSolutionStep() override
This function will be executed at every time step BEFORE performing the solve phase.
Definition: dam_fix_temperature_condition_process.hpp:124
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: dam_fix_temperature_condition_process.hpp:196
double mTimeUnitConverter
Definition: dam_fix_temperature_condition_process.hpp:216
KRATOS_CLASS_POINTER_DEFINITION(DamFixTemperatureConditionProcess)
void ExecuteFinalizeSolutionStep() override
This function will be executed at every time step AFTER performing the solve phase.
Definition: dam_fix_temperature_condition_process.hpp:164
DamFixTemperatureConditionProcess(ModelPart &rModelPart, Parameters &rParameters)
Constructor.
Definition: dam_fix_temperature_condition_process.hpp:41
double mTemperature
Definition: dam_fix_temperature_condition_process.hpp:215
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: dam_fix_temperature_condition_process.hpp:202
Table< double, double > TableType
Definition: dam_fix_temperature_condition_process.hpp:36
ModelPart & mrModelPart
Member Variables.
Definition: dam_fix_temperature_condition_process.hpp:211
std::string Info() const override
Turn back information as a string.
Definition: dam_fix_temperature_condition_process.hpp:190
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: dam_fix_temperature_condition_process.hpp:87
Definition: flags.h:58
KratosComponents class encapsulates a lookup table for a family of classes in a generic way.
Definition: kratos_components.h:49
NodesContainerType & Nodes()
Definition: mesh.h:346
NodeIterator NodesBegin()
Definition: mesh.h:326
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
ProcessInfo & GetProcessInfo()
Definition: model_part.h:1746
TableType::Pointer pGetTable(IndexType TableId)
Definition: model_part.h:595
MeshType & GetMesh(IndexType ThisIndex=0)
Definition: model_part.h:1791
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
double GetDouble() const
This method returns the double contained in the current Parameter.
Definition: kratos_parameters.cpp:657
int GetInt() const
This method returns the integer contained in the current Parameter.
Definition: kratos_parameters.cpp:666
void ValidateAndAssignDefaults(const Parameters &rDefaultParameters)
This function is designed to verify that the parameters under testing match the form prescribed by th...
Definition: kratos_parameters.cpp:1306
std::string GetString() const
This method returns the string contained in the current Parameter.
Definition: kratos_parameters.cpp:684
bool GetBool() const
This method returns the boolean contained in the current Parameter.
Definition: kratos_parameters.cpp:675
size_type size() const
Returns the number of elements in the container.
Definition: pointer_vector_set.h:502
The base class for all processes in Kratos.
Definition: process.h:49
Definition: table.h:435
#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
time
Definition: face_heat.py:85
int nnodes
Definition: sensitivityMatrix.py:24
integer i
Definition: TensorModule.f:17