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.
quadrilateral_2d_4.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: Riccardo Rossi
11 // Janosch Stascheit
12 // Felix Nagel
13 // contributors: Hoang Giang Bui
14 // Josep Maria Carbonell
15 //
16 
17 #pragma once
18 
19 // System includes
20 
21 // External includes
22 
23 // Project includes
24 #include "geometries/line_2d_2.h"
29 
30 namespace Kratos
31 {
34 
38 
42 
46 
50 
71 template<class TPointType> class Quadrilateral2D4
72  : public Geometry<TPointType>
73 {
74 public:
78 
83 
88 
93 
98 
104 
108  typedef TPointType PointType;
109 
116  typedef typename BaseType::IndexType IndexType;
117 
123  typedef typename BaseType::SizeType SizeType;
124 
130 
135 
141 
150 
157 
163 
169 
176 
183 
191 
199 
204 
208 
209 // Quadrilateral2D4( const PointType& FirstPoint,
210 // const PointType& SecondPoint,
211 // const PointType& ThirdPoint,
212 // const PointType& FourthPoint )
213 // : BaseType( PointsArrayType(), &msGeometryData )
214 // {
215 // this->Points().push_back( typename PointType::Pointer( new PointType( FirstPoint ) ) );
216 // this->Points().push_back( typename PointType::Pointer( new PointType( SecondPoint ) ) );
217 // this->Points().push_back( typename PointType::Pointer( new PointType( ThirdPoint ) ) );
218 // this->Points().push_back( typename PointType::Pointer( new PointType( FourthPoint ) ) );
219 // }
220 
221  Quadrilateral2D4( typename PointType::Pointer pFirstPoint,
222  typename PointType::Pointer pSecondPoint,
223  typename PointType::Pointer pThirdPoint,
224  typename PointType::Pointer pFourthPoint )
225  : BaseType( PointsArrayType(), &msGeometryData )
226  {
227  this->Points().push_back( pFirstPoint );
228  this->Points().push_back( pSecondPoint );
229  this->Points().push_back( pThirdPoint );
230  this->Points().push_back( pFourthPoint );
231  }
232 
233  explicit Quadrilateral2D4( const PointsArrayType& rThisPoints )
234  : BaseType( rThisPoints, &msGeometryData )
235  {
236  if ( this->PointsNumber() != 4 )
237  KRATOS_ERROR << "Invalid points number. Expected 4, given " << this->PointsNumber() << std::endl;
238  }
239 
242  const IndexType GeometryId,
243  const PointsArrayType& rThisPoints
244  ) : BaseType(GeometryId, rThisPoints, &msGeometryData)
245  {
246  KRATOS_ERROR_IF( this->PointsNumber() != 4 ) << "Invalid points number. Expected 4, given " << this->PointsNumber() << std::endl;
247  }
248 
251  const std::string& rGeometryName,
252  const PointsArrayType& rThisPoints
253  ) : BaseType(rGeometryName, rThisPoints, &msGeometryData)
254  {
255  KRATOS_ERROR_IF(this->PointsNumber() != 4) << "Invalid points number. Expected 4, given " << this->PointsNumber() << std::endl;
256  }
257 
268  : BaseType( rOther )
269  {
270  }
271 
284  template<class TOtherPointType> explicit Quadrilateral2D4( Quadrilateral2D4<TOtherPointType> const& rOther )
285  : BaseType( rOther )
286  {
287  }
288 
292  ~Quadrilateral2D4() override {}
293 
295  {
297  }
298 
300  {
302  }
303 
307 
320  {
321  BaseType::operator=( rOther );
322  return *this;
323  }
324 
336  template<class TOtherPointType>
338  {
339  BaseType::operator=( rOther );
340  return *this;
341  }
342 
346 
353  typename BaseType::Pointer Create(
354  const IndexType NewGeometryId,
355  PointsArrayType const& rThisPoints
356  ) const override
357  {
358  return typename BaseType::Pointer( new Quadrilateral2D4( NewGeometryId, rThisPoints ) );
359  }
360 
367  typename BaseType::Pointer Create(
368  const IndexType NewGeometryId,
369  const BaseType& rGeometry
370  ) const override
371  {
372  auto p_geometry = typename BaseType::Pointer( new Quadrilateral2D4( NewGeometryId, rGeometry.Points() ) );
373  p_geometry->SetData(rGeometry.GetData());
374  return p_geometry;
375  }
376 
378  SizeType PointsNumberInDirection(IndexType LocalDirectionIndex) const override
379  {
380  if ((LocalDirectionIndex == 0) || (LocalDirectionIndex == 1)) {
381  return 2;
382  }
383  KRATOS_ERROR << "Possible direction index reaches from 0-1. Given direction index: "
384  << LocalDirectionIndex << std::endl;
385  }
386 
392  Matrix& PointsLocalCoordinates( Matrix& rResult ) const override
393  {
394  if (rResult.size1() != 4 || rResult.size2() != 2)
395  rResult.resize(4, 2, false);
396  rResult(0, 0) = -1.0;
397  rResult(0, 1) = -1.0;
398  rResult(1, 0) = 1.0;
399  rResult(1, 1) = -1.0;
400  rResult(2, 0) = 1.0;
401  rResult(2, 1) = 1.0;
402  rResult(3, 0) = -1.0;
403  rResult(3, 1) = 1.0;
404  return rResult;
405  }
406 
410 
430  double Length() const override
431  {
432  //return sqrt(fabs( DeterminantOfJacobian(PointType())));
433  double length = 0.000;
434  length = sqrt( fabs( Area() ) );
435  return length;
436 
437  }
438 
449  double Area() const override
450  {
452  }
453 
461  double Volume() const override
462  {
463  KRATOS_WARNING("Quadrilateral2D4") << "Method not well defined. Replace with DomainSize() instead. This method preserves current behaviour but will be changed in June 2023 (returning error instead)" << std::endl;
464  return Area();
465  // TODO: Replace in June 2023
466  // KRATOS_ERROR << "Quadrilateral2D4:: Method not well defined. Replace with DomainSize() instead." << std::endl;
467  // return 0.0;
468  }
469 
479  double DomainSize() const override
480  {
481  return Area();
482  }
483 
492  bool IsInside(
493  const CoordinatesArrayType& rPoint,
494  CoordinatesArrayType& rResult,
495  const double Tolerance = std::numeric_limits<double>::epsilon()
496  ) const override
497  {
498  this->PointLocalCoordinates( rResult, rPoint );
499 
500  if ( std::abs(rResult[0]) <= (1.0+Tolerance) )
501  {
502  if ( std::abs(rResult[1]) <= (1.0+Tolerance) )
503  {
504  return true;
505  }
506  }
507 
508  return false;
509  }
510 
514 
526  SizeType EdgesNumber() const override
527  {
528  return 4;
529  }
530 
540  {
542  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 0 ), this->pGetPoint( 1 ) ) );
543  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 1 ), this->pGetPoint( 2 ) ) );
544  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 2 ), this->pGetPoint( 3 ) ) );
545  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 3 ), this->pGetPoint( 0 ) ) );
546  return edges;
547  }
548 
557  bool HasIntersection( const Point& rLowPoint, const Point& rHighPoint ) const override
558  {
559  Triangle2D3<PointType> triangle_0 (this->pGetPoint( 0 ),
560  this->pGetPoint( 1 ),
561  this->pGetPoint( 2 )
562  );
563  Triangle2D3<PointType> triangle_1 (this->pGetPoint( 2 ),
564  this->pGetPoint( 3 ),
565  this->pGetPoint( 0 )
566  );
567 
568  if ( triangle_0.HasIntersection(rLowPoint, rHighPoint) ) return true;
569  else if ( triangle_1.HasIntersection(rLowPoint, rHighPoint) ) return true;
570  else return false;
571  }
572 
576 
586  double ShapeFunctionValue( IndexType ShapeFunctionIndex,
587  const CoordinatesArrayType& rPoint ) const override
588  {
589  switch ( ShapeFunctionIndex )
590  {
591  case 0:
592  return( 0.25*( 1.0 - rPoint[0] )*( 1.0 - rPoint[1] ) );
593  case 1:
594  return( 0.25*( 1.0 + rPoint[0] )*( 1.0 - rPoint[1] ) );
595  case 2:
596  return( 0.25*( 1.0 + rPoint[0] )*( 1.0 + rPoint[1] ) );
597  case 3:
598  return( 0.25*( 1.0 - rPoint[0] )*( 1.0 + rPoint[1] ) );
599  default:
600  KRATOS_ERROR << "Wrong index of shape function!" << *this << std::endl;
601  }
602 
603  return 0;
604  }
605 
617  Vector& ShapeFunctionsValues (Vector &rResult, const CoordinatesArrayType& rCoordinates) const override
618  {
619  if(rResult.size() != 4) rResult.resize(4,false);
620  rResult[0] = 0.25*( 1.0 - rCoordinates[0] )*( 1.0 - rCoordinates[1] );
621  rResult[1] = 0.25*( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] );
622  rResult[2] = 0.25*( 1.0 + rCoordinates[0] )*( 1.0 + rCoordinates[1] );
623  rResult[3] = 0.25*( 1.0 - rCoordinates[0] )*( 1.0 + rCoordinates[1] );
624 
625  return rResult;
626  }
627 
628 
629 
633 
641  std::string Info() const override
642  {
643  return "2 dimensional quadrilateral with four nodes in 2D space";
644  }
645 
652  void PrintInfo( std::ostream& rOStream ) const override
653  {
654  rOStream << "2 dimensional quadrilateral with four nodes in 2D space";
655  }
656 
671  void PrintData( std::ostream& rOStream ) const override
672  {
673  // Base Geometry class PrintData call
674  BaseType::PrintData( rOStream );
675  std::cout << std::endl;
676 
677  // If the geometry has valid points, calculate and output its data
678  if (this->AllPointsAreValid()) {
679  Matrix jacobian;
680  this->Jacobian( jacobian, PointType() );
681  rOStream << " Jacobian in the origin\t : " << jacobian;
682  }
683  }
684 
690  IntegrationMethod ThisMethod )
691  {
692  ShapeFunctionsGradientsType localGradients
693  = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
694  const int integration_points_number
695  = msGeometryData.IntegrationPointsNumber( ThisMethod );
696  ShapeFunctionsGradientsType Result( integration_points_number );
697 
698  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
699  {
700  Result[pnt] = localGradients[pnt];
701  }
702 
703  return Result;
704  }
705 
711  {
712  IntegrationMethod ThisMethod = msGeometryData.DefaultIntegrationMethod();
713  ShapeFunctionsGradientsType localGradients
714  = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
715  const int integration_points_number
716  = msGeometryData.IntegrationPointsNumber( ThisMethod );
717  ShapeFunctionsGradientsType Result( integration_points_number );
718 
719  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
720  {
721  Result[pnt] = localGradients[pnt];
722  }
723 
724  return Result;
725  }
726 
737  const CoordinatesArrayType& rPoint ) const override
738  {
739  rResult.resize( 4, 2 , false);
740  noalias( rResult ) = ZeroMatrix( 4, 2 );
741  rResult( 0, 0 ) = -0.25 * ( 1.0 - rPoint[1] );
742  rResult( 0, 1 ) = -0.25 * ( 1.0 - rPoint[0] );
743  rResult( 1, 0 ) = 0.25 * ( 1.0 - rPoint[1] );
744  rResult( 1, 1 ) = -0.25 * ( 1.0 + rPoint[0] );
745  rResult( 2, 0 ) = 0.25 * ( 1.0 + rPoint[1] );
746  rResult( 2, 1 ) = 0.25 * ( 1.0 + rPoint[0] );
747  rResult( 3, 0 ) = -0.25 * ( 1.0 + rPoint[1] );
748  rResult( 3, 1 ) = 0.25 * ( 1.0 - rPoint[0] );
749  return rResult;
750  }
751 
752 
753 
763  virtual Matrix& ShapeFunctionsGradients( Matrix& rResult, PointType& rPoint )
764  {
765  rResult.resize( 4, 2 , false);
766  rResult( 0, 0 ) = -0.25 * ( 1.0 - rPoint.Y() );
767  rResult( 0, 1 ) = -0.25 * ( 1.0 - rPoint.X() );
768  rResult( 1, 0 ) = 0.25 * ( 1.0 - rPoint.Y() );
769  rResult( 1, 1 ) = -0.25 * ( 1.0 + rPoint.X() );
770  rResult( 2, 0 ) = 0.25 * ( 1.0 + rPoint.Y() );
771  rResult( 2, 1 ) = 0.25 * ( 1.0 + rPoint.X() );
772  rResult( 3, 0 ) = -0.25 * ( 1.0 + rPoint.Y() );
773  rResult( 3, 1 ) = 0.25 * ( 1.0 - rPoint.X() );
774  return rResult;
775  }
776 
784  {
785  if ( rResult.size() != this->PointsNumber() )
786  {
787  // KLUDGE: While there is a bug in
788  // ublas vector resize, I have to put this beside resizing!!
790  rResult.swap( temp );
791  }
792 
793  rResult[0].resize( 2, 2 , false);
794  rResult[1].resize( 2, 2 , false);
795  rResult[2].resize( 2, 2 , false);
796  rResult[3].resize( 2, 2 , false);
797 
798  rResult[0]( 0, 0 ) = 0.0;
799  rResult[0]( 0, 1 ) = 0.25;
800  rResult[0]( 1, 0 ) = 0.25;
801  rResult[0]( 1, 1 ) = 0.0;
802 
803  rResult[1]( 0, 0 ) = 0.0;
804  rResult[1]( 0, 1 ) = -0.25;
805  rResult[1]( 1, 0 ) = -0.25;
806  rResult[1]( 1, 1 ) = 0.0;
807 
808  rResult[2]( 0, 0 ) = 0.0;
809  rResult[2]( 0, 1 ) = 0.25;
810  rResult[2]( 1, 0 ) = 0.25;
811  rResult[2]( 1, 1 ) = 0.0;
812 
813  rResult[3]( 0, 0 ) = 0.0;
814  rResult[3]( 0, 1 ) = -0.25;
815  rResult[3]( 1, 0 ) = -0.25;
816  rResult[3]( 1, 1 ) = 0.0;
817 
818  return rResult;
819  }
820 
828  {
829  if ( rResult.size() != this->PointsNumber() )
830  {
831  // KLUDGE: While there is a bug in
832  // ublas vector resize, I have to put this beside resizing!!
833 // ShapeFunctionsGradientsType
835  rResult.swap( temp );
836  }
837 
838  for ( IndexType i = 0; i < rResult.size(); i++ )
839  {
841  rResult[i].swap( temp );
842  }
843 
844  for ( unsigned int i = 0; i < this->PointsNumber(); i++ )
845  {
846  for ( unsigned int j = 0; j < 2; j++ )
847  {
848  rResult[i][j].resize( 2, 2 , false);
849  noalias( rResult[i][j] ) = ZeroMatrix( 2, 2 );
850  }
851  }
852 
853 // rResult(0,0).resize( 2, 2);
854 // rResult(0,1).resize( 2, 2);
855 // rResult(1,0).resize( 2, 2);
856 // rResult(1,1).resize( 2, 2);
857 // rResult(2,0).resize( 2, 2);
858 // rResult(2,1).resize( 2, 2);
859 // rResult(3,0).resize( 2, 2);
860 // rResult(3,0).resize( 2, 2);
861 
862  for ( int i = 0; i < 4; i++ )
863  {
864  rResult[i][0]( 0, 0 ) = 0.0;
865  rResult[i][0]( 0, 1 ) = 0.0;
866  rResult[i][0]( 1, 0 ) = 0.0;
867  rResult[i][0]( 1, 1 ) = 0.0;
868  rResult[i][1]( 0, 0 ) = 0.0;
869  rResult[i][1]( 0, 1 ) = 0.0;
870  rResult[i][1]( 1, 0 ) = 0.0;
871  rResult[i][1]( 1, 1 ) = 0.0;
872  }
873 
874  return rResult;
875  }
876 
880 
882 
883 protected:
888 private:
891 
892  static const GeometryData msGeometryData;
893 
894  static const GeometryDimension msGeometryDimension;
895 
899 
900 
904 
905  friend class Serializer;
906 
907  void save( Serializer& rSerializer ) const override
908  {
910  }
911 
912  void load( Serializer& rSerializer ) override
913  {
915  }
916 
917  Quadrilateral2D4(): BaseType( PointsArrayType(), &msGeometryData ) {}
918 
919 
923 
924 
928 
936  static Matrix CalculateShapeFunctionsIntegrationPointsValues(
937  typename BaseType::IntegrationMethod ThisMethod )
938  {
939  IntegrationPointsContainerType all_integration_points =
940  AllIntegrationPoints();
941  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
942  //number of integration points
943  const int integration_points_number = integration_points.size();
944  //number of nodes in current geometry
945  const int points_number = 4;
946  //setting up return matrix
947  Matrix shape_function_values( integration_points_number, points_number );
948  //loop over all integration points
949 
950  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
951  {
952  shape_function_values( pnt, 0 ) =
953  0.25 * ( 1.0 - integration_points[pnt].X() )
954  * ( 1.0 - integration_points[pnt].Y() );
955  shape_function_values( pnt, 1 ) =
956  0.25 * ( 1.0 + integration_points[pnt].X() )
957  * ( 1.0 - integration_points[pnt].Y() );
958  shape_function_values( pnt, 2 ) =
959  0.25 * ( 1.0 + integration_points[pnt].X() )
960  * ( 1.0 + integration_points[pnt].Y() );
961  shape_function_values( pnt, 3 ) =
962  0.25 * ( 1.0 - integration_points[pnt].X() )
963  * ( 1.0 + integration_points[pnt].Y() );
964  }
965 
966  return shape_function_values;
967  }
968 
979  CalculateShapeFunctionsIntegrationPointsLocalGradients(
980  typename BaseType::IntegrationMethod ThisMethod )
981  {
982  IntegrationPointsContainerType all_integration_points =
983  AllIntegrationPoints();
984  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
985  //number of integration points
986  const int integration_points_number = integration_points.size();
987  ShapeFunctionsGradientsType d_shape_f_values( integration_points_number );
988  //initialising container
989  //std::fill(d_shape_f_values.begin(), d_shape_f_values.end(), Matrix(4,2));
990  //loop over all integration points
991 
992  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
993  {
994  Matrix result( 4, 2 );
995  result( 0, 0 ) = -0.25 * ( 1.0 - integration_points[pnt].Y() );
996  result( 0, 1 ) = -0.25 * ( 1.0 - integration_points[pnt].X() );
997  result( 1, 0 ) = 0.25 * ( 1.0 - integration_points[pnt].Y() );
998  result( 1, 1 ) = -0.25 * ( 1.0 + integration_points[pnt].X() );
999  result( 2, 0 ) = 0.25 * ( 1.0 + integration_points[pnt].Y() );
1000  result( 2, 1 ) = 0.25 * ( 1.0 + integration_points[pnt].X() );
1001  result( 3, 0 ) = -0.25 * ( 1.0 + integration_points[pnt].Y() );
1002  result( 3, 1 ) = 0.25 * ( 1.0 - integration_points[pnt].X() );
1003  d_shape_f_values[pnt] = result;
1004  }
1005 
1006  return d_shape_f_values;
1007  }
1008 
1012  static const IntegrationPointsContainerType AllIntegrationPoints()
1013  {
1014  IntegrationPointsContainerType integration_points =
1015  {
1016  {
1017  Quadrature < QuadrilateralGaussLegendreIntegrationPoints1,
1018  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1019  Quadrature < QuadrilateralGaussLegendreIntegrationPoints2,
1020  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1021  Quadrature < QuadrilateralGaussLegendreIntegrationPoints3,
1022  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1023  Quadrature < QuadrilateralGaussLegendreIntegrationPoints4,
1024  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1025  Quadrature < QuadrilateralGaussLegendreIntegrationPoints5,
1026  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1027  Quadrature < QuadrilateralCollocationIntegrationPoints1,
1028  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1029  Quadrature < QuadrilateralCollocationIntegrationPoints2,
1030  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1031  Quadrature < QuadrilateralCollocationIntegrationPoints3,
1032  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1033  Quadrature < QuadrilateralCollocationIntegrationPoints4,
1034  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1035  Quadrature < QuadrilateralCollocationIntegrationPoints5,
1036  2, IntegrationPoint<3> >::GenerateIntegrationPoints()
1037  }
1038  };
1039  return integration_points;
1040  }
1041 
1045  static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
1046  {
1047  ShapeFunctionsValuesContainerType shape_functions_values =
1048  {
1049  {
1050  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1052  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1054  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1056  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1058  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1060  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1062  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1064  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1066  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1068  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1070  }
1071  };
1072  return shape_functions_values;
1073  }
1074 
1079  AllShapeFunctionsLocalGradients()
1080  {
1081  ShapeFunctionsLocalGradientsContainerType shape_functions_local_gradients =
1082  {
1083  {
1084  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_1 ),
1085  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_2 ),
1086  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_3 ),
1087  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_4 ),
1088  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_5 ),
1089  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_1 ),
1090  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_2 ),
1091  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_3 ),
1092  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_4 ),
1093  Quadrilateral2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_5 ),
1094  }
1095  };
1096  return shape_functions_local_gradients;
1097  }
1098 
1102 
1103 
1107 
1108 
1112 
1113  template<class TOtherPointType> friend class Quadrilateral2D4;
1114 
1118 
1119 
1121 }; // Class Geometry
1122 
1126 
1127 
1131 
1134 template<class TPointType> inline std::istream& operator >> (
1135  std::istream& rIStream,
1140 template<class TPointType> inline std::ostream& operator << (
1141  std::ostream& rOStream,
1142  const Quadrilateral2D4<TPointType>& rThis )
1143 {
1144  rThis.PrintInfo( rOStream );
1145  rOStream << std::endl;
1146  rThis.PrintData( rOStream );
1147  return rOStream;
1148 }
1149 
1151 
1152 template<class TPointType> const
1153 GeometryData Quadrilateral2D4<TPointType>::msGeometryData(
1156  Quadrilateral2D4<TPointType>::AllIntegrationPoints(),
1157  Quadrilateral2D4<TPointType>::AllShapeFunctionsValues(),
1158  AllShapeFunctionsLocalGradients()
1159 );
1160 
1161 template<class TPointType> const
1163 
1164 }// namespace Kratos.
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
SizeType IntegrationPointsNumber() const
Definition: geometry_data.h:430
IntegrationMethod
Definition: geometry_data.h:76
KratosGeometryFamily
Definition: geometry_data.h:91
IntegrationMethod DefaultIntegrationMethod() const
Definition: geometry_data.h:425
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
DataValueContainer & GetData()
Definition: geometry.h:591
std::size_t SizeType
Definition: geometry.h:144
const TPointType::Pointer pGetPoint(const int Index) const
Definition: geometry.h:1790
std::size_t IndexType
Definition: geometry.h:137
std::array< Matrix, static_cast< int >GeometryData::IntegrationMethod::NumberOfIntegrationMethods)> ShapeFunctionsValuesContainerType
Definition: geometry.h:172
virtual CoordinatesArrayType & PointLocalCoordinates(CoordinatesArrayType &rResult, const CoordinatesArrayType &rPoint) const
Returns the local coordinates of a given arbitrary point.
Definition: geometry.h:1854
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
JacobiansType & Jacobian(JacobiansType &rResult) const
Definition: geometry.h:2901
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
GeometryData::IntegrationMethod IntegrationMethod
Definition: geometry.h:122
Short class definition.
Definition: integration_point.h:52
static double ComputeArea2DGeometry(const Geometry< TPointType > &rGeometry)
This method calculates and returns the volume of the geometry from a 3D geometry.
Definition: integration_utilities.h:107
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
An two node 2D line geometry with linear shape functions.
Definition: line_2d_2.h:65
This class defines the node.
Definition: node.h:65
Point class.
Definition: point.h:59
double Y() const
Definition: point.h:187
double X() const
Definition: point.h:181
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
A four node 2D quadrilateral geometry with bi-linear shape functions.
Definition: quadrilateral_2d_4.h:73
bool IsInside(const CoordinatesArrayType &rPoint, CoordinatesArrayType &rResult, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Returns whether given arbitrary point is inside the Geometry and the respective local point for the g...
Definition: quadrilateral_2d_4.h:492
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_4.h:586
BaseType::IntegrationPointType IntegrationPointType
Definition: quadrilateral_2d_4.h:140
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: quadrilateral_2d_4.h:392
Quadrilateral2D4(Quadrilateral2D4 const &rOther)
Definition: quadrilateral_2d_4.h:267
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: quadrilateral_2d_4.h:168
virtual Matrix & ShapeFunctionsGradients(Matrix &rResult, PointType &rPoint)
Definition: quadrilateral_2d_4.h:763
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: quadrilateral_2d_4.h:182
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: quadrilateral_2d_4.h:134
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients(IntegrationMethod ThisMethod)
Definition: quadrilateral_2d_4.h:689
double Volume() const override
This method calculates and returns the volume of this geometry.
Definition: quadrilateral_2d_4.h:461
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: quadrilateral_2d_4.h:294
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients()
Definition: quadrilateral_2d_4.h:710
BaseType::IndexType IndexType
Definition: quadrilateral_2d_4.h:116
Quadrilateral2D4(const std::string &rGeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: quadrilateral_2d_4.h:250
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: quadrilateral_2d_4.h:299
void PrintInfo(std::ostream &rOStream) const override
Definition: quadrilateral_2d_4.h:652
friend class Quadrilateral2D4
Definition: quadrilateral_2d_4.h:1113
double Area() const override
This method calculates and returns area or surface area of this geometry depending to it's dimension.
Definition: quadrilateral_2d_4.h:449
BaseType::Pointer Create(const IndexType NewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_4.h:353
ShapeFunctionsThirdDerivativesType & ShapeFunctionsThirdDerivatives(ShapeFunctionsThirdDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_4.h:827
Quadrilateral2D4(typename PointType::Pointer pFirstPoint, typename PointType::Pointer pSecondPoint, typename PointType::Pointer pThirdPoint, typename PointType::Pointer pFourthPoint)
Definition: quadrilateral_2d_4.h:221
~Quadrilateral2D4() override
Definition: quadrilateral_2d_4.h:292
bool HasIntersection(const Point &rLowPoint, const Point &rHighPoint) const override
Definition: quadrilateral_2d_4.h:557
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: quadrilateral_2d_4.h:156
BaseType::NormalType NormalType
Definition: quadrilateral_2d_4.h:203
std::string Info() const override
Definition: quadrilateral_2d_4.h:641
KRATOS_CLASS_POINTER_DEFINITION(Quadrilateral2D4)
ShapeFunctionsSecondDerivativesType & ShapeFunctionsSecondDerivatives(ShapeFunctionsSecondDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_4.h:783
Line2D2< TPointType > EdgeType
Definition: quadrilateral_2d_4.h:87
BaseType::SizeType SizeType
Definition: quadrilateral_2d_4.h:123
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_4.h:367
BaseType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: quadrilateral_2d_4.h:190
Quadrilateral2D4 & operator=(const Quadrilateral2D4 &rOther)
Definition: quadrilateral_2d_4.h:319
Quadrilateral2D4(const IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: quadrilateral_2d_4.h:241
BaseType::PointsArrayType PointsArrayType
Definition: quadrilateral_2d_4.h:129
BaseType::ShapeFunctionsThirdDerivativesType ShapeFunctionsThirdDerivativesType
Definition: quadrilateral_2d_4.h:198
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: quadrilateral_2d_4.h:162
void PrintData(std::ostream &rOStream) const override
Definition: quadrilateral_2d_4.h:671
GeometryData::IntegrationMethod IntegrationMethod
Definition: quadrilateral_2d_4.h:97
Quadrilateral2D4(const PointsArrayType &rThisPoints)
Definition: quadrilateral_2d_4.h:233
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: quadrilateral_2d_4.h:526
TPointType PointType
Definition: quadrilateral_2d_4.h:108
SizeType PointsNumberInDirection(IndexType LocalDirectionIndex) const override
Returns number of points per direction.
Definition: quadrilateral_2d_4.h:378
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: quadrilateral_2d_4.h:617
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: quadrilateral_2d_4.h:149
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: quadrilateral_2d_4.h:539
Geometry< TPointType > BaseType
Definition: quadrilateral_2d_4.h:82
Quadrilateral2D4 & operator=(Quadrilateral2D4< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_4.h:337
BaseType::JacobiansType JacobiansType
Definition: quadrilateral_2d_4.h:175
Quadrilateral2D4(Quadrilateral2D4< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_4.h:284
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_4.h:736
double Length() const override
Definition: quadrilateral_2d_4.h:430
double DomainSize() const override
This method calculates and returns length, area or volume of this geometry depending to it's dimensio...
Definition: quadrilateral_2d_4.h:479
BaseType::GeometriesArrayType GeometriesArrayType
Definition: quadrilateral_2d_4.h:103
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
A three node 2D triangle geometry with linear shape functions.
Definition: triangle_2d_3.h:74
bool HasIntersection(const BaseType &rThisGeometry) const override
Detect if this triangle is intersected with another geometry.
Definition: triangle_2d_3.h:518
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
#define KRATOS_ERROR_IF(conditional)
Definition: exception.h:162
#define KRATOS_WARNING(label)
Definition: logger.h:265
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
const GeometryData Quadrilateral2D4< TPointType >::msGeometryData & msGeometryDimension(), Quadrilateral2D4< 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
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
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
int j
Definition: quadrature.py:648
float temp
Definition: rotating_cone.py:85
integer i
Definition: TensorModule.f:17