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.
line_gauss_lobatto_3d_2.h
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: Josep Maria Carbonell
11 //
12 // contributors: Hoang Giang Bui
13 // Riccardo Rossi
14 // Janosch Stascheit
15 // Felix Nagel
16 
17 #if !defined(KRATOS_LINE_GAUSS_LOBATTO_3D_2_H_INCLUDED )
18 #define KRATOS_LINE_GAUSS_LOBATTO_3D_2_H_INCLUDED
19 
20 // System includes
21 
22 // External includes
23 
24 // Project includes
25 #include "geometries/geometry.h"
27 
28 
29 namespace Kratos
30 {
31 
34 
38 
42 
46 
50 
53 template<class TPointType>
54 
55 class LineGaussLobatto3D2 : public Geometry<TPointType>
56 {
57 
58 public:
62 
65 
68 
72 
77 
80  typedef TPointType PointType;
81 
86  typedef typename BaseType::IndexType IndexType;
87 
88 
93  typedef typename BaseType::SizeType SizeType;
94 
99 
105 
112 
118 
123 
128 
134 
140 
144 
149 
153 
154 // LineGaussLobatto3D2( const PointType& FirstPoint, const PointType& SecondPoint )
155 // : BaseType( PointsArrayType(), &msGeometryData )
156 // {
157 // BaseType::Points().push_back( typename PointType::Pointer( new PointType( FirstPoint ) ) );
158 // BaseType::Points().push_back( typename PointType::Pointer( new PointType( SecondPoint ) ) );
159 // }
160 
161  LineGaussLobatto3D2( typename PointType::Pointer pFirstPoint, typename PointType::Pointer pSecondPoint )
162  : BaseType( PointsArrayType(), &msGeometryData )
163  {
164  BaseType::Points().push_back( pFirstPoint );
165  BaseType::Points().push_back( pSecondPoint );
166  }
167 
169  : BaseType( ThisPoints, &msGeometryData )
170  {
171  if ( BaseType::PointsNumber() != 2 )
172  KRATOS_ERROR << "Invalid points number. Expected 2, given " << BaseType::PointsNumber() << std::endl;
173  }
174 
184  : BaseType( rOther )
185  {
186  }
187 
188 
200  template<class TOtherPointType> LineGaussLobatto3D2( LineGaussLobatto3D2<TOtherPointType> const& rOther )
201  : BaseType( rOther )
202  {
203  }
204 
206  ~LineGaussLobatto3D2() override {}
207 
209  {
211  }
212 
214  {
216  }
217 
221 
233  {
234  BaseType::operator=( rOther );
235 
236  return *this;
237  }
238 
249  template<class TOtherPointType>
251  {
252  BaseType::operator=( rOther );
253 
254  return *this;
255  }
256 
260 
261  typename BaseType::Pointer Create( PointsArrayType const& ThisPoints ) const override
262  {
263  return typename BaseType::Pointer( new LineGaussLobatto3D2( ThisPoints ) );
264  }
265 
275  Vector& rResult,
276  const typename BaseType::LumpingMethods LumpingMethod = BaseType::LumpingMethods::ROW_SUM
277  ) const override
278  {
279  if(rResult.size() != 2)
280  rResult.resize( 2, false );
281  rResult[0] = 0.5;
282  rResult[1] = 0.5;
283  return rResult;
284  }
285 
289 
302  double Length() const override
303  {
304  const TPointType& point0 = BaseType::GetPoint(0);
305  const TPointType& point1 = BaseType::GetPoint(1);
306  const double lx = point0.X() - point1.X();
307  const double ly = point0.Y() - point1.Y();
308  const double lz = point0.Z() - point1.Z();
309 
310  const double length = lx * lx + ly * ly + lz * lz;
311 
312  return sqrt( length );
313  }
314 
326  double Area() const override
327  {
328  return Length();
329  }
330 
331 
342  double DomainSize() const override
343  {
344  const TPointType& point0 = BaseType::GetPoint(0);
345  const TPointType& point1 = BaseType::GetPoint(1);
346  const double lx = point0.X() - point1.X();
347  const double ly = point0.Y() - point1.Y();
348  const double lz = point0.Z() - point1.Z();
349 
350  const double length = lx * lx + ly * ly + lz * lz;
351 
352  return sqrt( length );
353  }
354 
355 
359 
360 
375  JacobiansType& Jacobian( JacobiansType& rResult, IntegrationMethod ThisMethod ) const override
376  {
377  Matrix jacobian( 3, 1 );
378  jacobian( 0, 0 ) = ( this->GetPoint( 1 ).X() - this->GetPoint( 0 ).X() ) * 0.5; //on the Gauss points (J is constant at each element)
379  jacobian( 1, 0 ) = ( this->GetPoint( 1 ).Y() - this->GetPoint( 0 ).Y() ) * 0.5;
380  jacobian( 2, 0 ) = ( this->GetPoint( 1 ).Z() - this->GetPoint( 0 ).Z() ) * 0.5;
381 
382  if ( rResult.size() != BaseType::IntegrationPointsNumber( ThisMethod ) )
383  {
384  // KLUDGE: While there is a bug in ublas vector resize, I have to put this beside resizing!!
386  rResult.swap( temp );
387  }
388 
389  std::fill( rResult.begin(), rResult.end(), jacobian );
390 
391  return rResult;
392  }
393 
411  JacobiansType& Jacobian( JacobiansType& rResult, IntegrationMethod ThisMethod, Matrix & DeltaPosition ) const override
412  {
413  Matrix jacobian( 3, 1 );
414  jacobian( 0, 0 ) = ( (this->GetPoint( 1 ).X() - DeltaPosition(1,0)) - (this->GetPoint( 0 ).X() - DeltaPosition(0,0)) ) * 0.5; //on the Gauss points (J is constant at each element)
415  jacobian( 1, 0 ) = ( (this->GetPoint( 1 ).Y() - DeltaPosition(1,1)) - (this->GetPoint( 0 ).Y() - DeltaPosition(0,1)) ) * 0.5;
416  jacobian( 2, 0 ) = ( (this->GetPoint( 1 ).Z() - DeltaPosition(1,2)) - (this->GetPoint( 0 ).Z() - DeltaPosition(0,2)) ) * 0.5;
417 
418  if ( rResult.size() != BaseType::IntegrationPointsNumber( ThisMethod ) )
419  {
420  // KLUDGE: While there is a bug in ublas vector resize, I have to put this beside resizing!!
422  rResult.swap( temp );
423  }
424 
425  std::fill( rResult.begin(), rResult.end(), jacobian );
426 
427  return rResult;
428  }
429 
447  Matrix& Jacobian( Matrix& rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod ) const override
448  {
449  rResult.resize( 3, 1, false );
450  //on the Gauss points (J is constant at each element)
451  rResult( 0, 0 ) = ( this->GetPoint( 1 ).X() - this->GetPoint( 0 ).X() ) * 0.5;
452  rResult( 1, 0 ) = ( this->GetPoint( 1 ).Y() - this->GetPoint( 0 ).Y() ) * 0.5;
453  rResult( 2, 0 ) = ( this->GetPoint( 1 ).Z() - this->GetPoint( 0 ).Z() ) * 0.5;
454  return rResult;
455  }
456 
468  Matrix& Jacobian( Matrix& rResult, const CoordinatesArrayType& rPoint ) const override
469  {
470  rResult.resize( 3, 1, false );
471  //on the Gauss points (J is constant at each element)
472  rResult( 0, 0 ) = ( this->GetPoint( 1 ).X() - this->GetPoint( 0 ).X() ) * 0.5;
473  rResult( 1, 0 ) = ( this->GetPoint( 1 ).Y() - this->GetPoint( 0 ).Y() ) * 0.5;
474  rResult( 2, 0 ) = ( this->GetPoint( 1 ).Z() - this->GetPoint( 0 ).Z() ) * 0.5;
475  return rResult;
476  }
477 
489  Vector& DeterminantOfJacobian( Vector& rResult, IntegrationMethod ThisMethod ) const override
490  {
491  rResult = ZeroVector( 1 );
492  rResult[0] = 0.5 * MathUtils<double>::Norm3(( this->GetPoint( 1 ) ) - ( this->GetPoint( 0 ) ) );
493  return rResult;
494  }
495 
514  double DeterminantOfJacobian( IndexType IntegrationPointIndex, IntegrationMethod ThisMethod ) const override
515  {
516  return( 0.5*MathUtils<double>::Norm3(( this->GetPoint( 1 ) ) - ( this->GetPoint( 0 ) ) ) );
517  }
518 
531  double DeterminantOfJacobian( const CoordinatesArrayType& rPoint ) const override
532  {
533  return( 0.5*MathUtils<double>::Norm3(( this->GetPoint( 1 ) ) - ( this->GetPoint( 0 ) ) ) );
534  }
535 
536 
551  JacobiansType& InverseOfJacobian( JacobiansType& rResult, IntegrationMethod ThisMethod ) const override
552  {
553  rResult[0] = ZeroMatrix( 1, 1 );
554  rResult[0]( 0, 0 ) = 2.0 * MathUtils<double>::Norm3(( this->GetPoint( 1 ) ) - ( this->GetPoint( 0 ) ) );
555  return rResult;
556  }
557 
575  Matrix& InverseOfJacobian( Matrix& rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod ) const override
576  {
577  rResult = ZeroMatrix( 1, 1 );
578  rResult( 0, 0 ) = 2.0 * MathUtils<double>::Norm3(( this->GetPoint( 1 ) ) - ( this->GetPoint( 0 ) ) );
579  return( rResult );
580  }
581 
593  Matrix& InverseOfJacobian( Matrix& rResult, const CoordinatesArrayType& rPoint ) const override
594  {
595  rResult = ZeroMatrix( 1, 1 );
596  rResult( 0, 0 ) = 2.0 * MathUtils<double>::Norm3(( this->GetPoint( 1 ) ) - ( this->GetPoint( 0 ) ) );
597  return( rResult );
598  }
599 
600 
604  SizeType EdgesNumber() const override
605  {
606  return 2;
607  }
608 
609 
610  SizeType FacesNumber() const override
611  {
612  return 2;
613  }
614 
615 
619 
620  double ShapeFunctionValue( IndexType ShapeFunctionIndex,
621  const CoordinatesArrayType& rPoint ) const override
622  {
623  switch ( ShapeFunctionIndex )
624  {
625 
626  case 0:
627  return( 0.5*( 1.0 - rPoint[0] ) );
628 
629  case 1:
630  return( 0.5*( 1.0 + rPoint[0] ) );
631 
632  default:
633  KRATOS_ERROR << "Wrong index of shape function!" << *this << std::endl;
634  }
635 
636  return 0;
637  }
638 
640  const CoordinatesArrayType& rPoint ) const override
641  {
642  rResult = ZeroMatrix( 2, 1 );
643  rResult( 0, 0 ) = -0.5;
644  rResult( 1, 0 ) = 0.5;
645  return( rResult );
646  }
647 
654  std::string Info() const override
655  {
656  return "1 dimensional line with 2 nodes in 3D space";
657  }
658 
665  void PrintInfo( std::ostream& rOStream ) const override
666  {
667  rOStream << "1 dimensional line with 2 nodes in 3D space";
668  }
669 
678  void PrintData( std::ostream& rOStream ) const override
679  {
680  // Base Geometry class PrintData call
681  BaseType::PrintData( rOStream );
682  std::cout << std::endl;
683 
684  // If the geometry has valid points, calculate and output its data
685  if (this->AllPointsAreValid()) {
686  Matrix jacobian;
687  this->Jacobian( jacobian, PointType() );
688  rOStream << " Jacobian\t : " << jacobian;
689  }
690  }
691 
695 
696 
698 
699 protected:
702 
703 
707 
708 
712 
713 
717 
718 
722 
723 
727 
728 
732 
733 
735 
736 private:
739 
740  static const GeometryData msGeometryData;
741 
742  static const GeometryDimension msGeometryDimension;
743 
747 
751 
752  friend class Serializer;
753 
754  void save( Serializer& rSerializer ) const override
755  {
757  }
758 
759  void load( Serializer& rSerializer ) override
760  {
762  }
763 
764  LineGaussLobatto3D2(): BaseType( PointsArrayType(), &msGeometryData ) {}
765 
766 
770 
771 
775 
776  static Matrix CalculateShapeFunctionsIntegrationPointsValues( typename BaseType::IntegrationMethod ThisMethod )
777  {
778  const IntegrationPointsContainerType& all_integration_points = AllIntegrationPoints();
779  const IntegrationPointsArrayType& IntegrationPoints = all_integration_points[static_cast<int>(ThisMethod)];
780  int integration_points_number = IntegrationPoints.size();
781  Matrix N( integration_points_number, 2 );
782 
783  for ( int it_gp = 0; it_gp < integration_points_number; it_gp++ )
784  {
785  double e = IntegrationPoints[it_gp].X();
786  N( it_gp, 0 ) = 0.5 * ( 1 - e );
787  N( it_gp, 1 ) = 0.5 * ( 1 + e );
788  }
789 
790  return N;
791  }
792 
793  static ShapeFunctionsGradientsType CalculateShapeFunctionsIntegrationPointsLocalGradients( typename BaseType::IntegrationMethod ThisMethod )
794  {
795  const IntegrationPointsContainerType& all_integration_points = AllIntegrationPoints();
796  const IntegrationPointsArrayType& IntegrationPoints = all_integration_points[static_cast<int>(ThisMethod)];
798 
799  for ( unsigned int it_gp = 0; it_gp < IntegrationPoints.size(); it_gp++ )
800  {
801  Matrix aux_mat = ZeroMatrix(2,1);
802  aux_mat(0,0) = -0.5;
803  aux_mat(1,0) = 0.5;
804  DN_De[it_gp] = aux_mat;
805  }
806 
807  return DN_De;
808 
809  }
810 
811  static const IntegrationPointsContainerType AllIntegrationPoints()
812  {
813  IntegrationPointsContainerType integration_points = {{
814  Quadrature<LineGaussLobattoIntegrationPoints1, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
815  Quadrature<LineGaussLobattoIntegrationPoints2, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
816  Quadrature<LineGaussLobattoIntegrationPoints3, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
817  Quadrature<LineGaussLobattoIntegrationPoints4, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
818  Quadrature<LineGaussLobattoIntegrationPoints5, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
819  Quadrature<LineGaussLobattoIntegrationPoints6, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
820  Quadrature<LineGaussLobattoIntegrationPoints7, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
821  Quadrature<LineGaussLobattoIntegrationPoints8, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
822  Quadrature<LineGaussLobattoIntegrationPoints9, 1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
823  Quadrature<LineGaussLobattoIntegrationPoints10,1, IntegrationPoint<3> >::GenerateIntegrationPoints(),
824  }
825  };
826  return integration_points;
827  }
828 
829  static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
830  {
831  ShapeFunctionsValuesContainerType shape_functions_values = {{
832  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_1 ),
833  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_2 ),
834  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_3 ),
835  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_4 ),
836  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_5 ),
837  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_1 ),
838  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_2 ),
839  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_3 ),
840  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_4 ),
841  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_5 )
842  }
843  };
844  return shape_functions_values;
845  }
846 
847  static const ShapeFunctionsLocalGradientsContainerType AllShapeFunctionsLocalGradients()
848  {
849  ShapeFunctionsLocalGradientsContainerType shape_functions_local_gradients = {{
850  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_1 ),
851  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_2 ),
852  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_3 ),
853  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_4 ),
854  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_5 ),
855  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_1 ),
856  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_2 ),
857  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_3 ),
858  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_4 ),
859  LineGaussLobatto3D2<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_5 )
860 
861  }
862  };
863  return shape_functions_local_gradients;
864  }
865 
869 
870 
874 
875 
879 
880  template<class TOtherPointType> friend class LineGaussLobatto3D2;
881 
885 
887 
888 }; // Class Geometry
889 
891 
894 
895 
899 
900 
902 template<class TPointType>
903 inline std::istream& operator >> ( std::istream& rIStream,
905 
907 template<class TPointType>
908 inline std::ostream& operator << ( std::ostream& rOStream,
909  const LineGaussLobatto3D2<TPointType>& rThis )
910 {
911  rThis.PrintInfo( rOStream );
912  rOStream << std::endl;
913  rThis.PrintData( rOStream );
914 
915  return rOStream;
916 }
917 
919 
920 
921 template<class TPointType>
922 const GeometryData LineGaussLobatto3D2<TPointType>::msGeometryData(
925  LineGaussLobatto3D2<TPointType>::AllIntegrationPoints(),
926  LineGaussLobatto3D2<TPointType>::AllShapeFunctionsValues(),
927  AllShapeFunctionsLocalGradients() );
928 
929 template<class TPointType>
930 const GeometryDimension LineGaussLobatto3D2<TPointType>::msGeometryDimension(3, 1);
931 
932 } // namespace Kratos.
933 
934 #endif // KRATOS_LINE_GAUSS_LOBATTO_3D_2_H_INCLUDED defined
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
IntegrationMethod
Definition: geometry_data.h:76
KratosGeometryFamily
Definition: geometry_data.h:91
Definition: geometry_dimension.h:42
Geometry base class.
Definition: geometry.h:71
SizeType PointsNumber() const
Definition: geometry.h:528
virtual void PrintData(std::ostream &rOStream) const
Print object's data.
Definition: geometry.h:3834
Geometry & operator=(const Geometry &rOther)
Definition: geometry.h:400
std::vector< IntegrationPointType > IntegrationPointsArrayType
Definition: geometry.h:161
std::size_t SizeType
Definition: geometry.h:144
std::size_t IndexType
Definition: geometry.h:137
std::array< Matrix, static_cast< int >GeometryData::IntegrationMethod::NumberOfIntegrationMethods)> ShapeFunctionsValuesContainerType
Definition: geometry.h:172
const PointsArrayType & Points() const
Definition: geometry.h:1768
bool AllPointsAreValid() const
Checks if the geometry points are valid Checks if the geometry points are valid from the pointer valu...
Definition: geometry.h:4093
std::array< IntegrationPointsArrayType, static_cast< int >GeometryData::IntegrationMethod::NumberOfIntegrationMethods)> IntegrationPointsContainerType
Definition: geometry.h:167
LumpingMethods
This defines the different methods to compute the lumping methods.
Definition: geometry.h:109
const IntegrationPointsArrayType & IntegrationPoints() const
Definition: geometry.h:2284
TPointType const & GetPoint(const int Index) const
Definition: geometry.h:1816
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
SizeType IntegrationPointsNumber() const
Definition: geometry.h:2257
GeometryData::IntegrationMethod IntegrationMethod
Definition: geometry.h:122
Short class definition.
Definition: integration_point.h:52
Definition: amatrix_interface.h:41
void swap(Matrix &Other)
Definition: amatrix_interface.h:289
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
iterator end()
Definition: amatrix_interface.h:243
iterator begin()
Definition: amatrix_interface.h:241
Definition: line_gauss_lobatto_3d_2.h:56
GeometryData::IntegrationMethod IntegrationMethod
Definition: line_gauss_lobatto_3d_2.h:71
Matrix & Jacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: line_gauss_lobatto_3d_2.h:447
BaseType::IntegrationPointType IntegrationPointType
Definition: line_gauss_lobatto_3d_2.h:104
Matrix & Jacobian(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: line_gauss_lobatto_3d_2.h:468
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: line_gauss_lobatto_3d_2.h:122
Vector & DeterminantOfJacobian(Vector &rResult, IntegrationMethod ThisMethod) const override
Definition: line_gauss_lobatto_3d_2.h:489
double Area() const override
Definition: line_gauss_lobatto_3d_2.h:326
SizeType FacesNumber() const override
Returns the number of faces of the current geometry.
Definition: line_gauss_lobatto_3d_2.h:610
double DeterminantOfJacobian(IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: line_gauss_lobatto_3d_2.h:514
JacobiansType & InverseOfJacobian(JacobiansType &rResult, IntegrationMethod ThisMethod) const override
Definition: line_gauss_lobatto_3d_2.h:551
double DomainSize() const override
Definition: line_gauss_lobatto_3d_2.h:342
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: line_gauss_lobatto_3d_2.h:639
LineGaussLobatto3D2(const PointsArrayType &ThisPoints)
Definition: line_gauss_lobatto_3d_2.h:168
BaseType::IndexType IndexType
Definition: line_gauss_lobatto_3d_2.h:86
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: line_gauss_lobatto_3d_2.h:139
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: line_gauss_lobatto_3d_2.h:117
void PrintData(std::ostream &rOStream) const override
Definition: line_gauss_lobatto_3d_2.h:678
SizeType EdgesNumber() const override
Definition: line_gauss_lobatto_3d_2.h:604
Vector & LumpingFactors(Vector &rResult, const typename BaseType::LumpingMethods LumpingMethod=BaseType::LumpingMethods::ROW_SUM) const override
Lumping factors for the calculation of the lumped mass matrix.
Definition: line_gauss_lobatto_3d_2.h:274
LineGaussLobatto3D2(LineGaussLobatto3D2 const &rOther)
Definition: line_gauss_lobatto_3d_2.h:183
Matrix & InverseOfJacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: line_gauss_lobatto_3d_2.h:575
BaseType::PointsArrayType PointsArrayType
Definition: line_gauss_lobatto_3d_2.h:98
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: line_gauss_lobatto_3d_2.h:213
BaseType::JacobiansType JacobiansType
Definition: line_gauss_lobatto_3d_2.h:133
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: line_gauss_lobatto_3d_2.h:208
LineGaussLobatto3D2 & operator=(const LineGaussLobatto3D2 &rOther)
Definition: line_gauss_lobatto_3d_2.h:232
LineGaussLobatto3D2(typename PointType::Pointer pFirstPoint, typename PointType::Pointer pSecondPoint)
Definition: line_gauss_lobatto_3d_2.h:161
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: line_gauss_lobatto_3d_2.h:111
double DeterminantOfJacobian(const CoordinatesArrayType &rPoint) const override
Definition: line_gauss_lobatto_3d_2.h:531
LineGaussLobatto3D2(LineGaussLobatto3D2< TOtherPointType > const &rOther)
Definition: line_gauss_lobatto_3d_2.h:200
JacobiansType & Jacobian(JacobiansType &rResult, IntegrationMethod ThisMethod, Matrix &DeltaPosition) const override
Definition: line_gauss_lobatto_3d_2.h:411
LineGaussLobatto3D2 & operator=(LineGaussLobatto3D2< TOtherPointType > const &rOther)
Definition: line_gauss_lobatto_3d_2.h:250
void PrintInfo(std::ostream &rOStream) const override
Definition: line_gauss_lobatto_3d_2.h:665
double Length() const override
Definition: line_gauss_lobatto_3d_2.h:302
Geometry< TPointType > BaseType
Geometry as base class.
Definition: line_gauss_lobatto_3d_2.h:64
TPointType PointType
Definition: line_gauss_lobatto_3d_2.h:80
std::string Info() const override
Definition: line_gauss_lobatto_3d_2.h:654
friend class LineGaussLobatto3D2
Definition: line_gauss_lobatto_3d_2.h:880
BaseType::Pointer Create(PointsArrayType const &ThisPoints) const override
Creates a new geometry pointer.
Definition: line_gauss_lobatto_3d_2.h:261
~LineGaussLobatto3D2() override
Destructor. Do nothing!!!
Definition: line_gauss_lobatto_3d_2.h:206
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: line_gauss_lobatto_3d_2.h:148
BaseType::GeometriesArrayType GeometriesArrayType
Definition: line_gauss_lobatto_3d_2.h:76
Matrix & InverseOfJacobian(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: line_gauss_lobatto_3d_2.h:593
JacobiansType & Jacobian(JacobiansType &rResult, IntegrationMethod ThisMethod) const override
Definition: line_gauss_lobatto_3d_2.h:375
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: line_gauss_lobatto_3d_2.h:127
KRATOS_CLASS_POINTER_DEFINITION(LineGaussLobatto3D2)
Pointer definition of LineGaussLobatto3D2.
BaseType::NormalType NormalType
Definition: line_gauss_lobatto_3d_2.h:143
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: line_gauss_lobatto_3d_2.h:620
BaseType::SizeType SizeType
Definition: line_gauss_lobatto_3d_2.h:93
Various mathematical utilitiy functions.
Definition: math_utils.h:62
static double Norm3(const TVectorType &a)
Calculates the norm of vector "a" which is assumed to be of size 3.
Definition: math_utils.h:691
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
void push_back(const TPointerType &x)
Definition: pointer_vector.h:270
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
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
#define KRATOS_ERROR
Definition: exception.h:161
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
const GeometryData LineGaussLobatto3D2< TPointType >::msGeometryData & msGeometryDimension(), LineGaussLobatto3D2< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
Internals::Matrix< double, AMatrix::dynamic, AMatrix::dynamic > Matrix
Definition: amatrix_interface.h:470
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
const GeometryData BrepCurve< TContainerPointType, TContainerPointEmbeddedType >::msGeometryData & msGeometryDimension
Definition: brep_curve.h:483
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
def load(f)
Definition: ode_solve.py:307
float temp
Definition: rotating_cone.py:85
N
Definition: sensitivityMatrix.py:29
e
Definition: run_cpp_mpi_tests.py:31