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_noorzai_heat_source_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_NOORZAI_HEAT_SOURCE_PROCESS)
15 #define KRATOS_DAM_NOORZAI_HEAT_SOURCE_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 
36  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
37 
40  Parameters &rParameters) : Process(Flags()), mrModelPart(rModelPart)
41  {
43 
44  //only include validation with c++11 since raw_literals do not exist in c++03
45  Parameters default_parameters(R"(
46  {
47  "model_part_name":"PLEASE_CHOOSE_MODEL_PART_NAME",
48  "variable_name": "PLEASE_PRESCRIBE_VARIABLE_NAME",
49  "density" : 0.0,
50  "specific_heat" : 0.0,
51  "t_max" : 0.0,
52  "alpha" : 0.0,
53  "interval":[
54  0.0,
55  0.0
56  ]
57  } )");
58 
59  // Some values need to be mandatorily prescribed since no meaningful default value exist. For this reason try accessing to them
60  // So that an error is thrown if they don't exist
61  rParameters["t_max"];
62  rParameters["alpha"];
63  rParameters["specific_heat"];
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  mDensity = rParameters["density"].GetDouble();
70  mSpecificHeat = rParameters["specific_heat"].GetDouble();
71  mTMax = rParameters["t_max"].GetDouble();
72  mAlpha = rParameters["alpha"].GetDouble();
73 
74  KRATOS_CATCH("");
75  }
76 
78 
81 
82  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
83 
84  void Execute() override
85  {
86  }
87 
88  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
89 
90  void ExecuteInitialize() override
91  {
92  KRATOS_TRY;
93 
94  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
96 
97  const double time = mrModelPart.GetProcessInfo()[TIME];
98  const double delta_time = mrModelPart.GetProcessInfo()[DELTA_TIME];
99 
100  double value = mDensity * mSpecificHeat * mAlpha * mTMax * (exp(-mAlpha * time + 0.5 * delta_time));
101 
102  if (nnodes != 0)
103  {
104  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
105 
106  #pragma omp parallel for
107  for (int i = 0; i < nnodes; i++)
108  {
109  ModelPart::NodesContainerType::iterator it = it_begin + i;
110  it->FastGetSolutionStepValue(var) = value;
111  }
112  }
113  KRATOS_CATCH("");
114  }
115 
116  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
117 
119  {
120  KRATOS_TRY;
121 
122  this->ExecuteInitialize();
123 
124  KRATOS_CATCH("");
125  }
126 
128 
130  std::string Info() const override
131  {
132  return "DamNoorzaiHeatFluxProcess";
133  }
134 
136  void PrintInfo(std::ostream &rOStream) const override
137  {
138  rOStream << "DamNoorzaiHeatFluxProcess";
139  }
140 
142  void PrintData(std::ostream &rOStream) const override
143  {
144  }
145 
147 
148  protected:
151  std::string mVariableName;
152  double mDensity;
154  double mAlpha;
155  double mTMax;
156 
157  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
158 
159  private:
161  DamNoorzaiHeatFluxProcess &operator=(DamNoorzaiHeatFluxProcess const &rOther);
162 
163 }; //Class
164 
166 inline std::istream &operator>>(std::istream &rIStream,
168 
170 inline std::ostream &operator<<(std::ostream &rOStream,
171  const DamNoorzaiHeatFluxProcess &rThis)
172 {
173  rThis.PrintInfo(rOStream);
174  rOStream << std::endl;
175  rThis.PrintData(rOStream);
176 
177  return rOStream;
178 }
179 
180 } /* namespace Kratos.*/
181 
182 #endif /* KRATOS_DAM_NOORZAI_HEAT_SOURCE_PROCESS defined */
Definition: dam_noorzai_heat_source_process.hpp:31
double mAlpha
Definition: dam_noorzai_heat_source_process.hpp:154
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: dam_noorzai_heat_source_process.hpp:84
virtual ~DamNoorzaiHeatFluxProcess()
Destructor.
Definition: dam_noorzai_heat_source_process.hpp:80
DamNoorzaiHeatFluxProcess(ModelPart &rModelPart, Parameters &rParameters)
Constructor.
Definition: dam_noorzai_heat_source_process.hpp:39
std::string mVariableName
Definition: dam_noorzai_heat_source_process.hpp:151
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: dam_noorzai_heat_source_process.hpp:136
std::string Info() const override
Turn back information as a string.
Definition: dam_noorzai_heat_source_process.hpp:130
double mSpecificHeat
Definition: dam_noorzai_heat_source_process.hpp:153
double mDensity
Definition: dam_noorzai_heat_source_process.hpp:152
void ExecuteInitialize() override
This function is designed for being called at the beginning of the computations right after reading t...
Definition: dam_noorzai_heat_source_process.hpp:90
KRATOS_CLASS_POINTER_DEFINITION(DamNoorzaiHeatFluxProcess)
void ExecuteInitializeSolutionStep() override
This function will be executed at every time step BEFORE performing the solve phase.
Definition: dam_noorzai_heat_source_process.hpp:118
double mTMax
Definition: dam_noorzai_heat_source_process.hpp:155
ModelPart & mrModelPart
Member Variables.
Definition: dam_noorzai_heat_source_process.hpp:150
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: dam_noorzai_heat_source_process.hpp:142
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
delta_time
Definition: generate_frictional_mortar_condition.py:130
int nnodes
Definition: sensitivityMatrix.py:24
integer i
Definition: TensorModule.f:17