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.
isotropic_shell_element.hpp
Go to the documentation of this file.
1 // KRATOS ___| | | |
2 // \___ \ __| __| | | __| __| | | __| _` | |
3 // | | | | | ( | | | | ( | |
4 // _____/ \__|_| \__,_|\___|\__|\__,_|_| \__,_|_| MECHANICS
5 //
6 // License: BSD License
7 // license: StructuralMechanicsApplication/license.txt
8 //
9 // Main authors: Riccardo Rossi
10 // Massimo Petracca
11 //
12 
13 #pragma once
14 
15 // System includes
16 
17 // External includes
18 
19 // Project includes
20 #include "includes/element.h"
21 
22 namespace Kratos
23 {
24 
27 
31 
35 
39 
43 
45 
47 class KRATOS_API(STRUCTURAL_MECHANICS_APPLICATION) IsotropicShellElement
48  : public Element
49 {
50 public:
53 
56 
60 
62  IsotropicShellElement(IndexType NewId, GeometryType::Pointer pGeometry);
63  IsotropicShellElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties);
64 
66  ~IsotropicShellElement() override;
67 
68 
72 
73 
77 
85  Element::Pointer Create(
86  IndexType NewId,
87  GeometryType::Pointer pGeom,
88  PropertiesType::Pointer pProperties
89  ) const override;
90 
98  Element::Pointer Create(
99  IndexType NewId,
100  NodesArrayType const& ThisNodes,
101  PropertiesType::Pointer pProperties
102  ) const override;
103 
104  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix, VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
105 
106  void CalculateRightHandSide(VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
107 
108  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
109 
110  void GetDofList(DofsVectorType& ElementalDofList, const ProcessInfo& CurrentProcessInfo) const override;
111 
112  void InitializeSolutionStep(const ProcessInfo& CurrentProcessInfo) override;
113 
114  void GetValuesVector(Vector& values, int Step) const override;
115  void GetFirstDerivativesVector(Vector& values, int Step = 0) const override;
116  void GetSecondDerivativesVector(Vector& values, int Step = 0) const override;
117 
118  void CalculateOnIntegrationPoints(const Variable<Matrix >& rVariable, std::vector< Matrix >& Output, const ProcessInfo& rCurrentProcessInfo) override;
119 
120  void CalculateOnIntegrationPoints(const Variable<double >& rVariable, std::vector<double>& Output, const ProcessInfo& rCurrentProcessInfo) override;
121 
122  void Calculate(const Variable<Matrix >& rVariable, Matrix& Output, const ProcessInfo& rCurrentProcessInfo) override;
123 
124  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
125 
126  void FinalizeNonLinearIteration(const ProcessInfo& CurrentProcessInfo) override;
127 
128  void CalculateMassMatrix(MatrixType& rMassMatrix, const ProcessInfo& rCurrentProcessInfo) override;
129 
130  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
131 
135 
136 
140 
141 
145 
147  // virtual String Info() const;
148 
150  // virtual void PrintInfo(std::ostream& rOStream) const;
151 
153  // virtual void PrintData(std::ostream& rOStream) const;
154 
155 
159 
160 
162 
163 protected:
166 
167 
171 
172 
176 
177 
181 
182 
186 
187 
191 
192 
196 
197 
199 
200 private:
203 
209 
210  array_1d< array_1d<double,3>, 3> rot_oldit;
211 
212  double mOrientationAngle;
213 
214 
218  void CalculateLocalGlobalTransformation(
219  double& x12,
220  double& x23,
221  double& x31,
222  double& y12,
223  double& y23,
224  double& y31,
225  array_1d<double,3>& v1,
226  array_1d<double,3>& v2,
227  array_1d<double,3>& v3,
228  double& area
229  );
230 
231  void CalculateMembraneB(
233  const double& beta0,
234  const double& loc1,
235  const double& loc2,
236  const double& loc3,
237  const double& x12,
238  const double& x23,
239  const double& x31,
240  const double& y12,
241  const double& y23,
242  const double& y31
243  );
244 
245 
246  void CalculateBendingB(
248  const double& loc2,
249  const double& loc3,
250  const double& x12,
251  const double& x23,
252  const double& x31,
253  const double& y12,
254  const double& y23,
255  const double& y31
256  );
257 
258  void CalculateMembraneContribution(
259  const BoundedMatrix<double,9,3>& Bm,
260  const BoundedMatrix<double,3,3>& Em,
262  );
263 
264 
265  void AssembleMembraneContribution(
266  const BoundedMatrix<double,9,9>& Km,
267  const double& coeff,
268  BoundedMatrix<double,18,18>& Kloc_system
269  );
270 
271  void CalculateBendingContribution(
272  const BoundedMatrix<double,9,3>& Bb,
273  const BoundedMatrix<double,3,3>& Eb,
275  );
276 
277  void AssembleBendingContribution(
278  const BoundedMatrix<double,9,9>& Kb,
279  const double& coeff,
280  BoundedMatrix<double,18,18>& Kloc_system
281  );
282 
283  void CalculateGaussPointContribution(
284  BoundedMatrix<double,18,18>& Kloc_system ,
285  const BoundedMatrix<double,3,3>& Em,
286  const BoundedMatrix<double,3,3>& Eb,
287  const double& weight,
288  const double& h, /*thickness*/
289  const double& loc1, /*local coords*/
290  const double& loc2,
291  const double& loc3,
292  const double& x12,
293  const double& x23,
294  const double& x31,
295  const double& y12,
296  const double& y23,
297  const double& y31
298  );
299 
300  double CalculateBeta(
301  const BoundedMatrix<double,3,3>& Em
302  );
303 
304  void CalculateMembraneElasticityTensor(
306  const double& h
307  );
308 
309  void CalculateBendingElasticityTensor(
311  const double& h );
312 
313  void CalculateAllMatrices(
314  MatrixType& rLeftHandSideMatrix,
315  VectorType& rRightHandSideVector,
316  const ProcessInfo& rCurrentProcessInfo
317  );
318 
319  Vector& CalculateVolumeForce(
320  Vector& rVolumeForce
321  );
322 
323  void AddBodyForce(
324  const double& h,
325  const double& Area,
326  const Vector& VolumeForce,
327  VectorType& rRightHandSideVector
328  );
329 
330  void RotateToGlobal(
331  const array_1d<double,3>& v1,
332  const array_1d<double,3>& v2,
333  const array_1d<double,3>& v3,
334  const BoundedMatrix<double,18,18>& Kloc_system,
335  Matrix& rLeftHandSideMatrix
336  );
337 
338  void RotateToGlobal(
339  const array_1d<double,3>& v1,
340  const array_1d<double,3>& v2,
341  const array_1d<double,3>& v3,
342  const BoundedMatrix<double,18,18>& Kloc_system,
343  Matrix& rLeftHandSideMatrix,
344  VectorType& rRightHandSideVector
345  );
346 
347 
348 
349  void NicePrint(const Matrix& A);
350 
351  void AddVoigtTensorComponents(
352  const double local_component,
354  const array_1d<double,3>& a,
355  const array_1d<double,3>& b
356  );
357 
358  void CalculateAndAddKg(
359  MatrixType& LHS,
360  BoundedMatrix<double,18,18>& rWorkMatrix,
361  const double& x12,
362  const double& x23,
363  const double& x31,
364  const double& y12,
365  const double& y23,
366  const double& y31,
367  const array_1d<double,3>& v1,
368  const array_1d<double,3>& v2,
369  const array_1d<double,3>& v3,
370  const double& A
371  );
372 
373  void CalculateKg_GaussPointContribution(
374  BoundedMatrix<double,18,18>& Kloc_system ,
375  const BoundedMatrix<double,3,3>& Em,
376  const double& weight,
377  const double& h, /*thickness*/
378  const double& loc1, /*local coords*/
379  const double& loc2,
380  const double& loc3,
381  const double& x12,
382  const double& x23,
383  const double& x31,
384  const double& y12,
385  const double& y23,
386  const double& y31,
387  const array_1d<double,9>& membrane_disp
388  );
389 
390  void CalculateLocalShapeDerivatives(
391  double alpha,
392  BoundedMatrix<double,2,9>& DNu_loc ,
393  BoundedMatrix<double,2,9>& DNv_loc ,
394  BoundedMatrix<double,2,9>& DNw_loc ,
395  const double& a, /*local coords*/ //loc1
396  const double& b, //loc2
397  const double& c, //loc3
398  const double& x12,
399  const double& x23,
400  const double& x31,
401  const double& y12,
402  const double& y23,
403  const double& y31
404  );
405 
406  void CalculateProjectionOperator(
407  BoundedMatrix<double,18,18>& rProjOperator,
408  const double& x12,
409  const double& x23,
410  const double& x31,
411  const double& y12,
412  const double& y23,
413  const double& y31
414  );
415 
416  void ApplyProjection(
417  BoundedMatrix<double,18,18>& rLeftHandSideMatrix,
418  VectorType& rRightHandSideVector,
419  BoundedMatrix<double,18,18>& rWorkMatrix,
420  array_1d<double,18>& rWorkArray,
421  const BoundedMatrix<double,18,18>& rProjOperator
422  );
423 
424  void UpdateNodalReferenceSystem(
425  const double& x12,
426  const double& x23,
427  const double& x31,
428  const double& y12,
429  const double& y23,
430  const double& y31
431  );
432 
433  void SaveOriginalReference(
434  const array_1d<double,3>& v1,
435  const array_1d<double,3>& v2,
436  const array_1d<double,3>& v3
437  );
438 
439  void CalculatePureDisplacement(
440  Vector& values,
441  const array_1d<double,3>& v1,
442  const array_1d<double,3>& v2,
443  const array_1d<double,3>& v3
444  );
445 
446  void CalculatePureMembraneDisplacement(
448  const array_1d<double,3>& v1,
449  const array_1d<double,3>& v2,
450  const array_1d<double,3>& v3
451  );
452 
453  void CalculatePureBendingDisplacement(
455  const array_1d<double,3>& v1,
456  const array_1d<double,3>& v2,
457  const array_1d<double,3>& v3
458  );
459 
460  void InvertMatrix(
461  const BoundedMatrix<double,3,3>& InputMatrix,
462  BoundedMatrix<double,3,3>& InvertedMatrix,
463  double& InputMatrixDet
464  );
465 
466  void SetupOrientationAngles();
467 
468 
472 
473 
477 
478 
482  friend class Serializer;
483 
484  // A private default constructor necessary for serialization
486 
487  void save(Serializer& rSerializer) const override
488  {
490  }
491 
492  void load(Serializer& rSerializer) override
493  {
494  KRATOS_SERIALIZE_LOAD_BASE_CLASS( rSerializer, Element )
495  }
496 
500 
501 
505 
507  //IsotropicShellElement& operator=(const IsotropicShellElement& rOther);
508 
510  //IsotropicShellElement(const IsotropicShellElement& rOther);
511 
512 
514 
515 }; // Class IsotropicShellElement
516 
518 /* inline std::istream& operator >> (std::istream& rIStream,
519  IsotropicShellElement& rThis);
520 */
522 /* inline std::ostream& operator << (std::ostream& rOStream,
523  const IsotropicShellElement& rThis)
524 {
525 rThis.PrintInfo(rOStream);
526 rOStream << std::endl;
527 rThis.PrintData(rOStream);
528 
529 return rOStream;
530 }*/
532 
533 } // namespace Kratos.
534 
535 
Base class for all Elements.
Definition: element.h:60
std::vector< DofType::Pointer > DofsVectorType
Definition: element.h:100
std::vector< std::size_t > EquationIdVectorType
Definition: element.h:98
std::size_t IndexType
Definition: flags.h:74
Short class definition.
Definition: isotropic_shell_element.hpp:49
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(IsotropicShellElement)
Counted pointer of IsotropicShellElement.
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
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
Short class definition.
Definition: array_1d.h:61
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
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
void InitializeSolutionStep(ConstructionUtility &rThisUtil, std::string ThermalSubModelPartName, std::string MechanicalSubModelPartName, std::string HeatFluxSubModelPartName, std::string HydraulicPressureSubModelPartName, bool thermal_conditions, bool mechanical_conditions, int phase)
Definition: add_custom_utilities_to_python.cpp:45
TDataType Calculate(GeometryType &dummy, const Variable< TDataType > &rVariable)
Definition: add_geometries_to_python.cpp:103
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
list coeff
Definition: bombardelli_test.py:41
list values
Definition: bombardelli_test.py:42
v
Definition: generate_convection_diffusion_explicit_element.py:114
alpha
Definition: generate_convection_diffusion_explicit_element.py:113
h
Definition: generate_droplet_dynamics.py:91
a
Definition: generate_stokes_twofluid_element.py:77
b
Definition: generate_total_lagrangian_mixed_volumetric_strain_element.py:31
c
Definition: generate_weakly_compressible_navier_stokes_element.py:108
def load(f)
Definition: ode_solve.py:307
A
Definition: sensitivityMatrix.py:70
B
Definition: sensitivityMatrix.py:76