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_hydro_condition_load_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_HYDRO_CONDITION_LOAD_PROCESS)
15 #define KRATOS_DAM_HYDRO_CONDITION_LOAD_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  "Modify" : true,
52  "Gravity_Direction" : "Y",
53  "Reservoir_Bottom_Coordinate_in_Gravity_Direction" : 0.0,
54  "Spe_weight" : 0.0,
55  "Water_level" : 0.0,
56  "Water_Table" : 0,
57  "interval":[
58  0.0,
59  0.0
60  ]
61  } )");
62 
63  // Some values need to be mandatorily prescribed since no meaningful default value exist. For this reason try accessing to them
64  // So that an error is thrown if they don't exist
65  rParameters["Reservoir_Bottom_Coordinate_in_Gravity_Direction"];
66  rParameters["variable_name"];
67  rParameters["model_part_name"];
68 
69  // Now validate agains defaults -- this also ensures no type mismatch
70  rParameters.ValidateAndAssignDefaults(default_parameters);
71 
72  mVariableName = rParameters["variable_name"].GetString();
73  mGravityDirection = rParameters["Gravity_Direction"].GetString();
74  mReferenceCoordinate = rParameters["Reservoir_Bottom_Coordinate_in_Gravity_Direction"].GetDouble();
75  mSpecific = rParameters["Spe_weight"].GetDouble();
76  mWaterLevel = rParameters["Water_level"].GetDouble();
77 
78  mTimeUnitConverter = mrModelPart.GetProcessInfo()[TIME_UNIT_CONVERTER];
79  mTableId = rParameters["Water_Table"].GetInt();
80 
81  if (mTableId != 0)
83 
84  KRATOS_CATCH("");
85  }
86 
88 
91 
92  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
93 
94  void Execute() override
95  {
96  }
97 
98  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
99 
100  void ExecuteInitialize() override
101  {
102  KRATOS_TRY;
103 
105  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
106  int direction;
107 
108  if (mGravityDirection == "X")
109  direction = 0;
110  else if (mGravityDirection == "Y")
111  direction = 1;
112  else
113  direction = 2;
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  double pressure = (mSpecific * (mWaterLevel - it->Coordinates()[direction]));
125 
126  if (pressure > 0.0)
127  {
128  it->FastGetSolutionStepValue(var) = pressure;
129  }
130  else
131  {
132  it->FastGetSolutionStepValue(var) = 0.0;
133  }
134  }
135  }
136 
137  KRATOS_CATCH("");
138  }
139 
140  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
141 
143  {
144 
145  KRATOS_TRY;
146 
148 
149  // Getting the values of table in case that it exist
150  if (mTableId != 0)
151  {
152  double time = mrModelPart.GetProcessInfo()[TIME];
154  mWaterLevel = mpTable->GetValue(time);
155  }
156 
157  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
158 
159  int direction;
160 
161  if (mGravityDirection == "X")
162  direction = 0;
163  else if (mGravityDirection == "Y")
164  direction = 1;
165  else
166  direction = 2;
167 
168  if (nnodes != 0)
169  {
170  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
171 
172 #pragma omp parallel for
173  for (int i = 0; i < nnodes; i++)
174  {
175  ModelPart::NodesContainerType::iterator it = it_begin + i;
176 
177  double pressure = (mSpecific * (mWaterLevel - it->Coordinates()[direction]));
178 
179  if (pressure > 0.0)
180  {
181  it->FastGetSolutionStepValue(var) = pressure;
182  }
183  else
184  {
185  it->FastGetSolutionStepValue(var) = 0.0;
186  }
187  }
188  }
189 
190  KRATOS_CATCH("");
191  }
192 
194  std::string Info() const override
195  {
196  return "DamHydroConditionLoadProcess";
197  }
198 
200  void PrintInfo(std::ostream &rOStream) const override
201  {
202  rOStream << "DamHydroConditionLoadProcess";
203  }
204 
206  void PrintData(std::ostream &rOStream) const override
207  {
208  }
209 
211 
212  protected:
214 
216  std::string mVariableName;
217  std::string mGravityDirection;
219  double mSpecific;
220  double mWaterLevel;
222  TableType::Pointer mpTable;
223  int mTableId;
224 
225  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
226 
227  private:
230 
231 }; //Class
232 
234 inline std::istream &operator>>(std::istream &rIStream,
236 
238 inline std::ostream &operator<<(std::ostream &rOStream,
239  const DamHydroConditionLoadProcess &rThis)
240 {
241  rThis.PrintInfo(rOStream);
242  rOStream << std::endl;
243  rThis.PrintData(rOStream);
244 
245  return rOStream;
246 }
247 
248 } /* namespace Kratos.*/
249 
250 #endif /* KRATOS_DAM_HYDRO_CONDITION_LOAD_PROCESS defined */
Definition: dam_hydro_condition_load_process.hpp:31
double mWaterLevel
Definition: dam_hydro_condition_load_process.hpp:220
double mTimeUnitConverter
Definition: dam_hydro_condition_load_process.hpp:221
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: dam_hydro_condition_load_process.hpp:200
Table< double, double > TableType
Definition: dam_hydro_condition_load_process.hpp:36
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: dam_hydro_condition_load_process.hpp:94
int mTableId
Definition: dam_hydro_condition_load_process.hpp:223
std::string Info() const override
Turn back information as a string.
Definition: dam_hydro_condition_load_process.hpp:194
ModelPart & mrModelPart
Member Variables.
Definition: dam_hydro_condition_load_process.hpp:215
TableType::Pointer mpTable
Definition: dam_hydro_condition_load_process.hpp:222
void ExecuteInitializeSolutionStep() override
This function will be executed at every time step BEFORE performing the solve phase.
Definition: dam_hydro_condition_load_process.hpp:142
std::string mVariableName
Definition: dam_hydro_condition_load_process.hpp:216
virtual ~DamHydroConditionLoadProcess()
Destructor.
Definition: dam_hydro_condition_load_process.hpp:90
DamHydroConditionLoadProcess(ModelPart &rModelPart, Parameters &rParameters)
Constructor.
Definition: dam_hydro_condition_load_process.hpp:41
std::string mGravityDirection
Definition: dam_hydro_condition_load_process.hpp:217
double mSpecific
Definition: dam_hydro_condition_load_process.hpp:219
void ExecuteInitialize() override
This function is designed for being called at the beginning of the computations right after reading t...
Definition: dam_hydro_condition_load_process.hpp:100
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: dam_hydro_condition_load_process.hpp:206
KRATOS_CLASS_POINTER_DEFINITION(DamHydroConditionLoadProcess)
double mReferenceCoordinate
Definition: dam_hydro_condition_load_process.hpp:218
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
float pressure
Definition: test_pureconvectionsolver_benchmarking.py:101