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.
truss_embedded_edge_element.h
Go to the documentation of this file.
1 #if !defined(KRATOS_TRUSS_EMBEDDED_EDGE_ELEMENT_H_INCLUDED )
2 #define KRATOS_TRUSS_EMBEDDED_EDGE_ELEMENT_H_INCLUDED
3 
4 
5 // System includes
6 
7 // External includes
8 
9 // Project includes
10 #include "includes/define.h"
11 #include "includes/element.h"
12 #include "utilities/math_utils.h"
13 
14 // Application includes
16 
17 namespace Kratos
18 {
22 
24 class KRATOS_API(IGA_APPLICATION) TrussEmbeddedEdgeElement
25  : public Element
26 {
27 protected:
29  enum class ConfigurationType {
30  Current,
31  Reference
32  };
33 
34 public:
37 
40 
42  typedef std::size_t SizeType;
43  typedef std::size_t IndexType;
44 
45  // GometryType
47 
51 
54  IndexType NewId,
55  GeometryType::Pointer pGeometry)
56  : Element(NewId, pGeometry)
57  {};
58 
61  IndexType NewId,
62  GeometryType::Pointer pGeometry,
63  PropertiesType::Pointer pProperties)
64  : Element(NewId, pGeometry, pProperties)
65  {};
66 
69  : Element()
70  {};
71 
73  virtual ~TrussEmbeddedEdgeElement() = default;
74 
78 
80  Element::Pointer Create(
81  IndexType NewId,
82  GeometryType::Pointer pGeom,
83  PropertiesType::Pointer pProperties
84  ) const override
85  {
86  return Kratos::make_intrusive<TrussEmbeddedEdgeElement>(
87  NewId, pGeom, pProperties);
88  };
89 
91  Element::Pointer Create(
92  IndexType NewId,
93  NodesArrayType const& ThisNodes,
94  PropertiesType::Pointer pProperties
95  ) const override
96  {
97  return Kratos::make_intrusive< TrussEmbeddedEdgeElement >(
98  NewId, GetGeometry().Create(ThisNodes), pProperties);
99  };
100 
104 
112  VectorType& rRightHandSideVector,
113  const ProcessInfo& rCurrentProcessInfo) override
114  {
115  const SizeType number_of_nodes = GetGeometry().size();
116  const SizeType mat_size = number_of_nodes * 3;
117 
118  if (rRightHandSideVector.size() != mat_size)
119  rRightHandSideVector.resize(mat_size);
120  noalias(rRightHandSideVector) = ZeroVector(mat_size);
121 
122  MatrixType left_hand_side_matrix;
123 
124  CalculateAll(left_hand_side_matrix, rRightHandSideVector,
125  rCurrentProcessInfo, false, true);
126  }
127 
135  MatrixType& rLeftHandSideMatrix,
136  const ProcessInfo& rCurrentProcessInfo) override
137  {
138  const SizeType number_of_nodes = GetGeometry().size();
139  const SizeType mat_size = number_of_nodes * 3;
140 
141  VectorType right_hand_side_vector;
142 
143  if (rLeftHandSideMatrix.size1() != mat_size)
144  rLeftHandSideMatrix.resize(mat_size, mat_size);
145  noalias(rLeftHandSideMatrix) = ZeroMatrix(mat_size, mat_size);
146 
147  CalculateAll(rLeftHandSideMatrix, right_hand_side_vector,
148  rCurrentProcessInfo, true, false);
149  }
150 
160  MatrixType& rLeftHandSideMatrix,
161  VectorType& rRightHandSideVector,
162  const ProcessInfo& rCurrentProcessInfo) override
163  {
164  const SizeType number_of_nodes = GetGeometry().size();
165  const SizeType mat_size = number_of_nodes * 3;
166 
167  if (rRightHandSideVector.size() != mat_size)
168  rRightHandSideVector.resize(mat_size);
169  noalias(rRightHandSideVector) = ZeroVector(mat_size);
170 
171  if (rLeftHandSideMatrix.size1() != mat_size)
172  rLeftHandSideMatrix.resize(mat_size, mat_size);
173  noalias(rLeftHandSideMatrix) = ZeroMatrix(mat_size, mat_size);
174 
175  CalculateAll(rLeftHandSideMatrix, rRightHandSideVector,
176  rCurrentProcessInfo, true, true);
177  }
178 
184  void CalculateMassMatrix(
185  MatrixType& rMassMatrix,
186  const ProcessInfo& rCurrentProcessInfo
187  ) override;
188 
194  void CalculateDampingMatrix(
195  MatrixType& rDampingMatrix,
196  const ProcessInfo& rCurrentProcessInfo
197  ) override;
198 
206  const Variable<double>& rVariable,
207  std::vector<double>& rValues,
208  const ProcessInfo& rCurrentProcessInfo
209  ) override;
210 
216  void EquationIdVector(
217  EquationIdVectorType& rResult,
218  const ProcessInfo& rCurrentProcessInfo
219  ) const override;
220 
226  void GetDofList(
227  DofsVectorType& rElementalDofList,
228  const ProcessInfo& rCurrentProcessInfo
229  ) const override;
230 
234 
235  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
236 
237  void GetValuesVector(
238  Vector& rValues,
239  int Step = 0) const override;
240 
241  void GetFirstDerivativesVector(
242  Vector& rValues,
243  int Step = 0) const override;
244 
245  void GetSecondDerivativesVector(
246  Vector& rValues,
247  int Step = 0) const override;
248 
252 
260  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
261 
265 
267  std::string Info() const override
268  {
269  std::stringstream buffer;
270  buffer << "TrussEmbeddedEdgeElement #" << Id();
271  return buffer.str();
272  }
273 
275  void PrintInfo(std::ostream& rOStream) const override
276  {
277  rOStream << "TrussEmbeddedEdgeElement #" << Id();
278  }
279 
281  void PrintData(std::ostream& rOStream) const override
282  {
283  pGetGeometry()->PrintData(rOStream);
284  }
285 
287 
288 private:
291 
292  std::vector<array_1d<double, 3>> mReferenceBaseVector;
293 
294  array_1d<double, 3> GetActualBaseVector(const Matrix& r_DN_De, const ConfigurationType& rConfiguration);
295 
299 
301  void CalculateAll(
302  MatrixType& rLeftHandSideMatrix,
303  VectorType& rRightHandSideVector,
304  const ProcessInfo& rCurrentProcessInfo,
305  const bool CalculateStiffnessMatrixFlag,
306  const bool CalculateResidualVectorFlag
307  );
308 
312 
313  friend class Serializer;
314 
315  void save(Serializer& rSerializer) const override
316  {
318  }
319 
320  void load(Serializer& rSerializer) override
321  {
323  }
324 
326 
327 }; // Class TrussEmbeddedEdgeElement
329 
330 } // namespace Kratos.
331 
332 #endif // KRATOS_TRUSS_EMBEDDED_EDGE_ELEMENT_H_INCLUDED defined
Base class for all Elements.
Definition: element.h:60
std::size_t SizeType
Definition: element.h:94
std::size_t IndexType
Definition: flags.h:74
Geometry base class.
Definition: geometry.h:71
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
Definition: truss_embedded_edge_element.h:26
std::size_t IndexType
Definition: truss_embedded_edge_element.h:43
TrussEmbeddedEdgeElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using an array of nodes.
Definition: truss_embedded_edge_element.h:53
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
This function provides a more general interface to the element.
Definition: truss_embedded_edge_element.h:159
ConfigurationType
Internal flags used for calculate reference or current kinematic.
Definition: truss_embedded_edge_element.h:29
Geometry< Node > GeometryType
Definition: truss_embedded_edge_element.h:46
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: truss_embedded_edge_element.h:275
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: truss_embedded_edge_element.h:281
TrussEmbeddedEdgeElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constructor using an array of nodes with properties.
Definition: truss_embedded_edge_element.h:60
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const override
Create with Id, pointer to geometry and pointer to property.
Definition: truss_embedded_edge_element.h:91
virtual ~TrussEmbeddedEdgeElement()=default
Destructor.
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TrussEmbeddedEdgeElement)
Counted pointer of TrussEmbeddedEdgeElement.
std::string Info() const override
Turn back information as a string.
Definition: truss_embedded_edge_element.h:267
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
This is called during the assembling process in order to calculate the condition right hand side matr...
Definition: truss_embedded_edge_element.h:111
TrussEmbeddedEdgeElement()
Default constructor necessary for serialization.
Definition: truss_embedded_edge_element.h:68
Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override
Create with Id, pointer to geometry and pointer to property.
Definition: truss_embedded_edge_element.h:80
std::size_t SizeType
Size types.
Definition: truss_embedded_edge_element.h:42
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
This is called during the assembling process in order to calculate the condition left hand side matri...
Definition: truss_embedded_edge_element.h:134
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
Matrix MatrixType
Definition: geometrical_transformation_utilities.h:55
Modeler::Pointer Create(const std::string &ModelerName, Model &rModel, const Parameters ModelParameters)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:30
pybind11::list CalculateOnIntegrationPoints(TObject &dummy, const Variable< TDataType > &rVariable, const ProcessInfo &rProcessInfo)
Definition: add_mesh_to_python.cpp:142
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
@ Current
Definition: structural_mechanics_math_utilities.hpp:30
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
def load(f)
Definition: ode_solve.py:307