28 #include "utilities/geometry_utilities.h"
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 )
222 :
BaseType( rThisPoints, &msGeometryData )
231 ) :
BaseType( GeometryId, rThisPoints, &msGeometryData )
238 const std::string& rGeometryName,
240 ) :
BaseType(rGeometryName, rThisPoints, &msGeometryData)
321 template<
class TOtherPo
intType>
345 return typename BaseType::Pointer(
new Prism3D6( NewGeometryId, rThisPoints ) );
359 auto p_geometry =
typename BaseType::Pointer(
new Prism3D6( NewGeometryId, rGeometry.
Points() ) );
360 p_geometry->SetData(rGeometry.
GetData());
387 return std::pow(
volume, 1.0/3.0)/3.0;
448 if ( rResult.size1() != 6 || rResult.size2() != 3 )
449 rResult.
resize( 6, 3 ,
false);
451 rResult( 0, 0 ) = 0.0;
452 rResult( 0, 1 ) = 0.0;
453 rResult( 0, 2 ) = 0.0;
455 rResult( 1, 0 ) = 1.0;
456 rResult( 1, 1 ) = 0.0;
457 rResult( 1, 2 ) = 0.0;
459 rResult( 2, 0 ) = 0.0;
460 rResult( 2, 1 ) = 1.0;
461 rResult( 2, 2 ) = 0.0;
463 rResult( 3, 0 ) = 0.0;
464 rResult( 3, 1 ) = 0.0;
465 rResult( 3, 2 ) = 1.0;
467 rResult( 4, 0 ) = 1.0;
468 rResult( 4, 1 ) = 0.0;
469 rResult( 4, 2 ) = 1.0;
471 rResult( 5, 0 ) = 0.0;
472 rResult( 5, 1 ) = 1.0;
473 rResult( 5, 2 ) = 1.0;
489 const double Tolerance = std::numeric_limits<double>::epsilon()
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) ) )
635 if(
IsInside(rLowPoint,local_coordinates))
658 switch ( ShapeFunctionIndex )
661 return( 1.0 -( rPoint[0] + rPoint[1] + rPoint[2] - ( rPoint[0]*rPoint[2] ) - ( rPoint[1]*rPoint[2] ) ) );
663 return( rPoint[0] - ( rPoint[0]*rPoint[2] ) );
665 return( rPoint[1] - ( rPoint[1]*rPoint[2] ) );
667 return( rPoint[2] - ( rPoint[0]*rPoint[2] ) - ( rPoint[1]*rPoint[2] ) );
669 return( rPoint[0]*rPoint[2] );
671 return( rPoint[1]*rPoint[2] );
673 KRATOS_ERROR <<
"Wrong index of shape function!" << *
this << std::endl;
692 if(rResult.size() != 6)
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];
717 if(rResult.size1() != this->PointsNumber() || rResult.size2() != this->LocalSpaceDimension())
720 CalculateShapeFunctionsLocalGradients(rResult, rPoint);
743 const double Tolerance = std::numeric_limits<double>::epsilon()
747 const Point point(rPointGlobalCoordinates);
751 if (this->
IsInside(rPointGlobalCoordinates, aux_coordinates, Tolerance)) {
756 std::array<double, 5> distances;
762 const auto min = std::min_element(distances.begin(), distances.end());
777 std::string
Info()
const override
779 return "3 dimensional prism with six nodes in 3D space";
791 rOStream <<
"3 dimensional prism with six nodes in 3D space";
807 std::cout << std::endl;
813 rOStream <<
" Jacobian in the origin\t : " << jacobian;
837 void save(
Serializer& rSerializer )
const override
865 static Matrix& CalculateShapeFunctionsLocalGradients(
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];
907 const int integration_points_number = integration_points.size();
909 const int points_number = 6;
911 Matrix shape_function_values( integration_points_number, points_number );
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() );
932 return shape_function_values;
948 CalculateShapeFunctionsIntegrationPointsLocalGradients(
952 AllIntegrationPoints();
955 const int integration_points_number = integration_points.size();
960 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
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;
984 return d_shape_f_values;
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(),
1014 return integration_points;
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(
1044 return shape_functions_values;
1051 AllShapeFunctionsLocalGradients()
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(
1078 return shape_functions_local_gradients;
1109 rOStream << std::endl;
1115 template<
class TPo
intType>
const
1116 GeometryData Prism3D6<TPointType>::msGeometryData(
1119 Prism3D6<TPointType>::AllIntegrationPoints(),
1120 Prism3D6<TPointType>::AllShapeFunctionsValues(),
1121 AllShapeFunctionsLocalGradients()
1124 template<
class TPo
intType>
const
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