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.
free_scalar_pfem_dof_process.hpp
Go to the documentation of this file.
1 //
2 // Project Name: KratosPfemFluidDynamicsApplication $
3 // Created by: $Author: JMCarbonell $
4 // Last modified by: $Co-Author: $
5 // Date: $Date: August 2016 $
6 // Revision: $Revision: 0.0 $
7 //
8 //
9 
10 #if !defined(KRATOS_FREE_SCALAR_PFEM_DOF_PROCESS_H_INCLUDED)
11 #define KRATOS_FREE_SCALAR_PFEM_DOF_PROCESS_H_INCLUDED
12 
13 // System includes
14 
15 // External includes
16 
17 // Project includes
18 #include "includes/model_part.h"
20 #include "processes/process.h"
21 
22 namespace Kratos
23 {
24 
27 
29 
32 {
33 public:
36 
39 
44  Parameters rParameters) : Process(), mrModelPart(model_part)
45  {
47 
48  Parameters default_parameters(R"(
49  {
50  "model_part_name":"PLEASE_CHOOSE_MODEL_PART_NAME",
51  "variable_name": "PLEASE_PRESCRIBE_VARIABLE_NAME"
52  } )");
53 
54  // Validate against defaults -- this ensures no type mismatch
55  rParameters.ValidateAndAssignDefaults(default_parameters);
56 
57  mvariable_name = rParameters["variable_name"].GetString();
58 
59  if (KratosComponents<Variable<double>>::Has(mvariable_name)) //case of double variable
60  {
61  if (model_part.GetNodalSolutionStepVariablesList().Has(KratosComponents<Variable<double>>::Get(mvariable_name)) == false)
62  {
63  KRATOS_THROW_ERROR(std::runtime_error, "trying to set a variable that is not in the model_part - variable name is ", mvariable_name);
64  }
65  }
66  else if (KratosComponents<Variable<int>>::Has(mvariable_name)) //case of int variable
67  {
68  if (model_part.GetNodalSolutionStepVariablesList().Has(KratosComponents<Variable<int>>::Get(mvariable_name)) == false)
69  {
70  KRATOS_THROW_ERROR(std::runtime_error, "trying to set a variable that is not in the model_part - variable name is ", mvariable_name);
71  }
72  }
73  else if (KratosComponents<Variable<bool>>::Has(mvariable_name)) //case of bool variable
74  {
75  if (model_part.GetNodalSolutionStepVariablesList().Has(KratosComponents<Variable<bool>>::Get(mvariable_name)) == false)
76  {
77  KRATOS_THROW_ERROR(std::runtime_error, "trying to set a variable that is not in the model_part - variable name is ", mvariable_name);
78  }
79  }
80 
81  KRATOS_CATCH("");
82  }
83 
85  const Variable<double> &rVariable) : Process(), mrModelPart(model_part)
86  {
87  KRATOS_TRY;
88 
89  if (model_part.GetNodalSolutionStepVariablesList().Has(rVariable) == false)
90  {
91  KRATOS_THROW_ERROR(std::runtime_error, "trying to set a variable that is not in the model_part - variable name is ", rVariable);
92  }
93 
94  mvariable_name = rVariable.Name();
95 
96  KRATOS_CATCH("");
97  }
98 
100  const Variable<int> &rVariable) : Process(), mrModelPart(model_part)
101  {
102  KRATOS_TRY;
103 
104  if (model_part.GetNodalSolutionStepVariablesList().Has(rVariable) == false)
105  {
106  KRATOS_THROW_ERROR(std::runtime_error, "Trying to set a variable that is not in the model_part - variable name is ", rVariable);
107  }
108 
109  mvariable_name = rVariable.Name();
110 
111  KRATOS_CATCH("");
112  }
113 
115  const Variable<bool> &rVariable) : Process(), mrModelPart(model_part)
116  {
117  KRATOS_TRY;
118 
119  if (model_part.GetNodalSolutionStepVariablesList().Has(rVariable) == false)
120  {
121  KRATOS_THROW_ERROR(std::runtime_error, "Trying to set a variable that is not in the model_part - variable name is ", rVariable);
122  }
123 
124  mvariable_name = rVariable.Name();
125 
126  KRATOS_CATCH("");
127  }
128 
131 
135 
137  void operator()()
138  {
139  Execute();
140  }
141 
145 
147  void Execute() override
148  {
149 
150  KRATOS_TRY;
151 
152  if (KratosComponents<Variable<double>>::Has(mvariable_name)) //case of double variable
153  {
154  InternalFreeDof<>(KratosComponents<Variable<double>>::Get(mvariable_name));
155  }
156  else
157  {
158  KRATOS_THROW_ERROR(std::logic_error, "Not able to set the variable. Attempting to set variable:", mvariable_name);
159  }
160 
161  KRATOS_CATCH("");
162  }
163 
166  void ExecuteInitialize() override
167  {
168  }
169 
173  {
174  }
175 
178  {
179  }
180 
183  {
184  }
185 
187  void ExecuteBeforeOutputStep() override
188  {
189  }
190 
192  void ExecuteAfterOutputStep() override
193  {
194  }
195 
198  void ExecuteFinalize() override
199  {
200  }
201 
205 
209 
213 
215  std::string Info() const override
216  {
217  return "FreeScalarPfemDofProcess";
218  }
219 
221  void PrintInfo(std::ostream &rOStream) const override
222  {
223  rOStream << "FreeScalarPfemDofProcess";
224  }
225 
227  void PrintData(std::ostream &rOStream) const override
228  {
229  }
230 
235 
236 protected:
245 
248 
262 
263 private:
269 
270  ModelPart &mrModelPart;
271  std::string mvariable_name;
272 
276 
277  template <class TVarType>
278  void InternalFreeDof(TVarType &rVar)
279  {
280  const int nnodes = mrModelPart.GetMesh().Nodes().size();
281 
282  if (nnodes != 0)
283  {
284  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.NodesBegin();
285 
286 #pragma omp parallel for
287  for (int i = 0; i < nnodes; i++)
288  {
289  ModelPart::NodesContainerType::iterator it = it_begin + i;
290  //it->pAddDof(rVar)->FreeDof();
291  it->Free(rVar);
292  }
293  }
294  }
295 
302 
304  FreeScalarPfemDofProcess &operator=(FreeScalarPfemDofProcess const &rOther);
305 
315 
316 }; // Class FreeScalarPfemDofProcess
317 
319 
322 
326 
328 inline std::istream &operator>>(std::istream &rIStream,
329  FreeScalarPfemDofProcess &rThis);
330 
332 inline std::ostream &operator<<(std::ostream &rOStream,
333  const FreeScalarPfemDofProcess &rThis)
334 {
335  rThis.PrintInfo(rOStream);
336  rOStream << std::endl;
337  rThis.PrintData(rOStream);
338 
339  return rOStream;
340 }
342 
343 } // namespace Kratos.
344 
345 #endif // KRATOS_FREE_SCALAR_PFEM_DOF_PROCESS_H_INCLUDED defined
The base class for freeing scalar variable Dof or array_1d component Dof processes in Kratos.
Definition: free_scalar_pfem_dof_process.hpp:32
void ExecuteFinalizeSolutionStep() override
this function will be executed at every time step AFTER performing the solve phase
Definition: free_scalar_pfem_dof_process.hpp:182
void ExecuteAfterOutputStep() override
this function will be executed at every time step AFTER writing the output
Definition: free_scalar_pfem_dof_process.hpp:192
std::string Info() const override
Turn back information as a string.
Definition: free_scalar_pfem_dof_process.hpp:215
FreeScalarPfemDofProcess(FreeScalarPfemDofProcess const &rOther)
Copy constructor.
FreeScalarPfemDofProcess(ModelPart &model_part, const Variable< bool > &rVariable)
Definition: free_scalar_pfem_dof_process.hpp:114
FreeScalarPfemDofProcess(ModelPart &model_part, const Variable< int > &rVariable)
Definition: free_scalar_pfem_dof_process.hpp:99
void ExecuteInitialize() override
Definition: free_scalar_pfem_dof_process.hpp:166
KRATOS_CLASS_POINTER_DEFINITION(FreeScalarPfemDofProcess)
Pointer definition of FreeScalarPfemDofProcess.
void ExecuteInitializeSolutionStep() override
this function will be executed at every time step BEFORE performing the solve phase
Definition: free_scalar_pfem_dof_process.hpp:177
void ExecuteBeforeSolutionLoop() override
Definition: free_scalar_pfem_dof_process.hpp:172
FreeScalarPfemDofProcess(ModelPart &model_part, Parameters rParameters)
Definition: free_scalar_pfem_dof_process.hpp:43
void Execute() override
Execute method is used to execute the FreeScalarPfemDofProcess algorithms.
Definition: free_scalar_pfem_dof_process.hpp:147
void ExecuteBeforeOutputStep() override
this function will be executed at every time step BEFORE writing the output
Definition: free_scalar_pfem_dof_process.hpp:187
void ExecuteFinalize() override
Definition: free_scalar_pfem_dof_process.hpp:198
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: free_scalar_pfem_dof_process.hpp:221
~FreeScalarPfemDofProcess() override
Destructor.
Definition: free_scalar_pfem_dof_process.hpp:130
FreeScalarPfemDofProcess(ModelPart &model_part, const Variable< double > &rVariable)
Definition: free_scalar_pfem_dof_process.hpp:84
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: free_scalar_pfem_dof_process.hpp:227
void operator()()
This operator is provided to call the process as a function and simply calls the Execute method.
Definition: free_scalar_pfem_dof_process.hpp:137
KratosComponents class encapsulates a lookup table for a family of classes in a generic way.
Definition: kratos_components.h:49
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
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
The base class for all processes in Kratos.
Definition: process.h:49
const std::string & Name() const
Definition: variable_data.h:201
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
bool Has(const std::string &ModelerName)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:24
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
model_part
Definition: face_heat.py:14
int nnodes
Definition: sensitivityMatrix.py:24
integer i
Definition: TensorModule.f:17