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.
triangle_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
24 #include "geometries/line_3d_3.h"
26 #include "utilities/geometry_utilities.h"
27 
28 namespace Kratos
29 {
32 
36 
40 
44 
48 
66 template<class TPointType> class Triangle3D6
67  : public Geometry<TPointType>
68 {
69 public:
73 
78 
83 
88 
93 
99 
103  typedef TPointType PointType;
104 
111  typedef typename BaseType::IndexType IndexType;
112 
118  typedef typename BaseType::SizeType SizeType;
119 
125 
130 
136 
145 
152 
158 
164 
171 
178 
186 
194 
199 
203 
204 // Triangle3D6( const PointType& FirstPoint,
205 // const PointType& SecondPoint,
206 // const PointType& ThirdPoint,
207 // const PointType& FourthPoint,
208 // const PointType& FifthPoint,
209 // const PointType& SixthPoint
210 // )
211 // : BaseType( PointsArrayType(), &msGeometryData )
212 // {
213 // this->Points().push_back( typename PointType::Pointer( new PointType( FirstPoint ) ) );
214 // this->Points().push_back( typename PointType::Pointer( new PointType( SecondPoint ) ) );
215 // this->Points().push_back( typename PointType::Pointer( new PointType( ThirdPoint ) ) );
216 // this->Points().push_back( typename PointType::Pointer( new PointType( FourthPoint ) ) );
217 // this->Points().push_back( typename PointType::Pointer( new PointType( FifthPoint ) ) );
218 // this->Points().push_back( typename PointType::Pointer( new PointType( SixthPoint ) ) );
219 // }
220 
221  Triangle3D6( typename PointType::Pointer pFirstPoint,
222  typename PointType::Pointer pSecondPoint,
223  typename PointType::Pointer pThirdPoint,
224  typename PointType::Pointer pFourthPoint,
225  typename PointType::Pointer pFifthPoint,
226  typename PointType::Pointer pSixthPoint
227  )
228  : BaseType( PointsArrayType(), &msGeometryData )
229  {
230  this->Points().push_back( pFirstPoint );
231  this->Points().push_back( pSecondPoint );
232  this->Points().push_back( pThirdPoint );
233  this->Points().push_back( pFourthPoint );
234  this->Points().push_back( pFifthPoint );
235  this->Points().push_back( pSixthPoint );
236  }
237 
238  explicit Triangle3D6(
239  const PointsArrayType& ThisPoints
240  ) : BaseType( ThisPoints, &msGeometryData )
241  {
242  KRATOS_ERROR_IF( this->PointsNumber() != 6 ) << "Invalid points number. Expected 6, given " << this->PointsNumber() << std::endl;
243  }
244 
246  explicit Triangle3D6(
247  const IndexType GeometryId,
248  const PointsArrayType& rThisPoints
249  ) : BaseType(GeometryId, rThisPoints, &msGeometryData)
250  {
251  KRATOS_ERROR_IF( this->PointsNumber() != 6 ) << "Invalid points number. Expected 6, given " << this->PointsNumber() << std::endl;
252  }
253 
255  explicit Triangle3D6(
256  const std::string& rGeometryName,
257  const PointsArrayType& rThisPoints
258  ) : BaseType(rGeometryName, rThisPoints, &msGeometryData)
259  {
260  KRATOS_ERROR_IF(this->PointsNumber() != 20) << "Invalid points number. Expected 20, given " << this->PointsNumber() << std::endl;
261  }
262 
272  Triangle3D6( Triangle3D6 const& rOther )
273  : BaseType( rOther )
274  {
275  }
276 
289  template<class TOtherPointType> Triangle3D6( Triangle3D6<TOtherPointType> const& rOther )
290  : BaseType( rOther )
291  {
292  }
293 
297  ~Triangle3D6() override {}
298 
300  {
302  }
303 
305  {
307  }
308 
312 
325  {
326  BaseType::operator=( rOther );
327  return *this;
328  }
329 
341  template<class TOtherPointType>
343  {
344  BaseType::operator=( rOther );
345  return *this;
346  }
347 
351 
357  typename BaseType::Pointer Create(
358  PointsArrayType const& rThisPoints
359  ) const override
360  {
361  return typename BaseType::Pointer( new Triangle3D6( rThisPoints ) );
362  }
363 
370  typename BaseType::Pointer Create(
371  const IndexType NewGeometryId,
372  PointsArrayType const& rThisPoints
373  ) const override
374  {
375  return typename BaseType::Pointer( new Triangle3D6( NewGeometryId, rThisPoints ) );
376  }
377 
383  typename BaseType::Pointer Create(
384  const BaseType& rGeometry
385  ) const override
386  {
387  auto p_geometry = typename BaseType::Pointer( new Triangle3D6( rGeometry.Points() ) );
388  p_geometry->SetData(rGeometry.GetData());
389  return p_geometry;
390  }
391 
398  typename BaseType::Pointer Create(
399  const IndexType NewGeometryId,
400  const BaseType& rGeometry
401  ) const override
402  {
403  auto p_geometry = typename BaseType::Pointer( new Triangle3D6( NewGeometryId, rGeometry.Points() ) );
404  p_geometry->SetData(rGeometry.GetData());
405  return p_geometry;
406  }
407 
413  Matrix& PointsLocalCoordinates( Matrix& rResult ) const override
414  {
415  rResult.resize( 6, 2,false );
416  noalias( rResult ) = ZeroMatrix( 6, 2 );
417  rResult( 0, 0 ) = 0.0;
418  rResult( 0, 1 ) = 0.0;
419  rResult( 1, 0 ) = 1.0;
420  rResult( 1, 1 ) = 0.0;
421  rResult( 2, 0 ) = 0.0;
422  rResult( 2, 1 ) = 1.0;
423  rResult( 3, 0 ) = 0.5;
424  rResult( 3, 1 ) = 0.0;
425  rResult( 4, 0 ) = 0.5;
426  rResult( 4, 1 ) = 0.5;
427  rResult( 5, 0 ) = 0.0;
428  rResult( 5, 1 ) = 0.5;
429  return rResult;
430  }
431 
435 
451  double Length() const override
452  {
453  return std::sqrt(std::abs(this->DeterminantOfJacobian( PointType() ) ) );
454  }
455 
465  double Area() const override
466  {
467  const IntegrationMethod integration_method = msGeometryData.DefaultIntegrationMethod();
468  return IntegrationUtilities::ComputeDomainSize(*this, integration_method);
469  }
470 
471  // TODO: Code activated in June 2023
472  // /**
473  // * @brief This method calculates and returns the volume of this geometry.
474  // * @return Error, the volume of a 2D geometry is not defined
475  // * @see Length()
476  // * @see Area()
477  // * @see Volume()
478  // */
479  // double Volume() const override
480  // {
481  // KRATOS_ERROR << "Triangle3D6:: Method not well defined. Replace with DomainSize() instead" << std::endl;
482  // return 0.0;
483  // }
484 
493  double DomainSize() const override
494  {
495  return Area();
496  }
497 
506  bool IsInside(
507  const CoordinatesArrayType& rPoint,
508  CoordinatesArrayType& rResult,
509  const double Tolerance = std::numeric_limits<double>::epsilon()
510  )const override
511  {
512  this->PointLocalCoordinates( rResult, rPoint );
513 
514  if ( (rResult[0] >= (0.0-Tolerance)) && (rResult[0] <= (1.0+Tolerance)) )
515  {
516  if ( (rResult[1] >= (0.0-Tolerance)) && (rResult[1] <= (1.0+Tolerance)) )
517  {
518  if ( (rResult[0] + rResult[1]) <= (1.0+Tolerance) )
519  {
520  return true;
521  }
522  }
523  }
524 
525  return false;
526  }
527 
531 
544  double ShapeFunctionValue( IndexType ShapeFunctionIndex,
545  const CoordinatesArrayType& rPoint ) const override
546  {
547 
548  double thirdCoord = 1 - rPoint[0] - rPoint[1];
549 
550  switch ( ShapeFunctionIndex )
551  {
552  case 0:
553  return( thirdCoord*( 2*thirdCoord - 1 ) );
554  case 1:
555  return( rPoint[0]*( 2*rPoint[0] - 1 ) );
556  case 2:
557  return( rPoint[1]*( 2*rPoint[1] - 1 ) );
558  case 3:
559  return( 4*thirdCoord*rPoint[0] );
560  case 4:
561  return( 4*rPoint[0]*rPoint[1] );
562  case 5:
563  return( 4*rPoint[1]*thirdCoord );
564 
565  default:
566  KRATOS_ERROR << "Wrong index of shape function!" << *this << std::endl;
567  }
568 
569  return 0;
570  }
571 
575 
589  const CoordinatesArrayType& rPointGlobalCoordinates,
590  const double Tolerance = std::numeric_limits<double>::epsilon()
591  ) const override
592  {
593  const Point point(rPointGlobalCoordinates);
594  return GeometryUtils::PointDistanceToTriangle3D(this->GetPoint(0), this->GetPoint(1), this->GetPoint(2), this->GetPoint(3), this->GetPoint(4), this->GetPoint(5), point);
595  }
596 
600 
608  std::string Info() const override
609  {
610  return "2 dimensional triangle with six nodes in 3D space";
611  }
612 
619  void PrintInfo( std::ostream& rOStream ) const override
620  {
621  rOStream << "2 dimensional triangle with six nodes in 3D space";
622  }
623 
638  void PrintData( std::ostream& rOStream ) const override
639  {
640  // Base Geometry class PrintData call
641  BaseType::PrintData( rOStream );
642  std::cout << std::endl;
643 
644  // If the geometry has valid points, calculate and output its data
645  if (this->AllPointsAreValid()) {
646  Matrix jacobian;
647  this->Jacobian( jacobian, PointType() );
648  rOStream << " Jacobian in the origin\t : " << jacobian;
649  }
650  }
651 
655 
667  SizeType EdgesNumber() const override
668  {
669  return 3;
670  }
671 
681  {
683 
684  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 0 ), this->pGetPoint( 1 ), this->pGetPoint( 3 ) ) );
685  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 1 ), this->pGetPoint( 2 ), this->pGetPoint( 4 ) ) );
686  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 2 ), this->pGetPoint( 0 ), this->pGetPoint( 5 ) ) );
687  return edges;
688  }
689 
695  IntegrationMethod ThisMethod )
696  {
697  ShapeFunctionsGradientsType localGradients
698  = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
699  const int integration_points_number
700  = msGeometryData.IntegrationPointsNumber( ThisMethod );
701  ShapeFunctionsGradientsType Result( integration_points_number );
702 
703  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
704  {
705  Result[pnt] = localGradients[pnt];
706  }
707 
708  return Result;
709  }
710 
716  {
717  IntegrationMethod ThisMethod = msGeometryData.DefaultIntegrationMethod();
718  ShapeFunctionsGradientsType localGradients
719  = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
720  const int integration_points_number
721  = msGeometryData.IntegrationPointsNumber( ThisMethod );
722  ShapeFunctionsGradientsType Result( integration_points_number );
723 
724  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
725  {
726  Result[pnt] = localGradients[pnt];
727  }
728 
729  return Result;
730  }
731 
742  const CoordinatesArrayType& rPoint ) const override
743  {
744  rResult.resize( 6, 2,false );
745  double thirdCoord = 1 - rPoint[0] - rPoint[1];
746  double thirdCoord_DX = -1;
747  double thirdCoord_DY = -1;
748 
749  noalias( rResult ) = ZeroMatrix( 6, 2 );
750  rResult( 0, 0 ) = ( 4 * thirdCoord - 1 ) * thirdCoord_DX;
751  rResult( 0, 1 ) = ( 4 * thirdCoord - 1 ) * thirdCoord_DY;
752  rResult( 1, 0 ) = 4 * rPoint[0] - 1;
753  rResult( 1, 1 ) = 0;
754  rResult( 2, 0 ) = 0;
755  rResult( 2, 1 ) = 4 * rPoint[1] - 1;
756  rResult( 3, 0 ) = 4 * thirdCoord_DX * rPoint[0] + 4 * thirdCoord;
757  rResult( 3, 1 ) = 4 * thirdCoord_DY * rPoint[0];
758  rResult( 4, 0 ) = 4 * rPoint[1];
759  rResult( 4, 1 ) = 4 * rPoint[0];
760  rResult( 5, 0 ) = 4 * rPoint[1] * thirdCoord_DX;
761  rResult( 5, 1 ) = 4 * rPoint[1] * thirdCoord_DY + 4 * thirdCoord;
762  return rResult;
763  }
764 
765 
766 
777  {
778  rResult.resize( 6, 2 ,false);
779  double thirdCoord = 1 - rPoint[0] - rPoint[1];
780  double thirdCoord_DX = -1;
781  double thirdCoord_DY = -1;
782 
783  noalias( rResult ) = ZeroMatrix( 6, 2 );
784  rResult( 0, 0 ) = ( 4 * thirdCoord - 1 ) * thirdCoord_DX;
785  rResult( 0, 1 ) = ( 4 * thirdCoord - 1 ) * thirdCoord_DY;
786  rResult( 1, 0 ) = 4 * rPoint[0] - 1;
787  rResult( 1, 1 ) = 0;
788  rResult( 2, 0 ) = 0;
789  rResult( 2, 1 ) = 4 * rPoint[1] - 1;
790  rResult( 3, 0 ) = 4 * thirdCoord_DX * rPoint[0] + 4 * thirdCoord;
791  rResult( 3, 1 ) = 4 * thirdCoord_DY * rPoint[0];
792  rResult( 4, 0 ) = 4 * rPoint[1];
793  rResult( 4, 1 ) = 4 * rPoint[0];
794  rResult( 5, 0 ) = 4 * rPoint[1] * thirdCoord_DX;
795  rResult( 5, 1 ) = 4 * rPoint[1] * thirdCoord_DY + 4 * thirdCoord;
796  return rResult;
797  }
798 
806  {
807  if ( rResult.size() != this->PointsNumber() )
808  {
809  // KLUDGE: While there is a bug in
810  // ublas vector resize, I have to put this beside resizing!!
812  rResult.swap( temp );
813  }
814 
815  rResult[0].resize( 2, 2 ,false);
816  rResult[1].resize( 2, 2 ,false);
817  rResult[2].resize( 2, 2 ,false);
818  rResult[3].resize( 2, 2 ,false);
819  rResult[4].resize( 2, 2 ,false);
820  rResult[5].resize( 2, 2 ,false);
821 
822  rResult[0]( 0, 0 ) = 4.0;
823  rResult[0]( 0, 1 ) = 4.0;
824  rResult[0]( 1, 0 ) = 4.0;
825  rResult[0]( 1, 1 ) = 4.0;
826  rResult[1]( 0, 0 ) = 4.0;
827  rResult[1]( 0, 1 ) = 0.0;
828  rResult[1]( 1, 0 ) = 0.0;
829  rResult[1]( 1, 1 ) = 0.0;
830  rResult[2]( 0, 0 ) = 0.0;
831  rResult[2]( 0, 1 ) = 0.0;
832  rResult[2]( 1, 0 ) = 0.0;
833  rResult[2]( 1, 1 ) = 4.0;
834  rResult[3]( 0, 0 ) = -8.0;
835  rResult[3]( 0, 1 ) = -4.0;
836  rResult[3]( 1, 0 ) = -4.0;
837  rResult[3]( 1, 1 ) = 0.0;
838  rResult[4]( 0, 0 ) = 0.0;
839  rResult[4]( 0, 1 ) = 4.0;
840  rResult[4]( 1, 0 ) = 4.0;
841  rResult[4]( 1, 1 ) = 0.0;
842  rResult[5]( 0, 0 ) = 0.0;
843  rResult[5]( 0, 1 ) = -4.0;
844  rResult[5]( 1, 0 ) = -4.0;
845  rResult[5]( 1, 1 ) = -8.0;
846 
847  return rResult;
848  }
849 
857  {
858  if ( rResult.size() != this->PointsNumber() )
859  {
860  rResult.resize( this->PointsNumber() );
861  }
862 
863  for ( IndexType i = 0; i < rResult.size(); i++ )
864  {
865  rResult[i].resize( this->PointsNumber() );
866  }
867 
868  rResult[0][0].resize( 2, 2,false );
869 
870  rResult[0][1].resize( 2, 2,false );
871  rResult[1][0].resize( 2, 2,false );
872  rResult[1][1].resize( 2, 2 ,false);
873  rResult[2][0].resize( 2, 2 ,false);
874  rResult[2][1].resize( 2, 2 ,false);
875  rResult[3][0].resize( 2, 2 ,false);
876  rResult[3][1].resize( 2, 2 ,false);
877  rResult[4][0].resize( 2, 2 ,false);
878  rResult[4][1].resize( 2, 2 ,false);
879  rResult[5][0].resize( 2, 2 ,false);
880  rResult[5][1].resize( 2, 2 ,false);
881 
882 
883  for ( int i = 0; i < 6; i++ )
884  {
885  rResult[i][0]( 0, 0 ) = 0.0;
886  rResult[i][0]( 0, 1 ) = 0.0;
887  rResult[i][0]( 1, 0 ) = 0.0;
888  rResult[i][0]( 1, 1 ) = 0.0;
889  rResult[i][1]( 0, 0 ) = 0.0;
890  rResult[i][1]( 0, 1 ) = 0.0;
891  rResult[i][1]( 1, 0 ) = 0.0;
892  rResult[i][1]( 1, 1 ) = 0.0;
893  }
894 
895  return rResult;
896  }
897 
901 
903 
904 protected:
909 private:
912 
913  static const GeometryData msGeometryData;
914 
915  static const GeometryDimension msGeometryDimension;
916 
920 
921  friend class Serializer;
922 
923  void save( Serializer& rSerializer ) const override
924  {
926  }
927 
928  void load( Serializer& rSerializer ) override
929  {
931  }
932 
933  Triangle3D6(): BaseType( PointsArrayType(), &msGeometryData ) {}
934 
938 
939 
943 
944 
948 
959  static Matrix CalculateShapeFunctionsIntegrationPointsValues(
960  typename BaseType::IntegrationMethod ThisMethod )
961  {
962  IntegrationPointsContainerType all_integration_points =
963  AllIntegrationPoints();
964  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
965  //number of integration points
966  const int integration_points_number = integration_points.size();
967  //number of nodes in current geometry
968  const int points_number = 6;
969  //setting up return matrix
970  Matrix shape_function_values( integration_points_number, points_number );
971  //loop over all integration points
972 
973  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
974  {
975  double thirdCoord = 1 - integration_points[pnt].X() - integration_points[pnt].Y();
976 
977  shape_function_values( pnt, 0 ) = thirdCoord * ( 2 * thirdCoord - 1 ) ;
978  shape_function_values( pnt, 1 ) = integration_points[pnt].X() * ( 2 * integration_points[pnt].X() - 1 ) ;
979  shape_function_values( pnt, 2 ) = integration_points[pnt].Y() * ( 2 * integration_points[pnt].Y() - 1 ) ;
980  shape_function_values( pnt, 3 ) = 4 * thirdCoord * integration_points[pnt].X();
981  shape_function_values( pnt, 4 ) = 4 * integration_points[pnt].X() * integration_points[pnt].Y();
982  shape_function_values( pnt, 5 ) = 4 * integration_points[pnt].Y() * thirdCoord;
983 
984  }
985 
986  return shape_function_values;
987  }
988 
1002  CalculateShapeFunctionsIntegrationPointsLocalGradients(
1003  typename BaseType::IntegrationMethod ThisMethod )
1004  {
1005  IntegrationPointsContainerType all_integration_points =
1006  AllIntegrationPoints();
1007  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
1008  //number of integration points
1009  const int integration_points_number = integration_points.size();
1010  ShapeFunctionsGradientsType d_shape_f_values( integration_points_number );
1011  //initialising container
1012  //std::fill(d_shape_f_values.begin(), d_shape_f_values.end(), Matrix(4,2));
1013  //loop over all integration points
1014 
1015  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
1016  {
1017  Matrix result( 6, 2 );
1018  double thirdCoord = 1 - integration_points[pnt].X() - integration_points[pnt].Y();
1019  double thirdCoord_DX = -1;
1020  double thirdCoord_DY = -1;
1021 
1022  noalias( result ) = ZeroMatrix( 6, 2 );
1023  result( 0, 0 ) = ( 4 * thirdCoord - 1 ) * thirdCoord_DX;
1024  result( 0, 1 ) = ( 4 * thirdCoord - 1 ) * thirdCoord_DY;
1025  result( 1, 0 ) = 4 * integration_points[pnt].X() - 1;
1026  result( 1, 1 ) = 0;
1027  result( 2, 0 ) = 0;
1028  result( 2, 1 ) = 4 * integration_points[pnt].Y() - 1;
1029  result( 3, 0 ) = 4 * thirdCoord_DX * integration_points[pnt].X() + 4 * thirdCoord;
1030  result( 3, 1 ) = 4 * thirdCoord_DY * integration_points[pnt].X();
1031  result( 4, 0 ) = 4 * integration_points[pnt].Y();
1032  result( 4, 1 ) = 4 * integration_points[pnt].X();
1033  result( 5, 0 ) = 4 * integration_points[pnt].Y() * thirdCoord_DX;
1034  result( 5, 1 ) = 4 * integration_points[pnt].Y() * thirdCoord_DY + 4 * thirdCoord;
1035 
1036  d_shape_f_values[pnt] = result;
1037  }
1038 
1039  return d_shape_f_values;
1040  }
1041 
1045  static const IntegrationPointsContainerType AllIntegrationPoints()
1046  {
1047  IntegrationPointsContainerType integration_points =
1048  {
1049  {
1050  Quadrature<TriangleGaussLegendreIntegrationPoints1, 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1051  Quadrature<TriangleGaussLegendreIntegrationPoints2, 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1052  Quadrature<TriangleGaussLegendreIntegrationPoints3, 2, IntegrationPoint<3> >::GenerateIntegrationPoints()
1053  }
1054  };
1055  return integration_points;
1056  }
1057 
1061  static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
1062  {
1063  ShapeFunctionsValuesContainerType shape_functions_values =
1064  {
1065  {
1066  Triangle3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1068  Triangle3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1070  Triangle3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1072  }
1073  };
1074  return shape_functions_values;
1075  }
1076 
1081  AllShapeFunctionsLocalGradients()
1082  {
1083  ShapeFunctionsLocalGradientsContainerType shape_functions_local_gradients =
1084  {
1085  {
1086  Triangle3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_1 ),
1087  Triangle3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_2 ),
1088  Triangle3D6<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients( GeometryData::IntegrationMethod::GI_GAUSS_3 )
1089  }
1090  };
1091  return shape_functions_local_gradients;
1092  }
1093 
1097 
1098 
1102 
1103 
1107 
1108  template<class TOtherPointType> friend class Triangle3D6;
1109 
1113 
1115 }; // Class Geometry
1116 
1120 
1121 
1125 
1128 template<class TPointType> inline std::istream& operator >> (
1129  std::istream& rIStream,
1130  Triangle3D6<TPointType>& rThis );
1134 template<class TPointType> inline std::ostream& operator << (
1135  std::ostream& rOStream,
1136  const Triangle3D6<TPointType>& rThis )
1137 {
1138  rThis.PrintInfo( rOStream );
1139  rOStream << std::endl;
1140  rThis.PrintData( rOStream );
1141  return rOStream;
1142 }
1143 
1145 
1146 template<class TPointType> const
1147 GeometryData Triangle3D6<TPointType>::msGeometryData(
1150  Triangle3D6<TPointType>::AllIntegrationPoints(),
1151  Triangle3D6<TPointType>::AllShapeFunctionsValues(),
1152  AllShapeFunctionsLocalGradients()
1153 );
1154 
1155 template<class TPointType> const
1156 GeometryDimension Triangle3D6<TPointType>::msGeometryDimension(3, 2);
1157 
1158 }// 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
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
Short class definition.
Definition: integration_point.h:52
static double ComputeDomainSize(const TGeometryType &rGeometry)
This method calculates and returns the domain size of the geometry from any geometry in a generic man...
Definition: integration_utilities.h:63
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 three node 3D line geometry with quadratic shape functions.
Definition: line_3d_3.h:66
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 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
A six node 3D triangular geometry with quadratic shape functions.
Definition: triangle_3d_6.h:68
Triangle3D6(const PointsArrayType &ThisPoints)
Definition: triangle_3d_6.h:238
Triangle3D6(const std::string &rGeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: triangle_3d_6.h:255
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: triangle_3d_6.h:304
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: triangle_3d_6.h:299
BaseType::Pointer Create(const IndexType NewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: triangle_3d_6.h:370
Triangle3D6(Triangle3D6 const &rOther)
Definition: triangle_3d_6.h:272
GeometryData::IntegrationMethod IntegrationMethod
Definition: triangle_3d_6.h:92
Line3D3< TPointType > EdgeType
Definition: triangle_3d_6.h:82
ShapeFunctionsSecondDerivativesType & ShapeFunctionsSecondDerivatives(ShapeFunctionsSecondDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: triangle_3d_6.h:805
BaseType::PointsArrayType PointsArrayType
Definition: triangle_3d_6.h:124
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: triangle_3d_6.h:544
BaseType::IntegrationPointType IntegrationPointType
Definition: triangle_3d_6.h:135
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: triangle_3d_6.h:506
BaseType::IndexType IndexType
Definition: triangle_3d_6.h:111
TPointType PointType
Definition: triangle_3d_6.h:103
Geometry< TPointType > BaseType
Definition: triangle_3d_6.h:77
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: triangle_3d_6.h:741
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients(IntegrationMethod ThisMethod)
Definition: triangle_3d_6.h:694
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients()
Definition: triangle_3d_6.h:715
Triangle3D6 & operator=(Triangle3D6< TOtherPointType > const &rOther)
Definition: triangle_3d_6.h:342
void PrintData(std::ostream &rOStream) const override
Definition: triangle_3d_6.h:638
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: triangle_3d_6.h:157
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: triangle_3d_6.h:129
std::string Info() const override
Definition: triangle_3d_6.h:608
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: triangle_3d_6.h:144
double DomainSize() const override
This method calculates and returns length, area or volume of this geometry depending to it's dimensio...
Definition: triangle_3d_6.h:493
BaseType::Pointer Create(PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: triangle_3d_6.h:357
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: triangle_3d_6.h:667
Triangle3D6(Triangle3D6< TOtherPointType > const &rOther)
Definition: triangle_3d_6.h:289
BaseType::SizeType SizeType
Definition: triangle_3d_6.h:118
Triangle3D6 & operator=(const Triangle3D6 &rOther)
Definition: triangle_3d_6.h:324
Triangle3D6(const IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: triangle_3d_6.h:246
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: triangle_3d_6.h:588
BaseType::Pointer Create(const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: triangle_3d_6.h:383
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: triangle_3d_6.h:163
Triangle3D6(typename PointType::Pointer pFirstPoint, typename PointType::Pointer pSecondPoint, typename PointType::Pointer pThirdPoint, typename PointType::Pointer pFourthPoint, typename PointType::Pointer pFifthPoint, typename PointType::Pointer pSixthPoint)
Definition: triangle_3d_6.h:221
BaseType::JacobiansType JacobiansType
Definition: triangle_3d_6.h:170
virtual Matrix & ShapeFunctionsGradients(Matrix &rResult, CoordinatesArrayType &rPoint)
Definition: triangle_3d_6.h:776
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: triangle_3d_6.h:151
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: triangle_3d_6.h:398
double Area() const override
This method calculates and returns area or surface area of this geometry depending to it's dimension.
Definition: triangle_3d_6.h:465
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: triangle_3d_6.h:413
BaseType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: triangle_3d_6.h:185
BaseType::GeometriesArrayType GeometriesArrayType
Definition: triangle_3d_6.h:98
BaseType::NormalType NormalType
Definition: triangle_3d_6.h:198
double Length() const override
Definition: triangle_3d_6.h:451
~Triangle3D6() override
Definition: triangle_3d_6.h:297
KRATOS_CLASS_POINTER_DEFINITION(Triangle3D6)
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: triangle_3d_6.h:177
void PrintInfo(std::ostream &rOStream) const override
Definition: triangle_3d_6.h:619
ShapeFunctionsThirdDerivativesType & ShapeFunctionsThirdDerivatives(ShapeFunctionsThirdDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: triangle_3d_6.h:856
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: triangle_3d_6.h:680
BaseType::ShapeFunctionsThirdDerivativesType ShapeFunctionsThirdDerivativesType
Definition: triangle_3d_6.h:193
friend class Triangle3D6
Definition: triangle_3d_6.h:1108
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
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
const GeometryData Triangle3D6< TPointType >::msGeometryData & msGeometryDimension(), Triangle3D6< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
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
float temp
Definition: rotating_cone.py:85
integer i
Definition: TensorModule.f:17