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_dof_process.hpp
Go to the documentation of this file.
1 //
2 // Project Name: KratosSolidMechanicsApplication $
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_DOF_PROCESS_H_INCLUDED)
11 #define KRATOS_FREE_SCALAR_DOF_PROCESS_H_INCLUDED
12 
13 
14 // System includes
15 
16 // External includes
17 
18 // Project includes
19 #include "includes/model_part.h"
21 #include "processes/process.h"
22 
23 
24 namespace Kratos
25 {
26 
29 
31 
33 class FreeScalarDofProcess : public Process
34 {
35 public:
38 
41 
46  Parameters rParameters) : Process() , mrModelPart(model_part)
47  {
49 
50  Parameters default_parameters( R"(
51  {
52  "model_part_name":"PLEASE_CHOOSE_MODEL_PART_NAME",
53  "variable_name": "PLEASE_PRESCRIBE_VARIABLE_NAME"
54  } )" );
55 
56 
57  // Validate against defaults -- this ensures no type mismatch
58  rParameters.ValidateAndAssignDefaults(default_parameters);
59 
60  mvariable_name = rParameters["variable_name"].GetString();
61 
62  if( KratosComponents< Variable<double> >::Has( mvariable_name ) ) //case of double variable
63  {
64  if( model_part.GetNodalSolutionStepVariablesList().Has( KratosComponents< Variable<double> >::Get( mvariable_name ) ) == false )
65  {
66  KRATOS_THROW_ERROR(std::runtime_error,"trying to set a variable that is not in the model_part - variable name is ",mvariable_name);
67  }
68 
69  }
70  else if( KratosComponents< Variable<int> >::Has( mvariable_name ) ) //case of int variable
71  {
72  if( model_part.GetNodalSolutionStepVariablesList().Has( KratosComponents< Variable<int> >::Get( mvariable_name ) ) == false )
73  {
74  KRATOS_THROW_ERROR(std::runtime_error,"trying to set a variable that is not in the model_part - variable name is ",mvariable_name);
75  }
76 
77  }
78  else if( KratosComponents< Variable<bool> >::Has( mvariable_name ) ) //case of bool variable
79  {
80  if( model_part.GetNodalSolutionStepVariablesList().Has(KratosComponents< Variable<bool> >::Get( mvariable_name ) ) == false )
81  {
82  KRATOS_THROW_ERROR(std::runtime_error,"trying to set a variable that is not in the model_part - variable name is ",mvariable_name);
83  }
84  }
85 
86  KRATOS_CATCH("");
87  }
88 
90  const Variable<double>& rVariable) : Process(), mrModelPart(model_part)
91  {
92  KRATOS_TRY;
93 
94 
95  if( model_part.GetNodalSolutionStepVariablesList().Has( rVariable ) == false )
96  {
97  KRATOS_THROW_ERROR(std::runtime_error,"trying to set a variable that is not in the model_part - variable name is ",rVariable);
98  }
99 
100  mvariable_name = rVariable.Name();
101 
102  KRATOS_CATCH("");
103  }
104 
105 
107  const Variable< int >& rVariable) : Process(), mrModelPart(model_part)
108  {
109  KRATOS_TRY;
110 
111  if( model_part.GetNodalSolutionStepVariablesList().Has( rVariable ) == false )
112  {
113  KRATOS_THROW_ERROR(std::runtime_error,"Trying to set a variable that is not in the model_part - variable name is ",rVariable);
114  }
115 
116  mvariable_name = rVariable.Name();
117 
118  KRATOS_CATCH("");
119  }
120 
122  const Variable< bool >& rVariable) : Process(), mrModelPart(model_part)
123  {
124  KRATOS_TRY;
125 
126 
127  if( model_part.GetNodalSolutionStepVariablesList().Has( rVariable ) == false )
128  {
129  KRATOS_THROW_ERROR(std::runtime_error,"Trying to set a variable that is not in the model_part - variable name is ",rVariable);
130  }
131 
132  mvariable_name = rVariable.Name();
133 
134  KRATOS_CATCH("");
135  }
136 
137 
139  ~FreeScalarDofProcess() override {}
140 
141 
145 
147  void operator()()
148  {
149  Execute();
150  }
151 
152 
156 
157 
159  void Execute() override
160  {
161 
162  KRATOS_TRY;
163 
164  if( KratosComponents< Variable<double> >::Has( mvariable_name ) ) //case of double variable
165  {
166  InternalFreeDof<>(KratosComponents< Variable<double> >::Get(mvariable_name));
167  }
168  else
169  {
170  KRATOS_THROW_ERROR(std::logic_error, "Not able to set the variable. Attempting to set variable:",mvariable_name);
171  }
172 
173  KRATOS_CATCH("");
174 
175  }
176 
179  void ExecuteInitialize() override
180  {
181  }
182 
186  {
187  }
188 
189 
192  {
193  }
194 
197  {
198  }
199 
200 
202  void ExecuteBeforeOutputStep() override
203  {
204  }
205 
206 
208  void ExecuteAfterOutputStep() override
209  {
210  }
211 
212 
215  void ExecuteFinalize() override
216  {
217  }
218 
219 
223 
224 
228 
229 
233 
235  std::string Info() const override
236  {
237  return "FreeScalarDofProcess";
238  }
239 
241  void PrintInfo(std::ostream& rOStream) const override
242  {
243  rOStream << "FreeScalarDofProcess";
244  }
245 
247  void PrintData(std::ostream& rOStream) const override
248  {
249  }
250 
251 
256 
257 protected:
258 
267 
270 
284 
285 private:
286 
292 
293  ModelPart& mrModelPart;
294  std::string mvariable_name;
295 
299 
300  template< class TVarType >
301  void InternalFreeDof(TVarType& rVar)
302  {
303  const int nnodes = mrModelPart.GetMesh().Nodes().size();
304 
305  if(nnodes != 0)
306  {
307  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.NodesBegin();
308 
309  #pragma omp parallel for
310  for(int i = 0; i<nnodes; i++)
311  {
312  ModelPart::NodesContainerType::iterator it = it_begin + i;
313  //it->pAddDof(rVar)->FreeDof();
314  it->Free(rVar);
315  }
316  }
317  }
318 
325 
327  FreeScalarDofProcess& operator=(FreeScalarDofProcess const& rOther);
328 
329 
339 
340 }; // Class FreeScalarDofProcess
341 
343 
346 
347 
351 
352 
354 inline std::istream& operator >> (std::istream& rIStream,
355  FreeScalarDofProcess& rThis);
356 
358 inline std::ostream& operator << (std::ostream& rOStream,
359  const FreeScalarDofProcess& rThis)
360 {
361  rThis.PrintInfo(rOStream);
362  rOStream << std::endl;
363  rThis.PrintData(rOStream);
364 
365  return rOStream;
366 }
368 
369 
370 } // namespace Kratos.
371 
372 #endif // KRATOS_FREE_SCALAR_DOF_PROCESS_H_INCLUDED defined
The base class for freeing scalar variable Dof or array_1d component Dof processes in Kratos.
Definition: free_scalar_dof_process.hpp:38
void operator()()
This operator is provided to call the process as a function and simply calls the Execute method.
Definition: free_scalar_dof_process.hpp:147
~FreeScalarDofProcess() override
Destructor.
Definition: free_scalar_dof_process.hpp:139
void ExecuteAfterOutputStep() override
this function will be executed at every time step AFTER writing the output
Definition: free_scalar_dof_process.hpp:208
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: free_scalar_dof_process.hpp:241
FreeScalarDofProcess(FreeScalarDofProcess const &rOther)
Copy constructor.
FreeScalarDofProcess(ModelPart &model_part, const Variable< double > &rVariable)
Definition: free_scalar_dof_process.hpp:89
FreeScalarDofProcess(ModelPart &model_part, const Variable< int > &rVariable)
Definition: free_scalar_dof_process.hpp:106
KRATOS_CLASS_POINTER_DEFINITION(FreeScalarDofProcess)
Pointer definition of FreeScalarDofProcess.
void Execute() override
Execute method is used to execute the FreeScalarDofProcess algorithms.
Definition: free_scalar_dof_process.hpp:163
FreeScalarDofProcess(ModelPart &model_part, Parameters rParameters)
Definition: free_scalar_dof_process.hpp:45
std::string Info() const override
Turn back information as a string.
Definition: free_scalar_dof_process.hpp:235
FreeScalarDofProcess(ModelPart &model_part, const Variable< bool > &rVariable)
Definition: free_scalar_dof_process.hpp:121
void ExecuteInitialize() override
Definition: free_scalar_dof_process.hpp:179
void ExecuteFinalize() override
Definition: free_scalar_dof_process.hpp:215
void ExecuteFinalizeSolutionStep() override
this function will be executed at every time step AFTER performing the solve phase
Definition: free_scalar_dof_process.hpp:196
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: free_scalar_dof_process.hpp:247
void ExecuteInitializeSolutionStep() override
this function will be executed at every time step BEFORE performing the solve phase
Definition: free_scalar_dof_process.hpp:191
void ExecuteBeforeSolutionLoop() override
Definition: free_scalar_dof_process.hpp:185
void ExecuteBeforeOutputStep() override
this function will be executed at every time step BEFORE writing the output
Definition: free_scalar_dof_process.hpp:202
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