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