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.
prism_3d_6.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
28 #include "utilities/geometry_utilities.h"
29 
30 namespace Kratos
31 {
59 template<class TPointType> class Prism3D6 : public Geometry<TPointType>
60 {
61 public:
70 
77 
82 
87 
93 
97  typedef TPointType PointType;
98 
104  typedef typename BaseType::IndexType IndexType;
105 
106 
112  typedef typename BaseType::SizeType SizeType;
113 
119 
126 
134 
141 
148 
155 
162 
169 
174 
179 
184 
185 
190 // Prism3D6( const PointType& Point1, const PointType& Point2,
191 // const PointType& Point3, const PointType& Point4,
192 // const PointType& Point5, const PointType& Point6 )
193 // : BaseType( PointsArrayType(), &msGeometryData )
194 // {
195 // this->Points().reserve( 6 );
196 // this->Points().push_back( typename PointType::Pointer( new PointType( Point1 ) ) );
197 // this->Points().push_back( typename PointType::Pointer( new PointType( Point2 ) ) );
198 // this->Points().push_back( typename PointType::Pointer( new PointType( Point3 ) ) );
199 // this->Points().push_back( typename PointType::Pointer( new PointType( Point4 ) ) );
200 // this->Points().push_back( typename PointType::Pointer( new PointType( Point5 ) ) );
201 // this->Points().push_back( typename PointType::Pointer( new PointType( Point6 ) ) );
202 // }
203 
204  Prism3D6( typename PointType::Pointer pPoint1,
205  typename PointType::Pointer pPoint2,
206  typename PointType::Pointer pPoint3,
207  typename PointType::Pointer pPoint4,
208  typename PointType::Pointer pPoint5,
209  typename PointType::Pointer pPoint6 )
210  : BaseType( PointsArrayType(), &msGeometryData )
211  {
212  this->Points().reserve( 6 );
213  this->Points().push_back( pPoint1 );
214  this->Points().push_back( pPoint2 );
215  this->Points().push_back( pPoint3 );
216  this->Points().push_back( pPoint4 );
217  this->Points().push_back( pPoint5 );
218  this->Points().push_back( pPoint6 );
219  }
220 
221  explicit Prism3D6( const PointsArrayType& rThisPoints )
222  : BaseType( rThisPoints, &msGeometryData )
223  {
224  KRATOS_ERROR_IF( this->PointsNumber() != 6 ) << "Invalid points number. Expected 6, given " << this->PointsNumber() << std::endl;
225  }
226 
228  explicit Prism3D6(
229  const IndexType GeometryId,
230  const PointsArrayType& rThisPoints
231  ) : BaseType( GeometryId, rThisPoints, &msGeometryData )
232  {
233  KRATOS_ERROR_IF( this->PointsNumber() != 6 ) << "Invalid points number. Expected 6, given " << this->PointsNumber() << std::endl;
234  }
235 
237  explicit Prism3D6(
238  const std::string& rGeometryName,
239  const PointsArrayType& rThisPoints
240  ) : BaseType(rGeometryName, rThisPoints, &msGeometryData)
241  {
242  KRATOS_ERROR_IF(this->PointsNumber() != 6) << "Invalid points number. Expected 6, given " << this->PointsNumber() << std::endl;
243  }
244 
254  Prism3D6( Prism3D6 const& rOther )
255  : BaseType( rOther )
256  {
257  }
258 
271  template<class TOtherPointType> explicit Prism3D6( Prism3D6<TOtherPointType> const& rOther )
272  : BaseType( rOther )
273  {
274  }
275 
277  ~Prism3D6() override {}
278 
280  {
282  }
283 
285  {
287  }
288 
304  Prism3D6& operator=( const Prism3D6& rOther )
305  {
306  BaseType::operator=( rOther );
307  return *this;
308  }
309 
321  template<class TOtherPointType>
323  {
324  BaseType::operator=( rOther );
325 
326  return *this;
327  }
328 
329 
340  typename BaseType::Pointer Create(
341  const IndexType NewGeometryId,
342  PointsArrayType const& rThisPoints
343  ) const override
344  {
345  return typename BaseType::Pointer( new Prism3D6( NewGeometryId, rThisPoints ) );
346  }
347 
354  typename BaseType::Pointer Create(
355  const IndexType NewGeometryId,
356  const BaseType& rGeometry
357  ) const override
358  {
359  auto p_geometry = typename BaseType::Pointer( new Prism3D6( NewGeometryId, rGeometry.Points() ) );
360  p_geometry->SetData(rGeometry.GetData());
361  return p_geometry;
362  }
363 
383  double Length() const override
384  {
385  const double volume = Volume();
386 
387  return std::pow(volume, 1.0/3.0)/3.0;
388 // return std::sqrt( fabs( this->DeterminantOfJacobian( PointType() ) ) );
389  }
390 
405  double Area() const override
406  {
407  return std::abs( this->DeterminantOfJacobian( PointType() ) ) * 0.5;
408  }
409 
418  double Volume() const override
419  {
421  }
422 
436  double DomainSize() const override
437  {
438  return Volume();
439  }
440 
446  Matrix& PointsLocalCoordinates( Matrix& rResult ) const override
447  {
448  if ( rResult.size1() != 6 || rResult.size2() != 3 )
449  rResult.resize( 6, 3 ,false);
450 
451  rResult( 0, 0 ) = 0.0;
452  rResult( 0, 1 ) = 0.0;
453  rResult( 0, 2 ) = 0.0;
454 
455  rResult( 1, 0 ) = 1.0;
456  rResult( 1, 1 ) = 0.0;
457  rResult( 1, 2 ) = 0.0;
458 
459  rResult( 2, 0 ) = 0.0;
460  rResult( 2, 1 ) = 1.0;
461  rResult( 2, 2 ) = 0.0;
462 
463  rResult( 3, 0 ) = 0.0;
464  rResult( 3, 1 ) = 0.0;
465  rResult( 3, 2 ) = 1.0;
466 
467  rResult( 4, 0 ) = 1.0;
468  rResult( 4, 1 ) = 0.0;
469  rResult( 4, 2 ) = 1.0;
470 
471  rResult( 5, 0 ) = 0.0;
472  rResult( 5, 1 ) = 1.0;
473  rResult( 5, 2 ) = 1.0;
474 
475  return rResult;
476  }
477 
486  bool IsInside(
487  const CoordinatesArrayType& rPoint,
488  CoordinatesArrayType& rResult,
489  const double Tolerance = std::numeric_limits<double>::epsilon()
490  ) const override
491  {
492  this->PointLocalCoordinates( rResult, rPoint );
493 
494  if ( (rResult[0] >= (0.0 - Tolerance)) && (rResult[0] <= (1.0 + Tolerance)) )
495  if ( (rResult[1] >= (0.0 - Tolerance)) && (rResult[1] <= (1.0 + Tolerance)) )
496  if ( (rResult[2] >= (0.0 - Tolerance)) && (rResult[2] <= (1.0 + Tolerance)) )
497  if ((( 1.0 - ( rResult[0] + rResult[1] ) ) >= (0.0 - Tolerance) ) && (( 1.0 - ( rResult[0] + rResult[1] ) ) <= (1.0 + Tolerance) ) )
498  return true;
499 
500  return false;
501  }
502 
506 
518  SizeType EdgesNumber() const override
519  {
520  return 9;
521  }
522 
532  {
534  typedef typename Geometry<TPointType>::Pointer EdgePointerType;
535  edges.push_back( EdgePointerType( new EdgeType(
536  this->pGetPoint( 0 ),
537  this->pGetPoint( 1 ) ) ) );
538  edges.push_back( EdgePointerType( new EdgeType(
539  this->pGetPoint( 1 ),
540  this->pGetPoint( 2 ) ) ) );
541  edges.push_back( EdgePointerType( new EdgeType(
542  this->pGetPoint( 2 ),
543  this->pGetPoint( 0 ) ) ) );
544  edges.push_back( EdgePointerType( new EdgeType(
545  this->pGetPoint( 3 ),
546  this->pGetPoint( 4 ) ) ) );
547  edges.push_back( EdgePointerType( new EdgeType(
548  this->pGetPoint( 4 ),
549  this->pGetPoint( 5 ) ) ) );
550  edges.push_back( EdgePointerType( new EdgeType(
551  this->pGetPoint( 5 ),
552  this->pGetPoint( 3 ) ) ) );
553  edges.push_back( EdgePointerType( new EdgeType(
554  this->pGetPoint( 0 ),
555  this->pGetPoint( 3 ) ) ) );
556  edges.push_back( EdgePointerType( new EdgeType(
557  this->pGetPoint( 1 ),
558  this->pGetPoint( 4 ) ) ) );
559  edges.push_back( EdgePointerType( new EdgeType(
560  this->pGetPoint( 2 ),
561  this->pGetPoint( 5 ) ) ) );
562  return edges;
563  }
564 
568 
576  SizeType FacesNumber() const override
577  {
578  return 5;
579  }
580 
590  {
592  typedef typename Geometry<TPointType>::Pointer FacePointerType;
593  faces.push_back( FacePointerType( new FaceType1(
594  this->pGetPoint( 0 ),
595  this->pGetPoint( 2 ),
596  this->pGetPoint( 1 ) ) ) );
597  faces.push_back( FacePointerType( new FaceType1(
598  this->pGetPoint( 3 ),
599  this->pGetPoint( 4 ),
600  this->pGetPoint( 5 ) ) ) );
601  faces.push_back( FacePointerType( new FaceType2(
602  this->pGetPoint( 1 ),
603  this->pGetPoint( 2 ),
604  this->pGetPoint( 5 ),
605  this->pGetPoint( 4 ) ) ) );
606  faces.push_back( FacePointerType( new FaceType2(
607  this->pGetPoint( 0 ),
608  this->pGetPoint( 3 ),
609  this->pGetPoint( 5 ),
610  this->pGetPoint( 2 ) ) ) );
611  faces.push_back( FacePointerType( new FaceType2(
612  this->pGetPoint( 0 ),
613  this->pGetPoint( 1 ),
614  this->pGetPoint( 4 ),
615  this->pGetPoint( 3 ) ) ) );
616  return faces;
617  }
618 
619  bool HasIntersection( const Point& rLowPoint, const Point& rHighPoint ) const override
620  {
621  // Check if faces have intersection
622  if(FaceType1(this->pGetPoint(0),this->pGetPoint(2), this->pGetPoint(1)).HasIntersection(rLowPoint, rHighPoint))
623  return true;
624  if(FaceType1(this->pGetPoint(3),this->pGetPoint(4), this->pGetPoint(5)).HasIntersection(rLowPoint, rHighPoint))
625  return true;
626  if(FaceType2(this->pGetPoint(1),this->pGetPoint(2), this->pGetPoint(5), this->pGetPoint(4)).HasIntersection(rLowPoint, rHighPoint))
627  return true;
628  if(FaceType2(this->pGetPoint(0),this->pGetPoint(3), this->pGetPoint(5), this->pGetPoint(2)).HasIntersection(rLowPoint, rHighPoint))
629  return true;
630  if(FaceType2(this->pGetPoint(0),this->pGetPoint(1), this->pGetPoint(4), this->pGetPoint(3)).HasIntersection(rLowPoint, rHighPoint))
631  return true;
632 
633  CoordinatesArrayType local_coordinates;
634  // If there are no faces intersecting the box then or the box is inside the hexahedron or it does not have intersection
635  if(IsInside(rLowPoint,local_coordinates))
636  return true;
637 
638  return false;
639  }
640 
655  double ShapeFunctionValue( IndexType ShapeFunctionIndex,
656  const CoordinatesArrayType& rPoint ) const override
657  {
658  switch ( ShapeFunctionIndex )
659  {
660  case 0:
661  return( 1.0 -( rPoint[0] + rPoint[1] + rPoint[2] - ( rPoint[0]*rPoint[2] ) - ( rPoint[1]*rPoint[2] ) ) );
662  case 1:
663  return( rPoint[0] - ( rPoint[0]*rPoint[2] ) );
664  case 2:
665  return( rPoint[1] - ( rPoint[1]*rPoint[2] ) );
666  case 3:
667  return( rPoint[2] - ( rPoint[0]*rPoint[2] ) - ( rPoint[1]*rPoint[2] ) );
668  case 4:
669  return( rPoint[0]*rPoint[2] );
670  case 5:
671  return( rPoint[1]*rPoint[2] );
672  default:
673  KRATOS_ERROR << "Wrong index of shape function!" << *this << std::endl;
674  }
675 
676  return 0;
677  }
678 
688  Vector &rResult,
689  const CoordinatesArrayType& rCoordinates
690  ) const override
691  {
692  if(rResult.size() != 6)
693  rResult.resize(6,false);
694 
695  rResult[0] = 1.0 -( rCoordinates[0] + rCoordinates[1] + rCoordinates[2] - ( rCoordinates[0] * rCoordinates[2] ) - ( rCoordinates[1] * rCoordinates[2] ) );
696  rResult[1] = rCoordinates[0] - ( rCoordinates[0] * rCoordinates[2] );
697  rResult[2] = rCoordinates[1] - ( rCoordinates[1] * rCoordinates[2] );
698  rResult[3] = rCoordinates[2] - ( rCoordinates[0] * rCoordinates[2] ) - ( rCoordinates[1] * rCoordinates[2] );
699  rResult[4] = rCoordinates[0] * rCoordinates[2];
700  rResult[5] = rCoordinates[1] * rCoordinates[2];
701 
702  return rResult;
703  }
704 
713  Matrix& rResult,
714  const CoordinatesArrayType& rPoint
715  ) const override
716  {
717  if(rResult.size1() != this->PointsNumber() || rResult.size2() != this->LocalSpaceDimension())
718  rResult.resize(this->PointsNumber(),this->LocalSpaceDimension(),false);
719 
720  CalculateShapeFunctionsLocalGradients(rResult, rPoint);
721 
722  return rResult;
723  }
724 
728 
742  const CoordinatesArrayType& rPointGlobalCoordinates,
743  const double Tolerance = std::numeric_limits<double>::epsilon()
744  ) const override
745  {
746  // Point to compute distance to
747  const Point point(rPointGlobalCoordinates);
748 
749  // Check if point is inside
750  CoordinatesArrayType aux_coordinates;
751  if (this->IsInside(rPointGlobalCoordinates, aux_coordinates, Tolerance)) {
752  return 0.0;
753  }
754 
755  // Compute distance to faces
756  std::array<double, 5> distances;
757  distances[0] = GeometryUtils::PointDistanceToTriangle3D(this->GetPoint(0), this->GetPoint(2), this->GetPoint(1), point);
758  distances[1] = GeometryUtils::PointDistanceToTriangle3D(this->GetPoint(3), this->GetPoint(4), this->GetPoint(5), point);
759  distances[2] = GeometryUtils::PointDistanceToQuadrilateral3D(this->GetPoint(1), this->GetPoint(2), this->GetPoint(5), this->GetPoint(4), point);
760  distances[3] = GeometryUtils::PointDistanceToQuadrilateral3D(this->GetPoint(0), this->GetPoint(3), this->GetPoint(5), this->GetPoint(2), point);
761  distances[4] = GeometryUtils::PointDistanceToQuadrilateral3D(this->GetPoint(0), this->GetPoint(1), this->GetPoint(4), this->GetPoint(3), point);
762  const auto min = std::min_element(distances.begin(), distances.end());
763  return *min;
764  }
765 
769 
777  std::string Info() const override
778  {
779  return "3 dimensional prism with six nodes in 3D space";
780  }
781 
789  void PrintInfo( std::ostream& rOStream ) const override
790  {
791  rOStream << "3 dimensional prism with six nodes in 3D space";
792  }
793 
803  void PrintData( std::ostream& rOStream ) const override
804  {
805  // Base Geometry class PrintData call
806  BaseType::PrintData( rOStream );
807  std::cout << std::endl;
808 
809  // If the geometry has valid points, calculate and output its data
810  if (this->AllPointsAreValid()) {
811  Matrix jacobian;
812  this->Jacobian( jacobian, PointType() );
813  rOStream << " Jacobian in the origin\t : " << jacobian;
814  }
815  }
816 
817 protected:
818 
823 private:
826 
827  static const GeometryData msGeometryData;
828 
829  static const GeometryDimension msGeometryDimension;
830 
834 
835  friend class Serializer;
836 
837  void save( Serializer& rSerializer ) const override
838  {
840  }
841 
842  void load( Serializer& rSerializer ) override
843  {
845  }
846 
847  Prism3D6(): BaseType( PointsArrayType(), &msGeometryData ) {}
848 
849 
865  static Matrix& CalculateShapeFunctionsLocalGradients(
866  Matrix& rResult,
867  const CoordinatesArrayType& rPoint
868  )
869  {
870  rResult( 0, 0 ) = -1.0 + rPoint[2];
871  rResult( 0, 1 ) = -1.0 + rPoint[2];
872  rResult( 0, 2 ) = -1.0 + rPoint[0] + rPoint[1];
873  rResult( 1, 0 ) = 1.0 - rPoint[2];
874  rResult( 1, 1 ) = 0.0;
875  rResult( 1, 2 ) = -rPoint[0];
876  rResult( 2, 0 ) = 0.0;
877  rResult( 2, 1 ) = 1.0 - rPoint[2];
878  rResult( 2, 2 ) = -rPoint[1];
879  rResult( 3, 0 ) = -rPoint[2];
880  rResult( 3, 1 ) = -rPoint[2];
881  rResult( 3, 2 ) = 1.0 - rPoint[0] - rPoint[1];
882  rResult( 4, 0 ) = rPoint[2];
883  rResult( 4, 1 ) = 0.0;
884  rResult( 4, 2 ) = rPoint[0];
885  rResult( 5, 0 ) = 0.0;
886  rResult( 5, 1 ) = rPoint[2];
887  rResult( 5, 2 ) = rPoint[1];
888  return rResult;
889  }
890 
902  static Matrix CalculateShapeFunctionsIntegrationPointsValues(typename BaseType::IntegrationMethod ThisMethod )
903  {
904  IntegrationPointsContainerType all_integration_points = AllIntegrationPoints();
905  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
906  //number of integration points
907  const int integration_points_number = integration_points.size();
908  //number of nodes in current geometry
909  const int points_number = 6;
910  //setting up return matrix
911  Matrix shape_function_values( integration_points_number, points_number );
912  //loop over all integration points
913 
914  for ( int pnt = 0; pnt < integration_points_number; pnt++ ) {
915  shape_function_values( pnt, 0 ) = ( 1.0
916  - integration_points[pnt].X()
917  - integration_points[pnt].Y()
918  - integration_points[pnt].Z()
919  + ( integration_points[pnt].X() * integration_points[pnt].Z() )
920  + ( integration_points[pnt].Y() * integration_points[pnt].Z() ) );
921  shape_function_values( pnt, 1 ) = integration_points[pnt].X()
922  - ( integration_points[pnt].X() * integration_points[pnt].Z() );
923  shape_function_values( pnt, 2 ) = integration_points[pnt].Y()
924  - ( integration_points[pnt].Y() * integration_points[pnt].Z() );
925  shape_function_values( pnt, 3 ) = integration_points[pnt].Z()
926  - ( integration_points[pnt].X() * integration_points[pnt].Z() )
927  - ( integration_points[pnt].Y() * integration_points[pnt].Z() );
928  shape_function_values( pnt, 4 ) = ( integration_points[pnt].X() * integration_points[pnt].Z() );
929  shape_function_values( pnt, 5 ) = ( integration_points[pnt].Y() * integration_points[pnt].Z() );
930  }
931 
932  return shape_function_values;
933  }
934 
948  CalculateShapeFunctionsIntegrationPointsLocalGradients(
949  typename BaseType::IntegrationMethod ThisMethod )
950  {
951  IntegrationPointsContainerType all_integration_points =
952  AllIntegrationPoints();
953  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
954  //number of integration points
955  const int integration_points_number = integration_points.size();
956  ShapeFunctionsGradientsType d_shape_f_values( integration_points_number );
957  //initialising container
958  //loop over all integration points
959 
960  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
961  {
962  Matrix result = ZeroMatrix( 6, 3 );
963  result( 0, 0 ) = -1.0 + integration_points[pnt].Z();
964  result( 0, 1 ) = -1.0 + integration_points[pnt].Z();
965  result( 0, 2 ) = -1.0 + integration_points[pnt].X() + integration_points[pnt].Y();
966  result( 1, 0 ) = 1.0 - integration_points[pnt].Z();
967  result( 1, 1 ) = 0.0;
968  result( 1, 2 ) = -integration_points[pnt].X();
969  result( 2, 0 ) = 0.0;
970  result( 2, 1 ) = 1.0 - integration_points[pnt].Z();
971  result( 2, 2 ) = -integration_points[pnt].Y();
972  result( 3, 0 ) = -integration_points[pnt].Z();
973  result( 3, 1 ) = -integration_points[pnt].Z();
974  result( 3, 2 ) = 1.0 - integration_points[pnt].X() - integration_points[pnt].Y();
975  result( 4, 0 ) = integration_points[pnt].Z();
976  result( 4, 1 ) = 0.0;
977  result( 4, 2 ) = integration_points[pnt].X();
978  result( 5, 0 ) = 0.0;
979  result( 5, 1 ) = integration_points[pnt].Z();
980  result( 5, 2 ) = integration_points[pnt].Y();
981  d_shape_f_values[pnt] = result;
982  }
983 
984  return d_shape_f_values;
985  }
986 
987  static const IntegrationPointsContainerType AllIntegrationPoints()
988  {
989  IntegrationPointsContainerType integration_points =
990  {
991  {
992  Quadrature < PrismGaussLegendreIntegrationPoints1,
993  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
994  Quadrature < PrismGaussLegendreIntegrationPoints2,
995  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
996  Quadrature < PrismGaussLegendreIntegrationPoints3,
997  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
998  Quadrature < PrismGaussLegendreIntegrationPoints4,
999  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1000  Quadrature < PrismGaussLegendreIntegrationPoints5,
1001  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1002  Quadrature < PrismGaussLegendreIntegrationPointsExt1,
1003  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1004  Quadrature < PrismGaussLegendreIntegrationPointsExt2,
1005  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1006  Quadrature < PrismGaussLegendreIntegrationPointsExt3,
1007  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1008  Quadrature < PrismGaussLegendreIntegrationPointsExt4,
1009  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1010  Quadrature < PrismGaussLegendreIntegrationPointsExt5,
1011  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1012  }
1013  };
1014  return integration_points;
1015  }
1016 
1017  static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
1018  {
1019  ShapeFunctionsValuesContainerType shape_functions_values =
1020  {
1021  {
1022  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1024  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1026  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1028  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1030  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1032  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1034  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1036  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1038  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1040  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1042  }
1043  };
1044  return shape_functions_values;
1045  }
1046 
1051  AllShapeFunctionsLocalGradients()
1052  {
1053  ShapeFunctionsLocalGradientsContainerType shape_functions_local_gradients =
1054  {
1055  {
1056  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1058  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1060  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1062  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1064  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1066  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1068  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1070  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1072  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1074  Prism3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1076  }
1077  };
1078  return shape_functions_local_gradients;
1079  }
1080 
1081 
1086  template<class TOtherPointType> friend class Prism3D6;
1087 
1088 
1093 };// Class Prism3D6
1094 
1095 
1101 template<class TPointType> inline std::istream& operator >> (
1102  std::istream& rIStream, Prism3D6<TPointType>& rThis );
1103 
1105 template<class TPointType> inline std::ostream& operator << (
1106  std::ostream& rOStream, const Prism3D6<TPointType>& rThis )
1107 {
1108  rThis.PrintInfo( rOStream );
1109  rOStream << std::endl;
1110  rThis.PrintData( rOStream );
1111 
1112  return rOStream;
1113 }
1114 
1115 template<class TPointType> const
1116 GeometryData Prism3D6<TPointType>::msGeometryData(
1119  Prism3D6<TPointType>::AllIntegrationPoints(),
1120  Prism3D6<TPointType>::AllShapeFunctionsValues(),
1121  AllShapeFunctionsLocalGradients()
1122 );
1123 
1124 template<class TPointType> const
1125 GeometryDimension Prism3D6<TPointType>::msGeometryDimension(3, 3);
1126 
1127 }// namespace Kratos.
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
DataValueContainer & GetData()
Definition: geometry.h:591
SizeType LocalSpaceDimension() const
Definition: geometry.h:1300
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
Vector & DeterminantOfJacobian(Vector &rResult) const
Definition: geometry.h:3154
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
TPointType const & GetPoint(const int Index) const
Definition: geometry.h:1816
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
GeometryData::IntegrationMethod IntegrationMethod
Definition: geometry.h:122
static double PointDistanceToTriangle3D(const Point &rTrianglePoint1, const Point &rTrianglePoint2, const Point &rTrianglePoint3, const Point &rPoint)
This function calculates the distance of a 3D point to a 3D triangle.
Definition: geometry_utilities.cpp:172
static double PointDistanceToQuadrilateral3D(const Point &rQuadrilateralPoint1, const Point &rQuadrilateralPoint2, const Point &rQuadrilateralPoint3, const Point &rQuadrilateralPoint4, const Point &rPoint)
This function calculates the distance of a 3D point to a 3D quadrilateral.
Definition: geometry_utilities.cpp:376
Short class definition.
Definition: integration_point.h:52
static double ComputeVolume3DGeometry(const Geometry< TPointType > &rGeometry)
This method calculates and returns the volume of the geometry from a 3D geometry.
Definition: integration_utilities.h:168
Definition: amatrix_interface.h:41
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
An two node 3D line geometry with linear shape functions.
Definition: line_3d_2.h:64
Point class.
Definition: point.h:59
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
void reserve(size_type dim)
Definition: pointer_vector.h:319
void push_back(const TPointerType &x)
Definition: pointer_vector.h:270
A six node prism geometry with linear shape functions.
Definition: prism_3d_6.h:60
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: prism_3d_6.h:655
BaseType::GeometriesArrayType GeometriesArrayType
Definition: prism_3d_6.h:92
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: prism_3d_6.h:284
Triangle3D3< TPointType > FaceType1
Definition: prism_3d_6.h:75
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: prism_3d_6.h:140
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: prism_3d_6.h:178
BaseType::SizeType SizeType
Definition: prism_3d_6.h:112
Prism3D6(const IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: prism_3d_6.h:228
Line3D2< TPointType > EdgeType
Definition: prism_3d_6.h:74
BaseType::IndexType IndexType
Definition: prism_3d_6.h:104
bool IsInside(const CoordinatesArrayType &rPoint, CoordinatesArrayType &rResult, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Definition: prism_3d_6.h:486
double Length() const override
Definition: prism_3d_6.h:383
Prism3D6(const PointsArrayType &rThisPoints)
Definition: prism_3d_6.h:221
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: prism_3d_6.h:354
void PrintInfo(std::ostream &rOStream) const override
Definition: prism_3d_6.h:789
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: prism_3d_6.h:518
~Prism3D6() override
Destructor. Does nothing!!!
Definition: prism_3d_6.h:277
BaseType::JacobiansType JacobiansType
Definition: prism_3d_6.h:161
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: prism_3d_6.h:712
Prism3D6 & operator=(Prism3D6< TOtherPointType > const &rOther)
Definition: prism_3d_6.h:322
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: prism_3d_6.h:168
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: prism_3d_6.h:687
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: prism_3d_6.h:446
Geometry< TPointType > BaseType
Definition: prism_3d_6.h:69
void PrintData(std::ostream &rOStream) const override
Definition: prism_3d_6.h:803
Prism3D6 & operator=(const Prism3D6 &rOther)
Definition: prism_3d_6.h:304
Matrix MatrixType
Definition: prism_3d_6.h:183
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: prism_3d_6.h:133
BaseType::PointsArrayType PointsArrayType
Definition: prism_3d_6.h:118
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: prism_3d_6.h:531
bool HasIntersection(const Point &rLowPoint, const Point &rHighPoint) const override
Definition: prism_3d_6.h:619
std::string Info() const override
Definition: prism_3d_6.h:777
friend class Prism3D6
Definition: prism_3d_6.h:1086
Quadrilateral3D4< TPointType > FaceType2
Definition: prism_3d_6.h:76
double Area() const override
Definition: prism_3d_6.h:405
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: prism_3d_6.h:279
TPointType PointType
Definition: prism_3d_6.h:97
BaseType::IntegrationPointType IntegrationPointType
Definition: prism_3d_6.h:125
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: prism_3d_6.h:147
double DomainSize() const override
Definition: prism_3d_6.h:436
KRATOS_CLASS_POINTER_DEFINITION(Prism3D6)
Prism3D6(const std::string &rGeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: prism_3d_6.h:237
double CalculateDistance(const CoordinatesArrayType &rPointGlobalCoordinates, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Computes the distance between an point in global coordinates and the closest point of this geometry....
Definition: prism_3d_6.h:741
BaseType::Pointer Create(const IndexType NewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: prism_3d_6.h:340
Prism3D6(Prism3D6< TOtherPointType > const &rOther)
Definition: prism_3d_6.h:271
Prism3D6(typename PointType::Pointer pPoint1, typename PointType::Pointer pPoint2, typename PointType::Pointer pPoint3, typename PointType::Pointer pPoint4, typename PointType::Pointer pPoint5, typename PointType::Pointer pPoint6)
Definition: prism_3d_6.h:204
Prism3D6(Prism3D6 const &rOther)
Definition: prism_3d_6.h:254
SizeType FacesNumber() const override
Returns the number of faces of the current geometry.
Definition: prism_3d_6.h:576
GeometriesArrayType GenerateFaces() const override
Returns all faces of the current geometry.
Definition: prism_3d_6.h:589
GeometryData::IntegrationMethod IntegrationMethod
Definition: prism_3d_6.h:86
double Volume() const override
This method calculate and return volume of this geometry.
Definition: prism_3d_6.h:418
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: prism_3d_6.h:154
BaseType::NormalType NormalType
Definition: prism_3d_6.h:173
A four node 3D quadrilateral geometry with bi-linear shape functions.
Definition: quadrilateral_3d_4.h:76
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 3D triangle geometry with linear shape functions.
Definition: triangle_3d_3.h:77
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
static double min(double a, double b)
Definition: GeometryFunctions.h:71
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
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
const GeometryData Prism3D6< TPointType >::msGeometryData & msGeometryDimension(), Prism3D6< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
def load(f)
Definition: ode_solve.py:307
volume
Definition: sp_statistics.py:15