215 typename PointType::Pointer pPoint2,
216 typename PointType::Pointer pPoint3,
217 typename PointType::Pointer pPoint4,
218 typename PointType::Pointer pPoint5,
219 typename PointType::Pointer pPoint6,
220 typename PointType::Pointer pPoint7,
221 typename PointType::Pointer pPoint8,
222 typename PointType::Pointer pPoint9,
223 typename PointType::Pointer pPoint10,
224 typename PointType::Pointer pPoint11,
225 typename PointType::Pointer pPoint12,
226 typename PointType::Pointer pPoint13,
227 typename PointType::Pointer pPoint14,
228 typename PointType::Pointer pPoint15,
229 typename PointType::Pointer pPoint16,
230 typename PointType::Pointer pPoint17,
231 typename PointType::Pointer pPoint18,
232 typename PointType::Pointer pPoint19,
233 typename PointType::Pointer pPoint20 )
259 :
BaseType( ThisPoints, &msGeometryData )
268 ) :
BaseType(GeometryId, rThisPoints, &msGeometryData)
275 const std::string& GeometryName,
277 ) :
BaseType( GeometryName, rThisPoints, &msGeometryData)
359 template<
class TOtherPo
intType>
382 return typename BaseType::Pointer(
new Hexahedra3D20( NewGeometryId, rThisPoints ) );
396 auto p_geometry =
typename BaseType::Pointer(
new Hexahedra3D20( NewGeometryId, rGeometry.
Points() ) );
397 p_geometry->SetData(rGeometry.
GetData());
484 const CoordinatesArrayType& rPoint,
486 const double Tolerance = std::numeric_limits<double>::epsilon()
491 if ( std::abs( rResult[0] ) <= (1.0 + Tolerance) )
493 if ( std::abs( rResult[1] ) <= (1.0 + Tolerance) )
495 if ( std::abs( rResult[2] ) <= (1.0 + Tolerance) )
701 switch ( ShapeFunctionIndex )
704 return -(( 1.0 + rPoint[0] )
705 *( 1.0 - rPoint[1] )*( 2.0
706 - rPoint[0] + rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
708 return -(( 1.0 + rPoint[0] )
709 *( 1.0 + rPoint[1] )*( 2.0
710 - rPoint[0] - rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
712 return -(( 1.0 + rPoint[0] )
713 *( 1.0 + rPoint[1] )*( 1.0 - rPoint[2] )*( 2.0
714 - rPoint[0] - rPoint[1] + rPoint[2] ) ) / 8.0;
716 return -(( 1.0 + rPoint[0] )
717 *( 1.0 - rPoint[1] )*( 1.0 - rPoint[2] )*( 2.0
718 - rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
720 return -(( 1.0 - rPoint[0] )
721 *( 1.0 - rPoint[1] )*( 2.0
722 + rPoint[0] + rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
724 return -(( 1.0 - rPoint[0] )
725 *( 1.0 + rPoint[1] )*( 2.0
726 + rPoint[0] - rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
728 return -(( 1.0 - rPoint[0] )*( 1.0
729 + rPoint[1] )*( 1.0 - rPoint[2] )*( 2.0
730 + rPoint[0] - rPoint[1] + rPoint[2] ) ) / 8.0;
732 return -(( 1.0 - rPoint[0] )*( 1.0 - rPoint[1] )
733 *( 1.0 - rPoint[2] )*( 2.0 + rPoint[0]
734 + rPoint[1] + rPoint[2] ) ) / 8.0;
736 return (( 1.0 + rPoint[0] )
737 *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
739 return (( 1.0 + rPoint[0] )*( 1.0 + rPoint[1] )
740 *( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
742 return (( 1.0 + rPoint[0] )
743 *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0 ;
745 return (( 1.0 + rPoint[0] )*( 1.0 - rPoint[1] )
746 *( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
748 return (( 1.0 -rPoint[0]
749 *rPoint[0] )*( 1.0 - rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
751 return (( 1.0 -rPoint[0]
752 *rPoint[0] )*( 1.0 + rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
754 return (( 1.0 -rPoint[0]
755 *rPoint[0] )*( 1.0 + rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0 ;
757 return (( 1.0 -rPoint[0]*rPoint[0] )
758 *( 1.0 - rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0;
760 return (( 1.0 -rPoint[0] )
761 *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
763 return (( 1.0 -rPoint[0] )*( 1.0 + rPoint[1] )
764 *( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
766 return (( 1.0 -rPoint[0] )
767 *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0 ;
769 return (( 1.0 -rPoint[0] )
770 *( 1.0 - rPoint[1] )*( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
773 KRATOS_ERROR <<
"Wrong index of shape function!" << *
this << std::endl;
792 if(rResult.size() != 20) rResult.
resize(20,
false);
793 rResult[0] = -(( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] )*( 2.0
794 - rCoordinates[0] + rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
795 rResult[1] = -(( 1.0 + rCoordinates[0] )
796 *( 1.0 + rCoordinates[1] )*( 2.0
797 - rCoordinates[0] - rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
798 rResult[2] = -(( 1.0 + rCoordinates[0] )
799 *( 1.0 + rCoordinates[1] )*( 1.0 - rCoordinates[2] )*( 2.0
800 - rCoordinates[0] - rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
801 rResult[3] = -(( 1.0 + rCoordinates[0] )
802 *( 1.0 - rCoordinates[1] )*( 1.0 - rCoordinates[2] )*( 2.0
803 - rCoordinates[0] + rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
804 rResult[4] = -(( 1.0 - rCoordinates[0] )
805 *( 1.0 - rCoordinates[1] )*( 2.0
806 + rCoordinates[0] + rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
807 rResult[5] = -(( 1.0 - rCoordinates[0] )
808 *( 1.0 + rCoordinates[1] )*( 2.0
809 + rCoordinates[0] - rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
810 rResult[6] = -(( 1.0 - rCoordinates[0] )*( 1.0
811 + rCoordinates[1] )*( 1.0 - rCoordinates[2] )*( 2.0
812 + rCoordinates[0] - rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
813 rResult[7] = -(( 1.0 - rCoordinates[0] )*( 1.0 - rCoordinates[1] )
814 *( 1.0 - rCoordinates[2] )*( 2.0 + rCoordinates[0]
815 + rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
816 rResult[8] = (( 1.0 + rCoordinates[0] )
817 *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
818 rResult[9] = (( 1.0 + rCoordinates[0] )*( 1.0 + rCoordinates[1] )
819 *( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
820 rResult[10] = (( 1.0 + rCoordinates[0] )
821 *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0 ;
822 rResult[11] = (( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] )
823 *( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
824 rResult[12] = (( 1.0 -rCoordinates[0]
825 *rCoordinates[0] )*( 1.0 - rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
826 rResult[13] = (( 1.0 -rCoordinates[0]
827 *rCoordinates[0] )*( 1.0 + rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
828 rResult[14] = (( 1.0 -rCoordinates[0]
829 *rCoordinates[0] )*( 1.0 + rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0 ;
830 rResult[15] = (( 1.0 -rCoordinates[0]*rCoordinates[0] )
831 *( 1.0 - rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0;
832 rResult[16] = (( 1.0 -rCoordinates[0] )
833 *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
834 rResult[17] = (( 1.0 -rCoordinates[0] )*( 1.0 + rCoordinates[1] )
835 *( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
836 rResult[18] = (( 1.0 -rCoordinates[0] )
837 *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0 ;
838 rResult[19] = (( 1.0 -rCoordinates[0] )
839 *( 1.0 - rCoordinates[1] )*( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
853 std::string
Info()
const override
855 return "3 dimensional hexahedra with 20 nodes and quadratic shape functions in 3D space";
867 rOStream <<
"3 dimensional hexahedra with 20 nodes and quadratic shape functions in 3D space";
883 std::cout << std::endl;
889 rOStream <<
" Jacobian in the origin\t : " << jacobian;
908 if ( result.size1() != 20 || result.size2() != 3 )
909 result.
resize( 20, 3,
false );
911 result( 0, 0 ) = (( -1.0 + rPoint[1] ) * ( 1.0 - 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( 1.0
912 + rPoint[2] ) ) / 8.0;
914 result( 0, 1 ) = -(( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[2] ) * ( -1.0 + rPoint[0] - 2.0
915 * rPoint[1] + rPoint[2] ) ) / 8.0;
917 result( 0, 2 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] - rPoint[1] + 2.0
918 * rPoint[2] ) ) / 8.0;
920 result( 1, 0 ) = (( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) * ( -1.0 + 2.0
921 * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
923 result( 1, 1 ) = (( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[2] ) * ( -1.0 + rPoint[0] + 2.0
924 * rPoint[1] + rPoint[2] ) ) / 8.0;
926 result( 1, 2 ) = (( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] + rPoint[1] + 2.0
927 * rPoint[2] ) ) / 8.0;
929 result( 2, 0 ) = -(( 1.0 + rPoint[1] ) * ( -1.0 + 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( -1.0
930 + rPoint[2] ) ) / 8.0;
932 result( 2, 1 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[0] + 2.0 * rPoint[1] - rPoint[2] ) * ( -1.0
933 + rPoint[2] ) ) / 8.0;
935 result( 2, 2 ) = -(( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] + rPoint[1] - 2.0
936 * rPoint[2] ) ) / 8.0;
938 result( 3, 0 ) = -(( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) * ( 1.0 - 2.0
939 * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
941 result( 3, 1 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[0] - 2.0 * rPoint[1] - rPoint[2] ) * ( -1.0
942 + rPoint[2] ) ) / 8.0;
944 result( 3, 2 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] - rPoint[1] - 2.0
945 * rPoint[2] ) ) / 8.0;
947 result( 4, 0 ) = -(( -1.0 + rPoint[1] ) * ( 1.0 + 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( 1.0
948 + rPoint[2] ) ) / 8.0;
950 result( 4, 1 ) = -(( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[0] + 2.0 * rPoint[1] - rPoint[2] ) * ( 1.0
951 + rPoint[2] ) ) / 8.0;
953 result( 4, 2 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] + rPoint[1] - 2.0
954 * rPoint[2] ) ) / 8.0;
956 result( 5, 0 ) = -(( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) * ( -1.0 - 2.0
957 * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
959 result( 5, 1 ) = (( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[0] - 2.0 * rPoint[1] - rPoint[2] ) * ( 1.0
960 + rPoint[2] ) ) / 8.0;
962 result( 5, 2 ) = (( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] - rPoint[1] - 2.0
963 * rPoint[2] ) ) / 8.0;
965 result( 6, 0 ) = (( 1.0 + rPoint[1] ) * ( -1.0 - 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( -1.0
966 + rPoint[2] ) ) / 8.0;
968 result( 6, 1 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] ) * ( 1.0 + rPoint[0] - 2.0
969 * rPoint[1] + rPoint[2] ) ) / 8.0;
971 result( 6, 2 ) = -(( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] - rPoint[1] + 2.0
972 * rPoint[2] ) ) / 8.0;
974 result( 7, 0 ) = (( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) * ( 1.0 + 2.0
975 * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
977 result( 7, 1 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] ) * ( 1.0 + rPoint[0] + 2.0
978 * rPoint[1] + rPoint[2] ) ) / 8.0;
980 result( 7, 2 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] + rPoint[1] + 2.0
981 * rPoint[2] ) ) / 8.0;
983 result( 8, 0 ) = -(( -1.0 + rPoint[1] * rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
985 result( 8, 1 ) = -(( 1.0 + rPoint[0] ) * rPoint[1] * ( 1.0 + rPoint[2] ) ) / 2.0;
987 result( 8, 2 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
989 result( 9, 0 ) = -(( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
991 result( 9, 1 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
993 result( 9, 2 ) = -(( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
995 result( 10, 0 ) = (( -1.0 + rPoint[1] * rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
997 result( 10, 1 ) = (( 1.0 + rPoint[0] ) * rPoint[1] * ( -1.0 + rPoint[2] ) ) / 2.0;
999 result( 10, 2 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
1001 result( 11, 0 ) = (( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1003 result( 11, 1 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1005 result( 11, 2 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
1007 result( 12, 0 ) = ( rPoint[0] * ( -1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 2.0;
1009 result( 12, 1 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
1011 result( 12, 2 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[1] ) ) / 4.0;
1013 result( 13, 0 ) = -( rPoint[0] * ( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 2.0;
1015 result( 13, 1 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
1017 result( 13, 2 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[1] ) ) / 4.0;
1019 result( 14, 0 ) = ( rPoint[0] * ( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 2.0;
1021 result( 14, 1 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
1023 result( 14, 2 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[1] ) ) / 4.0;
1025 result( 15, 0 ) = -( rPoint[0] * ( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 2.0;
1027 result( 15, 1 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
1029 result( 15, 2 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[1] ) ) / 4.0;
1031 result( 16, 0 ) = (( -1.0 + rPoint[1] * rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
1033 result( 16, 1 ) = (( -1.0 + rPoint[0] ) * rPoint[1] * ( 1.0 + rPoint[2] ) ) / 2.0;
1035 result( 16, 2 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
1037 result( 17, 0 ) = (( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1039 result( 17, 1 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1041 result( 17, 2 ) = (( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
1043 result( 18, 0 ) = -(( -1.0 + rPoint[1] * rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
1045 result( 18, 1 ) = -(( -1.0 + rPoint[0] ) * rPoint[1] * ( -1.0 + rPoint[2] ) ) / 2.0;
1047 result( 18, 2 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
1049 result( 19, 0 ) = -(( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1051 result( 19, 1 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1053 result( 19, 2 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
1083 void save(
Serializer& rSerializer )
const override
1113 static Matrix CalculateShapeFunctionsIntegrationPointsValues(
1117 AllIntegrationPoints();
1120 const int integration_points_number = integration_points.size();
1122 const int points_number = 20;
1124 Matrix shape_function_values( integration_points_number, points_number );
1127 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1129 shape_function_values(pnt, 0 ) =
1130 -(( 1.0 + integration_points[pnt].X() )
1131 * ( 1.0 - integration_points[pnt].
Y() ) * ( 2.0
1132 - integration_points[pnt].
X() + integration_points[pnt].Y()
1133 - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].
Z() ) ) / 8.0;
1134 shape_function_values(pnt, 1 ) =
1135 -(( 1.0 + integration_points[pnt].X() )
1136 * ( 1.0 + integration_points[pnt].
Y() ) * ( 2.0
1137 - integration_points[pnt].
X() - integration_points[pnt].Y()
1138 - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].
Z() ) ) / 8.0;
1139 shape_function_values(pnt, 2 ) =
1140 -(( 1.0 + integration_points[pnt].X() )
1141 * ( 1.0 + integration_points[pnt].
Y() ) * ( 1.0 - integration_points[pnt].
Z() ) * ( 2.0
1142 - integration_points[pnt].
X() - integration_points[pnt].Y()
1143 + integration_points[pnt].Z() ) ) / 8.0;
1144 shape_function_values(pnt, 3 ) =
1145 -(( 1.0 + integration_points[pnt].X() )
1146 * ( 1.0 - integration_points[pnt].
Y() ) * ( 1.0 - integration_points[pnt].
Z() )
1147 * ( 2.0 - integration_points[pnt].
X()
1148 + integration_points[pnt].Y() + integration_points[pnt].Z() ) ) / 8.0;
1149 shape_function_values(pnt, 4 ) =
1150 -(( 1.0 - integration_points[pnt].X() )
1151 * ( 1.0 - integration_points[pnt].
Y() ) * ( 2.0
1152 + integration_points[pnt].
X() + integration_points[pnt].Y()
1153 - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].
Z() ) ) / 8.0;
1154 shape_function_values(pnt, 5 ) =
1155 -(( 1.0 - integration_points[pnt].X() )
1156 * ( 1.0 + integration_points[pnt].
Y() ) * ( 2.0
1157 + integration_points[pnt].
X() - integration_points[pnt].Y()
1158 - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].
Z() ) ) / 8.0;
1159 shape_function_values(pnt, 6 ) =
1160 -(( 1.0 - integration_points[pnt].X() ) * ( 1.0
1161 + integration_points[pnt].
Y() ) * ( 1.0 - integration_points[pnt].
Z() ) * ( 2.0
1162 + integration_points[pnt].
X() - integration_points[pnt].Y()
1163 + integration_points[pnt].Z() ) ) / 8.0;
1164 shape_function_values(pnt, 7 ) =
1165 -(( 1.0 - integration_points[pnt].X() ) * ( 1.0 - integration_points[pnt].
Y() )
1166 * ( 1.0 - integration_points[pnt].
Z() ) * ( 2.0 + integration_points[pnt].
X()
1167 + integration_points[pnt].Y() + integration_points[pnt].Z() ) ) / 8.0;
1168 shape_function_values(pnt, 8 ) =
1169 (( 1.0 + integration_points[pnt].X() )
1170 * ( 1.0 - integration_points[pnt].
Y() * integration_points[pnt].Y() )
1171 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0 ;
1172 shape_function_values(pnt, 9 ) =
1173 (( 1.0 + integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].
Y() )
1174 * ( 1.0 - integration_points[pnt].
Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1175 shape_function_values(pnt, 10 ) =
1176 (( 1.0 + integration_points[pnt].X() )
1177 * ( 1.0 - integration_points[pnt].
Y() * integration_points[pnt].Y() )
1178 * ( 1.0 - integration_points[pnt].
Z() ) ) / 4.0 ;
1179 shape_function_values(pnt, 11 ) =
1180 (( 1.0 + integration_points[pnt].X() ) * ( 1.0 - integration_points[pnt].
Y() )
1181 * ( 1.0 - integration_points[pnt].
Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1182 shape_function_values(pnt, 12 ) =
1183 (( 1.0 - integration_points[pnt].X()
1184 * integration_points[pnt].X() ) * ( 1.0 - integration_points[pnt].
Y() )
1185 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0 ;
1186 shape_function_values(pnt, 13 ) =
1187 (( 1.0 - integration_points[pnt].X()
1188 * integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].
Y() )
1189 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0 ;
1190 shape_function_values(pnt, 14 ) =
1191 (( 1.0 - integration_points[pnt].X()
1192 * integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].
Y() )
1193 * ( 1.0 - integration_points[pnt].
Z() ) ) / 4.0 ;
1194 shape_function_values(pnt, 15 ) =
1195 (( 1.0 - integration_points[pnt].X() * integration_points[pnt].X() )
1196 * ( 1.0 - integration_points[pnt].
Y() ) * ( 1.0 - integration_points[pnt].
Z() ) ) / 4.0;
1197 shape_function_values(pnt, 16 ) =
1198 (( 1.0 - integration_points[pnt].X() )
1199 * ( 1.0 - integration_points[pnt].
Y() * integration_points[pnt].Y() )
1200 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0 ;
1201 shape_function_values(pnt, 17 ) =
1202 (( 1.0 - integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].
Y() )
1203 * ( 1.0 - integration_points[pnt].
Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1204 shape_function_values(pnt, 18 ) =
1205 (( 1.0 - integration_points[pnt].X() )
1206 * ( 1.0 - integration_points[pnt].
Y() * integration_points[pnt].Y() )
1207 * ( 1.0 - integration_points[pnt].
Z() ) ) / 4.0 ;
1208 shape_function_values(pnt, 19 ) =
1209 (( 1.0 - integration_points[pnt].X() )
1210 * ( 1.0 - integration_points[pnt].
Y() ) * ( 1.0
1211 - integration_points[pnt].
Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1214 return shape_function_values;
1231 CalculateShapeFunctionsIntegrationPointsLocalGradients(
1235 AllIntegrationPoints();
1238 const int integration_points_number = integration_points.size();
1244 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1248 result( 0, 0 ) = (( -1.0 + integration_points[pnt].Y() )
1249 * ( 1.0 - 2.0 * integration_points[pnt].
X()
1250 + integration_points[pnt].Y()
1251 - integration_points[pnt].Z() ) * ( 1.0
1252 + integration_points[pnt].
Z() ) ) / 8.0;
1253 result( 0, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1254 * ( 1.0 + integration_points[pnt].
Z() ) * ( -1.0
1255 + integration_points[pnt].
X() - 2.0
1256 * integration_points[pnt].Y()
1257 + integration_points[pnt].Z() ) ) / 8.0;
1258 result( 0, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1259 * ( -1.0 + integration_points[pnt].
Y() ) * ( -1.0
1260 + integration_points[pnt].
X()
1261 - integration_points[pnt].Y() + 2.0
1262 * integration_points[pnt].Z() ) ) / 8.0;
1264 result( 1, 0 ) = (( 1.0 + integration_points[pnt].Y() )
1265 * ( 1.0 + integration_points[pnt].
Z() )
1266 * ( -1.0 + 2.0 * integration_points[pnt].
X()
1267 + integration_points[pnt].Y()
1268 + integration_points[pnt].Z() ) ) / 8.0;
1269 result( 1, 1 ) = (( 1.0 + integration_points[pnt].X() )
1270 * ( 1.0 + integration_points[pnt].
Z() ) * ( -1.0
1271 + integration_points[pnt].
X() + 2.0
1272 * integration_points[pnt].Y()
1273 + integration_points[pnt].Z() ) ) / 8.0;
1274 result( 1, 2 ) = (( 1.0 + integration_points[pnt].X() )
1275 * ( 1.0 + integration_points[pnt].
Y() ) * ( -1.0
1276 + integration_points[pnt].
X()
1277 + integration_points[pnt].Y() + 2.0
1278 * integration_points[pnt].Z() ) ) / 8.0;
1280 result( 2, 0 ) = -(( 1.0 + integration_points[pnt].Y() )
1281 * ( -1.0 + 2.0 * integration_points[pnt].
X()
1282 + integration_points[pnt].Y()
1283 - integration_points[pnt].Z() ) * ( -1.0 + integration_points[pnt].
Z() ) ) / 8.0;
1284 result( 2, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1285 * ( -1.0 + integration_points[pnt].
X() + 2.0
1286 * integration_points[pnt].Y()
1287 - integration_points[pnt].Z() ) * ( -1.0
1288 + integration_points[pnt].
Z() ) ) / 8.0;
1289 result( 2, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1290 * ( 1.0 + integration_points[pnt].
Y() )
1291 * ( -1.0 + integration_points[pnt].
X()
1292 + integration_points[pnt].Y() - 2.0
1293 * integration_points[pnt].Z() ) ) / 8.0;
1295 result( 3, 0 ) = -(( -1.0 + integration_points[pnt].Y() )
1296 * ( -1.0 + integration_points[pnt].
Z() ) * ( 1.0 - 2.0 * integration_points[pnt].
X()
1297 + integration_points[pnt].Y()
1298 + integration_points[pnt].Z() ) ) / 8.0;
1299 result( 3, 1 ) = (( 1.0 + integration_points[pnt].X() )
1300 * ( -1.0 + integration_points[pnt].
X() - 2.0
1301 * integration_points[pnt].Y()
1302 - integration_points[pnt].Z() ) * ( -1.0
1303 + integration_points[pnt].
Z() ) ) / 8.0;
1304 result( 3, 2 ) = (( 1.0 + integration_points[pnt].X() )
1305 * ( -1.0 + integration_points[pnt].
Y() )
1306 * ( -1.0 + integration_points[pnt].
X()
1307 - integration_points[pnt].Y() - 2.0
1308 * integration_points[pnt].Z() ) ) / 8.0;
1310 result( 4, 0 ) = -(( -1.0 + integration_points[pnt].Y() )
1311 * ( 1.0 + 2.0 * integration_points[pnt].
X()
1312 + integration_points[pnt].Y()
1313 - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].
Z() ) ) / 8.0;
1314 result( 4, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1315 * ( 1.0 + integration_points[pnt].
X() + 2.0
1316 * integration_points[pnt].Y()
1317 - integration_points[pnt].Z() ) * ( 1.0
1318 + integration_points[pnt].
Z() ) ) / 8.0;
1319 result( 4, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1320 * ( -1.0 + integration_points[pnt].
Y() ) * ( 1.0
1321 + integration_points[pnt].
X()
1322 + integration_points[pnt].Y() - 2.0
1323 * integration_points[pnt].Z() ) ) / 8.0;
1325 result( 5, 0 ) = -(( 1.0 + integration_points[pnt].Y() )
1326 * ( 1.0 + integration_points[pnt].
Z() ) * ( -1.0 - 2.0
1327 * integration_points[pnt].
X()
1328 + integration_points[pnt].Y()
1329 + integration_points[pnt].Z() ) ) / 8.0;
1330 result( 5, 1 ) = (( -1.0 + integration_points[pnt].X() )
1331 * ( 1.0 + integration_points[pnt].
X() - 2.0
1332 * integration_points[pnt].Y()
1333 - integration_points[pnt].Z() ) * ( 1.0
1334 + integration_points[pnt].
Z() ) ) / 8.0;
1335 result( 5, 2 ) = (( -1.0 + integration_points[pnt].X() )
1336 * ( 1.0 + integration_points[pnt].
Y() ) * ( 1.0
1337 + integration_points[pnt].
X()
1338 - integration_points[pnt].Y() - 2.0
1339 * integration_points[pnt].Z() ) ) / 8.0;
1341 result( 6, 0 ) = (( 1.0 + integration_points[pnt].Y() )
1342 * ( -1.0 - 2.0 * integration_points[pnt].
X()
1343 + integration_points[pnt].Y()
1344 - integration_points[pnt].Z() ) * ( -1.0
1345 + integration_points[pnt].
Z() ) ) / 8.0;
1346 result( 6, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1347 * ( -1.0 + integration_points[pnt].
Z() )
1348 * ( 1.0 + integration_points[pnt].
X()
1349 - 2.0 * integration_points[pnt].Y()
1350 + integration_points[pnt].Z() ) ) / 8.0;
1351 result( 6, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1352 * ( 1.0 + integration_points[pnt].
Y() ) * ( 1.0
1353 + integration_points[pnt].
X()
1354 - integration_points[pnt].Y() + 2.0
1355 * integration_points[pnt].Z() ) ) / 8.0;
1357 result( 7, 0 ) = (( -1.0 + integration_points[pnt].Y() )
1358 * ( -1.0 + integration_points[pnt].
Z() ) * ( 1.0 + 2.0
1359 * integration_points[pnt].
X() + integration_points[pnt].Y()
1360 + integration_points[pnt].Z() ) ) / 8.0;
1361 result( 7, 1 ) = (( -1.0 + integration_points[pnt].X() )
1362 * ( -1.0 + integration_points[pnt].
Z() )
1363 * ( 1.0 + integration_points[pnt].
X() + 2.0 * integration_points[pnt].Y()
1364 + integration_points[pnt].Z() ) ) / 8.0;
1365 result( 7, 2 ) = (( -1.0 + integration_points[pnt].X() )
1366 * ( -1.0 + integration_points[pnt].
Y() ) * ( 1.0 + integration_points[pnt].
X()
1367 + integration_points[pnt].Y() + 2.0
1368 * integration_points[pnt].Z() ) ) / 8.0;
1370 result( 8, 0 ) = -(( -1.0 + integration_points[pnt].Y()
1371 * integration_points[pnt].Y() )
1372 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1373 result( 8, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1374 * integration_points[pnt].
Y()
1375 * ( 1.0 + integration_points[pnt].Z() ) ) / 2.0;
1376 result( 8, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1377 * ( -1.0 + integration_points[pnt].
Y()
1378 * integration_points[pnt].Y() ) ) / 4.0;
1380 result( 9, 0 ) = -(( 1.0 + integration_points[pnt].Y() )
1381 * ( -1.0 + integration_points[pnt].
Z()
1382 * integration_points[pnt].Z() ) ) / 4.0;
1383 result( 9, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1384 * ( -1.0 + integration_points[pnt].
Z()
1385 * integration_points[pnt].Z() ) ) / 4.0;
1386 result( 9, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1387 * ( 1.0 + integration_points[pnt].
Y() )
1388 * integration_points[pnt].
Z() ) / 2.0;
1390 result( 10, 0 ) = (( -1.0 + integration_points[pnt].Y()
1391 * integration_points[pnt].Y() )
1392 * ( -1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1393 result( 10, 1 ) = (( 1.0 + integration_points[pnt].X() )
1394 * integration_points[pnt].
Y()
1395 * ( -1.0 + integration_points[pnt].Z() ) ) / 2.0;
1396 result( 10, 2 ) = (( 1.0 + integration_points[pnt].X() )
1397 * ( -1.0 + integration_points[pnt].
Y()
1398 * integration_points[pnt].Y() ) ) / 4.0;
1400 result( 11, 0 ) = (( -1.0 + integration_points[pnt].Y() )
1401 * ( -1.0 + integration_points[pnt].
Z()
1402 * integration_points[pnt].Z() ) ) / 4.0;
1403 result( 11, 1 ) = (( 1.0 + integration_points[pnt].X() )
1404 * ( -1.0 + integration_points[pnt].
Z()
1405 * integration_points[pnt].Z() ) ) / 4.0;
1406 result( 11, 2 ) = (( 1.0 + integration_points[pnt].X() )
1407 * ( -1.0 + integration_points[pnt].
Y() )
1408 * integration_points[pnt].
Z() ) / 2.0;
1410 result( 12, 0 ) = ( integration_points[pnt].X()
1411 * ( -1.0 + integration_points[pnt].Y() )
1412 * ( 1.0 + integration_points[pnt].
Z() ) ) / 2.0;
1413 result( 12, 1 ) = (( -1.0 + integration_points[pnt].X()
1414 * integration_points[pnt].X() )
1415 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1416 result( 12, 2 ) = (( -1.0 + integration_points[pnt].X()
1417 * integration_points[pnt].X() )
1418 * ( -1.0 + integration_points[pnt].
Y() ) ) / 4.0;
1420 result( 13, 0 ) = -( integration_points[pnt].X()
1421 * ( 1.0 + integration_points[pnt].Y() )
1422 * ( 1.0 + integration_points[pnt].
Z() ) ) / 2.0;
1423 result( 13, 1 ) = -(( -1.0 + integration_points[pnt].X()
1424 * integration_points[pnt].X() )
1425 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1426 result( 13, 2 ) = -(( -1.0 + integration_points[pnt].X()
1427 * integration_points[pnt].X() )
1428 * ( 1.0 + integration_points[pnt].
Y() ) ) / 4.0;
1430 result( 14, 0 ) = ( integration_points[pnt].X()
1431 * ( 1.0 + integration_points[pnt].Y() )
1432 * ( -1.0 + integration_points[pnt].
Z() ) ) / 2.0;
1433 result( 14, 1 ) = (( -1.0 + integration_points[pnt].X()
1434 * integration_points[pnt].X() )
1435 * ( -1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1436 result( 14, 2 ) = (( -1.0 + integration_points[pnt].X()
1437 * integration_points[pnt].X() )
1438 * ( 1.0 + integration_points[pnt].
Y() ) ) / 4.0;
1440 result( 15, 0 ) = -( integration_points[pnt].X() * ( -1.0
1441 + integration_points[pnt].Y() ) * ( -1.0
1442 + integration_points[pnt].
Z() ) ) / 2.0;
1443 result( 15, 1 ) = -(( -1.0 + integration_points[pnt].X()
1444 * integration_points[pnt].X() ) * ( -1.0
1445 + integration_points[pnt].
Z() ) ) / 4.0;
1446 result( 15, 2 ) = -(( -1.0 + integration_points[pnt].X()
1447 * integration_points[pnt].X() )
1448 * ( -1.0 + integration_points[pnt].
Y() ) ) / 4.0;
1450 result( 16, 0 ) = (( -1.0 + integration_points[pnt].Y()
1451 * integration_points[pnt].Y() )
1452 * ( 1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1453 result( 16, 1 ) = (( -1.0 + integration_points[pnt].X() )
1454 * integration_points[pnt].
Y()
1455 * ( 1.0 + integration_points[pnt].Z() ) ) / 2.0;
1456 result( 16, 2 ) = (( -1.0 + integration_points[pnt].X() )
1457 * ( -1.0 + integration_points[pnt].
Y()
1458 * integration_points[pnt].Y() ) ) / 4.0;
1460 result( 17, 0 ) = (( 1.0 + integration_points[pnt].Y() )
1461 * ( -1.0 + integration_points[pnt].
Z()
1462 * integration_points[pnt].Z() ) ) / 4.0;
1463 result( 17, 1 ) = (( -1.0 + integration_points[pnt].X() )
1464 * ( -1.0 + integration_points[pnt].
Z()
1465 * integration_points[pnt].Z() ) ) / 4.0;
1466 result( 17, 2 ) = (( -1.0 + integration_points[pnt].X() )
1467 * ( 1.0 + integration_points[pnt].
Y() )
1468 * integration_points[pnt].
Z() ) / 2.0;
1470 result( 18, 0 ) = -(( -1.0 + integration_points[pnt].Y()
1471 * integration_points[pnt].Y() )
1472 * ( -1.0 + integration_points[pnt].
Z() ) ) / 4.0;
1473 result( 18, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1474 * integration_points[pnt].
Y()
1475 * ( -1.0 + integration_points[pnt].Z() ) ) / 2.0;
1476 result( 18, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1477 * ( -1.0 + integration_points[pnt].
Y()
1478 * integration_points[pnt].Y() ) ) / 4.0;
1480 result( 19, 0 ) = -(( -1.0 + integration_points[pnt].Y() )
1481 * ( -1.0 + integration_points[pnt].
Z()
1482 * integration_points[pnt].Z() ) ) / 4.0;
1483 result( 19, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1484 * ( -1.0 + integration_points[pnt].
Z()
1485 * integration_points[pnt].Z() ) ) / 4.0;
1486 result( 19, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1487 * ( -1.0 + integration_points[pnt].
Y() )
1488 * integration_points[pnt].
Z() ) / 2.0;
1490 d_shape_f_values[pnt] = result;
1493 return d_shape_f_values;
1504 Quadrature < HexahedronGaussLegendreIntegrationPoints1,
1505 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1506 Quadrature < HexahedronGaussLegendreIntegrationPoints2,
1507 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1508 Quadrature < HexahedronGaussLegendreIntegrationPoints3,
1509 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1510 Quadrature < HexahedronGaussLegendreIntegrationPoints4,
1511 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1512 Quadrature < HexahedronGaussLegendreIntegrationPoints5,
1513 3, IntegrationPoint<3> >::GenerateIntegrationPoints()
1516 return integration_points;
1527 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1529 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1531 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1533 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1535 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1539 return shape_functions_values;
1546 AllShapeFunctionsLocalGradients()
1551 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1553 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1555 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1557 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1559 Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1563 return shape_functions_local_gradients;
1598 rOStream << std::endl;
1604 template<
class TPo
intType>
const
1605 GeometryData Hexahedra3D20<TPointType>::msGeometryData(
1608 Hexahedra3D20<TPointType>::AllIntegrationPoints(),
1609 Hexahedra3D20<TPointType>::AllShapeFunctionsValues(),
1610 AllShapeFunctionsLocalGradients()
1613 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
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
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
GeometryData::IntegrationMethod IntegrationMethod
Definition: geometry.h:122
A twenty node hexahedra geometry with serendipity shape functions.
Definition: hexahedra_3d_20.h:51
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: hexahedra_3d_20.h:698
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: hexahedra_3d_20.h:123
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: hexahedra_3d_20.h:168
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: hexahedra_3d_20.h:525
Hexahedra3D20 & operator=(const Hexahedra3D20 &rOther)
Definition: hexahedra_3d_20.h:342
friend class Hexahedra3D20
Definition: hexahedra_3d_20.h:1571
Quadrilateral3D8< TPointType > FaceType
Definition: hexahedra_3d_20.h:67
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: hexahedra_3d_20.h:790
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: hexahedra_3d_20.h:391
Hexahedra3D20(const std::string &GeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: hexahedra_3d_20.h:274
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: hexahedra_3d_20.h:130
Hexahedra3D20(const PointsArrayType &ThisPoints)
Definition: hexahedra_3d_20.h:258
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: hexahedra_3d_20.h:318
Matrix MatrixType
Definition: hexahedra_3d_20.h:173
GeometryData::IntegrationMethod IntegrationMethod
Definition: hexahedra_3d_20.h:77
BaseType::GeometriesArrayType GeometriesArrayType
Definition: hexahedra_3d_20.h:83
SizeType FacesNumber() const override
Returns the number of faces of the current geometry.
Definition: hexahedra_3d_20.h:609
double DomainSize() const override
Definition: hexahedra_3d_20.h:470
TPointType PointType
Definition: hexahedra_3d_20.h:88
KRATOS_CLASS_POINTER_DEFINITION(Hexahedra3D20)
void PrintData(std::ostream &rOStream) const override
Definition: hexahedra_3d_20.h:879
Matrix & ShapeFunctionsLocalGradients(Matrix &result, const CoordinatesArrayType &rPoint) const override
Definition: hexahedra_3d_20.h:904
Hexahedra3D20(const PointType &Point1, const PointType &Point2, const PointType &Point3, const PointType &Point4, const PointType &Point5, const PointType &Point6, const PointType &Point7, const PointType &Point8, const PointType &Point9, const PointType &Point10, const PointType &Point11, const PointType &Point12, const PointType &Point13, const PointType &Point14, const PointType &Point15, const PointType &Point16, const PointType &Point17, const PointType &Point18, const PointType &Point19, const PointType &Point20)
Definition: hexahedra_3d_20.h:179
double Length() const override
Definition: hexahedra_3d_20.h:420
BaseType::NormalType NormalType
Definition: hexahedra_3d_20.h:163
Hexahedra3D20 & operator=(Hexahedra3D20< TOtherPointType > const &rOther)
Definition: hexahedra_3d_20.h:360
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: hexahedra_3d_20.h:158
BaseType::IndexType IndexType
Definition: hexahedra_3d_20.h:95
Line3D3< TPointType > EdgeType
Definition: hexahedra_3d_20.h:66
GeometriesArrayType GenerateFaces() const override
Returns all faces of the current geometry.
Definition: hexahedra_3d_20.h:622
void PrintInfo(std::ostream &rOStream) const override
Definition: hexahedra_3d_20.h:865
double Area() const override
Definition: hexahedra_3d_20.h:438
Hexahedra3D20(const IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: hexahedra_3d_20.h:265
std::string Info() const override
Definition: hexahedra_3d_20.h:853
BaseType::JacobiansType JacobiansType
Definition: hexahedra_3d_20.h:151
double Volume() const override
This method calculate and return volume of this geometry.
Definition: hexahedra_3d_20.h:452
BaseType::Pointer Create(const IndexType NewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: hexahedra_3d_20.h:377
Hexahedra3D20(typename PointType::Pointer pPoint1, typename PointType::Pointer pPoint2, typename PointType::Pointer pPoint3, typename PointType::Pointer pPoint4, typename PointType::Pointer pPoint5, typename PointType::Pointer pPoint6, typename PointType::Pointer pPoint7, typename PointType::Pointer pPoint8, typename PointType::Pointer pPoint9, typename PointType::Pointer pPoint10, typename PointType::Pointer pPoint11, typename PointType::Pointer pPoint12, typename PointType::Pointer pPoint13, typename PointType::Pointer pPoint14, typename PointType::Pointer pPoint15, typename PointType::Pointer pPoint16, typename PointType::Pointer pPoint17, typename PointType::Pointer pPoint18, typename PointType::Pointer pPoint19, typename PointType::Pointer pPoint20)
Definition: hexahedra_3d_20.h:214
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: hexahedra_3d_20.h:323
Geometry< TPointType > BaseType
Definition: hexahedra_3d_20.h:61
Hexahedra3D20(Hexahedra3D20< TOtherPointType > const &rOther)
Definition: hexahedra_3d_20.h:307
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: hexahedra_3d_20.h:538
Hexahedra3D20(Hexahedra3D20 const &rOther)
Definition: hexahedra_3d_20.h:291
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: hexahedra_3d_20.h:137
BaseType::SizeType SizeType
Definition: hexahedra_3d_20.h:102
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: hexahedra_3d_20.h:483
BaseType::PointsArrayType PointsArrayType
Definition: hexahedra_3d_20.h:108
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: hexahedra_3d_20.h:144
~Hexahedra3D20() override
Definition: hexahedra_3d_20.h:315
BaseType::IntegrationPointType IntegrationPointType
Definition: hexahedra_3d_20.h:115
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 three node 3D line geometry with quadratic shape functions.
Definition: line_3d_3.h:66
This class defines the node.
Definition: node.h:65
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 eight node 3D quadrilateral geometry with quadratic shape functions.
Definition: quadrilateral_3d_8.h:49
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
Short class definition.
Definition: array_1d.h:61
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
#define KRATOS_ERROR
Definition: exception.h:161
#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 Hexahedra3D20< TPointType >::msGeometryData & msGeometryDimension(), Hexahedra3D20< 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
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