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_nodal_young_modulus_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_NODAL_YOUNG_MODULUS_PROCESS)
15 #define KRATOS_DAM_NODAL_YOUNG_MODULUS_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  "is_fixed" : false,
50  "Young_Modulus_1" : 10.0,
51  "Young_Modulus_2" : 60.0,
52  "Young_Modulus_3" : 50.0,
53  "Young_Modulus_4" : 70.0,
54  "interval":[
55  0.0,
56  0.0
57  ]
58  } )");
59 
60  // Some values need to be mandatorily prescribed since no meaningful default value exist. For this reason try accessing to them
61  // So that an error is thrown if they don't exist
62  rParameters["Young_Modulus_1"];
63  rParameters["variable_name"];
64  rParameters["model_part_name"];
65 
66  // Now validate agains defaults -- this also ensures no type mismatch
67  rParameters.ValidateAndAssignDefaults(default_parameters);
68 
69  mVariableName = rParameters["variable_name"].GetString();
70  mIsFixed = rParameters["is_fixed"].GetBool();
71  mYoung1 = rParameters["Young_Modulus_1"].GetDouble();
72  mYoung2 = rParameters["Young_Modulus_2"].GetDouble();
73  mYoung3 = rParameters["Young_Modulus_3"].GetDouble();
74  mYoung4 = rParameters["Young_Modulus_4"].GetDouble();
75 
76  KRATOS_CATCH("");
77  }
78 
80 
83 
84  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
85 
86  void Execute() override
87  {
88  }
89 
90  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
91 
92  void ExecuteInitialize() override
93  {
94 
95  KRATOS_TRY;
96 
98  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
99 
100  if (nnodes != 0)
101  {
102  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
103 
104 #pragma omp parallel for
105  for (int i = 0; i < nnodes; i++)
106  {
107  ModelPart::NodesContainerType::iterator it = it_begin + i;
108 
109  if (mIsFixed)
110  {
111  it->Fix(var);
112  }
113 
114  double Young = mYoung1 + (mYoung2 * it->X()) + (mYoung3 * it->Y()) + (mYoung4 * it->Z());
115 
116  if (Young <= 0.0)
117  {
118  it->FastGetSolutionStepValue(var) = 0.0;
119  }
120  else
121  it->FastGetSolutionStepValue(var) = Young;
122  }
123  }
124 
125  KRATOS_CATCH("");
126  }
127 
128  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
129 
131  {
132 
133  KRATOS_TRY;
134 
136  const int nnodes = mrModelPart.GetMesh(0).Nodes().size();
137 
138  if (nnodes != 0)
139  {
140  ModelPart::NodesContainerType::iterator it_begin = mrModelPart.GetMesh(0).NodesBegin();
141 
142 #pragma omp parallel for
143  for (int i = 0; i < nnodes; i++)
144  {
145  ModelPart::NodesContainerType::iterator it = it_begin + i;
146 
147  if (mIsFixed)
148  {
149  it->Fix(var);
150  }
151 
152  double Young = mYoung1 + (mYoung2 * it->X()) + (mYoung3 * it->Y()) + (mYoung4 * it->Z());
153 
154  if (Young <= 0.0)
155  {
156  it->FastGetSolutionStepValue(var) = 0.0;
157  }
158  else
159  it->FastGetSolutionStepValue(var) = Young;
160  }
161  }
162 
163  KRATOS_CATCH("");
164  }
165 
167  std::string Info() const override
168  {
169  return "DamNodalYoungModulusProcess";
170  }
171 
173  void PrintInfo(std::ostream &rOStream) const override
174  {
175  rOStream << "DamNodalYoungModulusProcess";
176  }
177 
179  void PrintData(std::ostream &rOStream) const override
180  {
181  }
182 
184 
185  protected:
187 
189  std::string mVariableName;
190  bool mIsFixed;
191  double mYoung1;
192  double mYoung2;
193  double mYoung3;
194  double mYoung4;
195 
196  //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
197 
198  private:
201 
202 }; //Class
203 
205 inline std::istream &operator>>(std::istream &rIStream,
207 
209 inline std::ostream &operator<<(std::ostream &rOStream,
210  const DamNodalYoungModulusProcess &rThis)
211 {
212  rThis.PrintInfo(rOStream);
213  rOStream << std::endl;
214  rThis.PrintData(rOStream);
215 
216  return rOStream;
217 }
218 
219 } /* namespace Kratos.*/
220 
221 #endif /* KRATOS_DAM_NODAL_YOUNG_MODULUS_PROCESS defined */
Definition: dam_nodal_young_modulus_process.hpp:31
bool mIsFixed
Definition: dam_nodal_young_modulus_process.hpp:190
DamNodalYoungModulusProcess(ModelPart &rModelPart, Parameters &rParameters)
Constructor.
Definition: dam_nodal_young_modulus_process.hpp:39
double mYoung2
Definition: dam_nodal_young_modulus_process.hpp:192
void Execute() override
Execute method is used to execute the Process algorithms.
Definition: dam_nodal_young_modulus_process.hpp:86
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: dam_nodal_young_modulus_process.hpp:179
double mYoung3
Definition: dam_nodal_young_modulus_process.hpp:193
std::string Info() const override
Turn back information as a string.
Definition: dam_nodal_young_modulus_process.hpp:167
void ExecuteInitializeSolutionStep() override
This function will be executed at every time step BEFORE performing the solve phase.
Definition: dam_nodal_young_modulus_process.hpp:130
virtual ~DamNodalYoungModulusProcess()
Destructor.
Definition: dam_nodal_young_modulus_process.hpp:82
std::string mVariableName
Definition: dam_nodal_young_modulus_process.hpp:189
KRATOS_CLASS_POINTER_DEFINITION(DamNodalYoungModulusProcess)
void ExecuteInitialize() override
This function is designed for being called at the beginning of the computations right after reading t...
Definition: dam_nodal_young_modulus_process.hpp:92
ModelPart & mrModelPart
Member Variables.
Definition: dam_nodal_young_modulus_process.hpp:188
double mYoung4
Definition: dam_nodal_young_modulus_process.hpp:194
double mYoung1
Definition: dam_nodal_young_modulus_process.hpp:191
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: dam_nodal_young_modulus_process.hpp:173
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
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
bool GetBool() const
This method returns the boolean contained in the current Parameter.
Definition: kratos_parameters.cpp:675
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
Young
Definition: isotropic_damage_automatic_differentiation.py:135
int nnodes
Definition: sensitivityMatrix.py:24
integer i
Definition: TensorModule.f:17