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_t_sol_air_heat_flux_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_T_SOL_AIR_HEAT_FLUX_PROCESS)
15 #define KRATOS_DAM_T_SOL_AIR_HEAT_FLUX_PROCESS
16 
17 #include <cmath>
18 
19 // Project includes
20 #include "includes/kratos_flags.h"
22 #include "processes/process.h"
23 
24 // Application include
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  "h_0" : 0.0,
52  "ambient_temperature" : 0.0,
53  "table_ambient_temperature" : 0,
54  "emisivity" : 0.0,
55  "delta_R" : 0.0,
56  "absorption_index" : 0.0,
57  "total_insolation" : 0.0,
58  "interval":[
59  0.0,
60  0.0
61  ]
62  } )");
63 
64  // Some values need to be mandatorily prescribed since no meaningful default value exist. For this reason try accessing to them
65  // So that an error is thrown if they don't exist
66  rParameters["h_0"];
67  rParameters["delta_R"];
68  rParameters["absorption_index"];
69 
70  // Now validate agains defaults -- this also ensures no type mismatch
71  rParameters.ValidateAndAssignDefaults(default_parameters);
72 
73  mVariableName = rParameters["variable_name"].GetString();
74  mH0 = rParameters["h_0"].GetDouble();
75  mAmbientTemperature = rParameters["ambient_temperature"].GetDouble();
76  mEmisivity = rParameters["emisivity"].GetDouble();
77  mDeltaR = rParameters["delta_R"].GetDouble();
78  mAbsorption_index = rParameters["absorption_index"].GetDouble();
79  mTotalInsolation = rParameters["total_insolation"].GetDouble();
80 
81  mTimeUnitConverter = mrModelPart.GetProcessInfo()[TIME_UNIT_CONVERTER];
82  mTableId = rParameters["table_ambient_temperature"].GetInt();
83 
84  if (mTableId != 0)
86 
87  KRATOS_CATCH("");
88  }
89 
91 
94 
95  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
96 
97  void Execute() override
98  {
99  }
100 
101  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
102 
103 
104  void ExecuteInitialize() override
105  {
106 
107  KRATOS_TRY;
108 
109  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
111 
112  // Computing the t_soil_air according to t_sol_air criteria
114 
115  if (nnodes != 0)
116  {
117  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
118 
119  #pragma omp parallel for
120  for (int i = 0; i < nnodes; i++)
121  {
122  ModelPart::NodesContainerType::iterator it = it_begin + i;
123 
124  const double temp_current = it->FastGetSolutionStepValue(TEMPERATURE);
125  const double heat_flux = mH0 * (t_sol_air - temp_current);
126  it->FastGetSolutionStepValue(var) = heat_flux;
127  }
128  }
129 
130  KRATOS_CATCH("");
131  }
132 
133  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
134 
136  {
137 
138  KRATOS_TRY;
139 
140  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
142 
143  // Getting the values of table in case that it exist
144  if (mTableId != 0)
145  {
146  double time = mrModelPart.GetProcessInfo()[TIME];
148  mAmbientTemperature = mpTable->GetValue(time);
149  }
150 
151  // Computing the t_soil_air according to t_sol_air criteria
153 
154  if (nnodes != 0)
155  {
156  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
157 
158  #pragma omp parallel for
159  for (int i = 0; i < nnodes; i++)
160  {
161  ModelPart::NodesContainerType::iterator it = it_begin + i;
162 
163  const double temp_current = it->FastGetSolutionStepValue(TEMPERATURE);
164  const double heat_flux = mH0 * (t_sol_air - temp_current);
165  it->FastGetSolutionStepValue(var) = heat_flux;
166  }
167  }
168 
169  KRATOS_CATCH("");
170  }
171 
173 
175  std::string Info() const override
176  {
177  return "DamTSolAirHeatFluxProcess";
178  }
179 
181  void PrintInfo(std::ostream &rOStream) const override
182  {
183  rOStream << "DamTSolAirHeatFluxProcess";
184  }
185 
187  void PrintData(std::ostream &rOStream) const override
188  {
189  }
190 
192 
193  protected:
196  std::string mVariableName;
197  double mH0;
199  double mEmisivity;
200  double mDeltaR;
204  TableType::Pointer mpTable;
205  int mTableId;
206 
207  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
208 
209  private:
211  DamTSolAirHeatFluxProcess &operator=(DamTSolAirHeatFluxProcess const &rOther);
212 
213 }; //Class
214 
216 inline std::istream &operator>>(std::istream &rIStream,
218 
220 inline std::ostream &operator<<(std::ostream &rOStream,
221  const DamTSolAirHeatFluxProcess &rThis)
222 {
223  rThis.PrintInfo(rOStream);
224  rOStream << std::endl;
225  rThis.PrintData(rOStream);
226 
227  return rOStream;
228 }
229 
230 } /* namespace Kratos.*/
231 
232 #endif /* KRATOS_DAM_T_SOL_AIR_HEAT_FLUX_PROCESS defined */
Definition: dam_t_sol_air_heat_flux_process.hpp:31
void ExecuteInitializeSolutionStep() override
This function will be executed at every time step BEFORE performing the solve phase.
Definition: dam_t_sol_air_heat_flux_process.hpp:135
int mTableId
Definition: dam_t_sol_air_heat_flux_process.hpp:205
std::string mVariableName
Definition: dam_t_sol_air_heat_flux_process.hpp:196
double mTotalInsolation
Definition: dam_t_sol_air_heat_flux_process.hpp:202
double mAbsorption_index
Definition: dam_t_sol_air_heat_flux_process.hpp:201
ModelPart & mrModelPart
Member Variables.
Definition: dam_t_sol_air_heat_flux_process.hpp:195
Table< double, double > TableType
Definition: dam_t_sol_air_heat_flux_process.hpp:36
std::string Info() const override
Turn back information as a string.
Definition: dam_t_sol_air_heat_flux_process.hpp:175
virtual ~DamTSolAirHeatFluxProcess()
Destructor.
Definition: dam_t_sol_air_heat_flux_process.hpp:93
double mH0
Definition: dam_t_sol_air_heat_flux_process.hpp:197
void ExecuteInitialize() override
This function is designed for being called at the beginning of the computations right after reading t...
Definition: dam_t_sol_air_heat_flux_process.hpp:104
double mEmisivity
Definition: dam_t_sol_air_heat_flux_process.hpp:199
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: dam_t_sol_air_heat_flux_process.hpp:187
KRATOS_CLASS_POINTER_DEFINITION(DamTSolAirHeatFluxProcess)
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: dam_t_sol_air_heat_flux_process.hpp:181
DamTSolAirHeatFluxProcess(ModelPart &rModelPart, Parameters &rParameters)
Constructor.
Definition: dam_t_sol_air_heat_flux_process.hpp:41
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: dam_t_sol_air_heat_flux_process.hpp:97
double mDeltaR
Definition: dam_t_sol_air_heat_flux_process.hpp:200
double mAmbientTemperature
Definition: dam_t_sol_air_heat_flux_process.hpp:198
TableType::Pointer mpTable
Definition: dam_t_sol_air_heat_flux_process.hpp:204
double mTimeUnitConverter
Definition: dam_t_sol_air_heat_flux_process.hpp:203
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
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