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_grouting_reference_temperature_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: D.J. Vicente
11 //
12 //
13 
14 #if !defined(KRATOS_DAM_GROUTING_REFERENCE_TEMPERATURE_PROCESS )
15 #define KRATOS_DAM_GROUTING_REFERENCE_TEMPERATURE_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:
34 
36 
37 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
38 
41  ) : Process(Flags()) , mrModelPart(rModelPart)
42  {
44 
45  //only include validation with c++11 since raw_literals do not exist in c++03
46  Parameters default_parameters( R"(
47  {
48  "model_part_name":"PLEASE_CHOOSE_MODEL_PART_NAME",
49  "variable_name" : "PLEASE_PRESCRIBE_VARIABLE_NAME",
50  "initial_value" : 0.0,
51  "time_grouting" : 0.0,
52  "interval":[
53  0.0,
54  0.0
55  ]
56  } )" );
57 
58  // Some values need to be mandatorily prescribed since no meaningful default value exist. For this reason try accessing to them
59  // So that an error is thrown if they don't exist
60  rParameters["variable_name"];
61  rParameters["model_part_name"];
62 
63  // Now validate agains defaults -- this also ensures no type mismatch
64  rParameters.ValidateAndAssignDefaults(default_parameters);
65 
66  mVariableName = rParameters["variable_name"].GetString();
67  mInitialValue = rParameters["initial_value"].GetDouble();
68  mTimeGrouting = rParameters["time_grouting"].GetDouble();
69  mTimeUnitConverter = mrModelPart.GetProcessInfo()[TIME_UNIT_CONVERTER];
70 
71  KRATOS_CATCH("");
72  }
73 
75 
78 
79 
80 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
81 
82  void Execute() override
83  {
84  }
85 
86  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
87 
88 
89  void ExecuteInitialize() override
90  {
91 
92  KRATOS_TRY;
93 
95 
96  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
97 
98  if(nnodes != 0)
99  {
100  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
101 
102  double time = mrModelPart.GetProcessInfo()[TIME];
104 
105  #pragma omp parallel for
106  for(int i = 0; i<nnodes; i++)
107  {
108  ModelPart::NodesContainerType::iterator it = it_begin + i;
109  it->FastGetSolutionStepValue(var) = mInitialValue;
110  }
111  }
112 
113  KRATOS_CATCH("");
114  }
115 
116 
117 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
118 
120  {
121 
122  KRATOS_TRY;
123 
125 
126  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
127 
128  if(nnodes != 0)
129  {
130  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
131 
132  double time = mrModelPart.GetProcessInfo()[TIME];
134 
135  if (time == mTimeGrouting)
136  {
137  #pragma omp parallel for
138  for(int i = 0; i<nnodes; i++)
139  {
140  ModelPart::NodesContainerType::iterator it = it_begin + i;
141  const double current_temp = it->FastGetSolutionStepValue(TEMPERATURE);
142  it->FastGetSolutionStepValue(var) = current_temp;
143 
144  }
145  }
146  }
147 
148  KRATOS_CATCH("");
149  }
150 
151 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
152 
154  std::string Info() const override
155  {
156  return "DamGroutingReferenceTemperatureProcess";
157  }
158 
160  void PrintInfo(std::ostream& rOStream) const override
161  {
162  rOStream << "DamGroutingReferenceTemperatureProcess";
163  }
164 
166  void PrintData(std::ostream& rOStream) const override
167  {
168  }
169 
171 
172 protected:
173 
175 
177  std::string mVariableName;
181 
182 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
183 
184 private:
185 
188 
189 };//Class
190 
191 
193 inline std::istream& operator >> (std::istream& rIStream,
195 
197 inline std::ostream& operator << (std::ostream& rOStream,
199 {
200  rThis.PrintInfo(rOStream);
201  rOStream << std::endl;
202  rThis.PrintData(rOStream);
203 
204  return rOStream;
205 }
206 
207 } /* namespace Kratos.*/
208 
209 #endif /* KRATOS_DAM_GROUTING_REFERENCE_TEMPERATURE_PROCESS defined */
210 
Definition: dam_grouting_reference_temperature_process.hpp:31
double mTimeUnitConverter
Definition: dam_grouting_reference_temperature_process.hpp:180
virtual ~DamGroutingReferenceTemperatureProcess()
Destructor.
Definition: dam_grouting_reference_temperature_process.hpp:77
DamGroutingReferenceTemperatureProcess(ModelPart &rModelPart, Parameters &rParameters)
Constructor.
Definition: dam_grouting_reference_temperature_process.hpp:40
void ExecuteInitialize() override
This function is designed for being called at the beginning of the computations right after reading t...
Definition: dam_grouting_reference_temperature_process.hpp:89
std::string mVariableName
Definition: dam_grouting_reference_temperature_process.hpp:177
double mInitialValue
Definition: dam_grouting_reference_temperature_process.hpp:178
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: dam_grouting_reference_temperature_process.hpp:160
void ExecuteFinalizeSolutionStep() override
This function will be executed at every time step AFTER performing the solve phase.
Definition: dam_grouting_reference_temperature_process.hpp:119
double mTimeGrouting
Definition: dam_grouting_reference_temperature_process.hpp:179
std::string Info() const override
Turn back information as a string.
Definition: dam_grouting_reference_temperature_process.hpp:154
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: dam_grouting_reference_temperature_process.hpp:166
KRATOS_CLASS_POINTER_DEFINITION(DamGroutingReferenceTemperatureProcess)
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: dam_grouting_reference_temperature_process.hpp:82
ModelPart & mrModelPart
Member Variables.
Definition: dam_grouting_reference_temperature_process.hpp:176
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
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
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
#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