208 typename PointType::Pointer pPoint2,
209 typename PointType::Pointer pPoint3,
210 typename PointType::Pointer pPoint4,
211 typename PointType::Pointer pPoint5,
212 typename PointType::Pointer pPoint6,
213 typename PointType::Pointer pPoint7,
214 typename PointType::Pointer pPoint8,
215 typename PointType::Pointer pPoint9 )
230 :
BaseType( ThisPoints, &msGeometryData )
240 ) :
BaseType(GeometryId, rThisPoints, &msGeometryData)
247 const std::string& rGeometryName,
249 ) :
BaseType(rGeometryName, rThisPoints, &msGeometryData)
334 template<
class TOtherPo
intType>
357 return typename BaseType::Pointer(
new Quadrilateral2D9( rNewGeometryId, rThisPoints ) );
371 auto p_geometry =
typename BaseType::Pointer(
new Quadrilateral2D9( NewGeometryId, rGeometry.
Points() ) );
372 p_geometry->SetData(rGeometry.
GetData());
383 if ((LocalDirectionIndex == 0) || (LocalDirectionIndex == 1)) {
386 KRATOS_ERROR <<
"Possible direction index reaches from 0-1. Given direction index: "
387 << LocalDirectionIndex << std::endl;
437 KRATOS_WARNING(
"Quadrilateral2D9") <<
"Method not well defined. Replace with DomainSize() instead. This method preserves current behaviour but will be changed in June 2023 (returning zero instead)" << std::endl;
471 const double Tolerance = std::numeric_limits<double>::epsilon()
476 if ( (rResult[0] >= (-1.0-Tolerance)) && (rResult[0] <= (1.0+Tolerance)) ) {
477 if ( (rResult[1] >= (-1.0-Tolerance)) && (rResult[1] <= (1.0+Tolerance)) ) {
538 double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
539 double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
540 double fx3 = 1 - rPoint[0] * rPoint[0];
541 double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
542 double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
543 double fy3 = 1 - rPoint[1] * rPoint[1];
545 switch ( ShapeFunctionIndex )
566 KRATOS_ERROR <<
"Wrong index of shape function!" << *
this << std::endl;
585 if(rResult.size() != 9) rResult.
resize(9,
false);
587 double fx1 = 0.5 * ( rCoordinates[0] - 1 ) * rCoordinates[0];
588 double fx2 = 0.5 * ( rCoordinates[0] + 1 ) * rCoordinates[0];
589 double fx3 = 1 - rCoordinates[0] * rCoordinates[0];
590 double fy1 = 0.5 * ( rCoordinates[1] - 1 ) * rCoordinates[1];
591 double fy2 = 0.5 * ( rCoordinates[1] + 1 ) * rCoordinates[1];
592 double fy3 = 1 - rCoordinates[1] * rCoordinates[1];
594 rResult[0] = fx1*fy1 ;
595 rResult[1] = fx2*fy1 ;
596 rResult[2] = fx2*fy2 ;
597 rResult[3] = fx1*fy2 ;
598 rResult[4] = fx3*fy1 ;
599 rResult[5] = fx2*fy3 ;
600 rResult[6] = fx3*fy2 ;
601 rResult[7] = fx1*fy3 ;
602 rResult[8] = fx3*fy3 ;
619 std::string
Info()
const override
621 return "2 dimensional quadrilateral with nine nodes in 2D space";
632 rOStream <<
"2 dimensional quadrilateral with nine nodes in 2D space";
648 std::cout << std::endl;
654 rOStream <<
" Jacobian in the origin\t : " << jacobian;
666 = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
667 const int integration_points_number
671 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
673 Result[pnt] = localGradients[pnt];
687 = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
688 const int integration_points_number
692 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
694 Result[pnt] = localGradients[pnt];
711 double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
712 double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
713 double fx3 = 1 - rPoint[0] * rPoint[0];
714 double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
715 double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
716 double fy3 = 1 - rPoint[1] * rPoint[1];
718 double gx1 = 0.5 * ( 2 * rPoint[0] - 1 );
719 double gx2 = 0.5 * ( 2 * rPoint[0] + 1 );
720 double gx3 = -2.0 * rPoint[0];
721 double gy1 = 0.5 * ( 2 * rPoint[1] - 1 );
722 double gy2 = 0.5 * ( 2 * rPoint[1] + 1 );
723 double gy3 = -2.0 * rPoint[1];
725 rResult.
resize( 9, 2,
false );
727 rResult( 0, 0 ) = gx1 * fy1;
728 rResult( 0, 1 ) = fx1 * gy1;
729 rResult( 1, 0 ) = gx2 * fy1;
730 rResult( 1, 1 ) = fx2 * gy1;
731 rResult( 2, 0 ) = gx2 * fy2;
732 rResult( 2, 1 ) = fx2 * gy2;
733 rResult( 3, 0 ) = gx1 * fy2;
734 rResult( 3, 1 ) = fx1 * gy2;
735 rResult( 4, 0 ) = gx3 * fy1;
736 rResult( 4, 1 ) = fx3 * gy1;
737 rResult( 5, 0 ) = gx2 * fy3;
738 rResult( 5, 1 ) = fx2 * gy3;
739 rResult( 6, 0 ) = gx3 * fy2;
740 rResult( 6, 1 ) = fx3 * gy2;
741 rResult( 7, 0 ) = gx1 * fy3;
742 rResult( 7, 1 ) = fx1 * gy3;
743 rResult( 8, 0 ) = gx3 * fy3;
744 rResult( 8, 1 ) = fx3 * gy3;
756 rResult.
resize( 9, 2,
false );
758 rResult( 0, 0 ) = -1.0;
759 rResult( 0, 1 ) = -1.0;
760 rResult( 1, 0 ) = 1.0;
761 rResult( 1, 1 ) = -1.0;
762 rResult( 2, 0 ) = 1.0;
763 rResult( 2, 1 ) = 1.0;
764 rResult( 3, 0 ) = -1.0;
765 rResult( 3, 1 ) = 1.0;
766 rResult( 4, 0 ) = 0.0;
767 rResult( 4, 1 ) = -1.0;
768 rResult( 5, 0 ) = 1.0;
769 rResult( 5, 1 ) = 0.0;
770 rResult( 6, 0 ) = 0.0;
771 rResult( 6, 1 ) = 1.0;
772 rResult( 7, 0 ) = -1.0;
773 rResult( 7, 1 ) = 0.0;
774 rResult( 8, 0 ) = 0.0;
775 rResult( 8, 1 ) = 0.0;
789 double fx1 = 0.5 * ( rPoint.
X() - 1 ) * rPoint.
X();
790 double fx2 = 0.5 * ( rPoint.
X() + 1 ) * rPoint.
X();
791 double fx3 = 1 - rPoint.
X() * rPoint.
X();
792 double fy1 = 0.5 * ( rPoint.
Y() - 1 ) * rPoint.
Y();
793 double fy2 = 0.5 * ( rPoint.
Y() + 1 ) * rPoint.
Y();
794 double fy3 = 1 - rPoint.
Y() * rPoint.
Y();
796 double gx1 = 0.5 * ( 2 * rPoint.
X() - 1 );
797 double gx2 = 0.5 * ( 2 * rPoint.
X() + 1 );
798 double gx3 = -2.0 * rPoint.
X();
799 double gy1 = 0.5 * ( 2 * rPoint.
Y() - 1 );
800 double gy2 = 0.5 * ( 2 * rPoint.
Y() + 1 );
801 double gy3 = -2.0 * rPoint.
Y();
803 rResult.
resize( 9, 2,
false );
805 rResult( 0, 0 ) = gx1 * fy1;
806 rResult( 0, 1 ) = fx1 * gy1;
807 rResult( 1, 0 ) = gx2 * fy1;
808 rResult( 1, 1 ) = fx2 * gy1;
809 rResult( 2, 0 ) = gx2 * fy2;
810 rResult( 2, 1 ) = fx2 * gy2;
811 rResult( 3, 0 ) = gx1 * fy2;
812 rResult( 3, 1 ) = fx1 * gy2;
813 rResult( 4, 0 ) = gx3 * fy1;
814 rResult( 4, 1 ) = fx3 * gy1;
815 rResult( 5, 0 ) = gx2 * fy3;
816 rResult( 5, 1 ) = fx2 * gy3;
817 rResult( 6, 0 ) = gx3 * fy2;
818 rResult( 6, 1 ) = fx3 * gy2;
819 rResult( 7, 0 ) = gx1 * fy3;
820 rResult( 7, 1 ) = fx1 * gy3;
821 rResult( 8, 0 ) = gx3 * fy3;
822 rResult( 8, 1 ) = fx3 * gy3;
835 if ( rResult.size() != this->PointsNumber() )
844 rResult[
i].
resize( 2, 2,
false );
848 double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
849 double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
850 double fx3 = 1 - rPoint[0] * rPoint[0];
851 double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
852 double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
853 double fy3 = 1 - rPoint[1] * rPoint[1];
855 double gx1 = 0.5 * ( 2 * rPoint[0] - 1 );
856 double gx2 = 0.5 * ( 2 * rPoint[0] + 1 );
857 double gx3 = -2.0 * rPoint[0];
858 double gy1 = 0.5 * ( 2 * rPoint[1] - 1 );
859 double gy2 = 0.5 * ( 2 * rPoint[1] + 1 );
860 double gy3 = -2.0 * rPoint[1];
869 rResult[0]( 0, 0 ) = hx1 * fy1;
870 rResult[0]( 0, 1 ) = gx1 * gy1;
871 rResult[0]( 1, 0 ) = gx1 * gy1;
872 rResult[0]( 1, 1 ) = fx1 * hy1;
874 rResult[1]( 0, 0 ) = hx2 * fy1;
875 rResult[1]( 0, 1 ) = gx2 * gy1;
876 rResult[1]( 1, 0 ) = gx2 * gy1;
877 rResult[1]( 1, 1 ) = fx2 * hy1;
879 rResult[2]( 0, 0 ) = hx2 * fy2;
880 rResult[2]( 0, 1 ) = gx2 * gy2;
881 rResult[2]( 1, 0 ) = gx2 * gy2;
882 rResult[2]( 1, 1 ) = fx2 * hy2;
884 rResult[3]( 0, 0 ) = hx1 * fy2;
885 rResult[3]( 0, 1 ) = gx1 * gy2;
886 rResult[3]( 1, 0 ) = gx1 * gy2;
887 rResult[3]( 1, 1 ) = fx1 * hy2;
889 rResult[4]( 0, 0 ) = hx3 * fy1;
890 rResult[4]( 0, 1 ) = gx3 * gy1;
891 rResult[4]( 1, 0 ) = gx3 * gy1;
892 rResult[4]( 1, 1 ) = fx3 * hy1;
894 rResult[5]( 0, 0 ) = hx2 * fy3;
895 rResult[5]( 0, 1 ) = gx2 * gy3;
896 rResult[5]( 1, 0 ) = gx2 * gy3;
897 rResult[5]( 1, 1 ) = fx2 * hy3;
899 rResult[6]( 0, 0 ) = hx3 * fy2;
900 rResult[6]( 0, 1 ) = gx3 * gy2;
901 rResult[6]( 1, 0 ) = gx3 * gy2;
902 rResult[6]( 1, 1 ) = fx3 * hy2;
904 rResult[7]( 0, 0 ) = hx1 * fy3;
905 rResult[7]( 0, 1 ) = gx1 * gy3;
906 rResult[7]( 1, 0 ) = gx1 * gy3;
907 rResult[7]( 1, 1 ) = fx1 * hy3;
909 rResult[8]( 0, 0 ) = hx3 * fy3;
910 rResult[8]( 0, 1 ) = gx3 * gy3;
911 rResult[8]( 1, 0 ) = gx3 * gy3;
912 rResult[8]( 1, 1 ) = fx3 * hy3;
927 if ( rResult.size() != this->PointsNumber() )
944 for (
unsigned int j = 0;
j < 2;
j++ )
946 rResult[
i][
j].
resize( 2, 2,
false );
959 double gx1 = 0.5 * ( 2 * rPoint[0] - 1 );
960 double gx2 = 0.5 * ( 2 * rPoint[0] + 1 );
961 double gx3 = -2.0 * rPoint[0];
962 double gy1 = 0.5 * ( 2 * rPoint[1] - 1 );
963 double gy2 = 0.5 * ( 2 * rPoint[1] + 1 );
964 double gy3 = -2.0 * rPoint[1];
973 rResult[0][0]( 0, 0 ) = 0.0;
974 rResult[0][0]( 0, 1 ) = hx1 * gy1;
975 rResult[0][0]( 1, 0 ) = hx1 * gy1;
976 rResult[0][0]( 1, 1 ) = gx1 * hy1;
977 rResult[0][1]( 0, 0 ) = hx1 * gy1;
978 rResult[0][1]( 0, 1 ) = gx1 * hy1;
979 rResult[0][1]( 1, 0 ) = gx1 * hy1;
980 rResult[0][1]( 1, 1 ) = 0.0;
982 rResult[1][0]( 0, 0 ) = 0.0;
983 rResult[1][0]( 0, 1 ) = hx2 * gy1;
984 rResult[1][0]( 1, 0 ) = hx2 * gy1;
985 rResult[1][0]( 1, 1 ) = gx2 * hy1;
986 rResult[1][1]( 0, 0 ) = hx2 * gy1;
987 rResult[1][1]( 0, 1 ) = gx2 * hy1;
988 rResult[1][1]( 1, 0 ) = gx2 * hy1;
989 rResult[1][1]( 1, 1 ) = 0.0;
991 rResult[2][0]( 0, 0 ) = 0.0;
992 rResult[2][0]( 0, 1 ) = hx2 * gy2;
993 rResult[2][0]( 1, 0 ) = hx2 * gy2;
994 rResult[2][0]( 1, 1 ) = gx2 * hy2;
995 rResult[2][1]( 0, 0 ) = hx2 * gy2;
996 rResult[2][1]( 0, 1 ) = gx2 * hy2;
997 rResult[2][1]( 1, 0 ) = gx2 * hy2;
998 rResult[2][1]( 1, 1 ) = 0.0;
1000 rResult[3][0]( 0, 0 ) = 0.0;
1001 rResult[3][0]( 0, 1 ) = hx1 * gy2;
1002 rResult[3][0]( 1, 0 ) = hx1 * gy2;
1003 rResult[3][0]( 1, 1 ) = gx1 * hy2;
1004 rResult[3][1]( 0, 0 ) = hx1 * gy2;
1005 rResult[3][1]( 0, 1 ) = gx1 * hy2;
1006 rResult[3][1]( 1, 0 ) = gx1 * hy2;
1007 rResult[3][1]( 1, 1 ) = 0.0;
1009 rResult[4][0]( 0, 0 ) = 0.0;
1010 rResult[4][0]( 0, 1 ) = hx3 * gy1;
1011 rResult[4][0]( 1, 0 ) = hx3 * gy1;
1012 rResult[4][0]( 1, 1 ) = gx3 * hy1;
1013 rResult[4][1]( 0, 0 ) = hx3 * gy1;
1014 rResult[4][1]( 0, 1 ) = gx3 * hy1;
1015 rResult[4][1]( 1, 0 ) = gx3 * hy1;
1016 rResult[4][1]( 1, 1 ) = 0.0;
1018 rResult[5][0]( 0, 0 ) = 0.0;
1019 rResult[5][0]( 0, 1 ) = hx2 * gy3;
1020 rResult[5][0]( 1, 0 ) = hx2 * gy3;
1021 rResult[5][0]( 1, 1 ) = gx2 * hy3;
1022 rResult[5][1]( 0, 0 ) = hx2 * gy3;
1023 rResult[5][1]( 0, 1 ) = gx2 * hy3;
1024 rResult[5][1]( 1, 0 ) = gx2 * hy3;
1025 rResult[5][1]( 1, 1 ) = 0.0;
1027 rResult[6][0]( 0, 0 ) = 0.0;
1028 rResult[6][0]( 0, 1 ) = hx3 * gy2;
1029 rResult[6][0]( 1, 0 ) = hx3 * gy2;
1030 rResult[6][0]( 1, 1 ) = gx3 * hy2;
1031 rResult[6][1]( 0, 0 ) = hx3 * gy2;
1032 rResult[6][1]( 0, 1 ) = gx3 * hy2;
1033 rResult[6][1]( 1, 0 ) = gx3 * hy2;
1034 rResult[6][1]( 1, 1 ) = 0.0;
1036 rResult[7][0]( 0, 0 ) = 0.0;
1037 rResult[7][0]( 0, 1 ) = hx1 * gy3;
1038 rResult[7][0]( 1, 0 ) = hx1 * gy3;
1039 rResult[7][0]( 1, 1 ) = gx1 * hy3;
1040 rResult[7][1]( 0, 0 ) = hx1 * gy3;
1041 rResult[7][1]( 0, 1 ) = gx1 * hy3;
1042 rResult[7][1]( 1, 0 ) = gx1 * hy3;
1043 rResult[7][1]( 1, 1 ) = 0.0;
1045 rResult[8][0]( 0, 0 ) = 0.0;
1046 rResult[8][0]( 0, 1 ) = hx3 * gy3;
1047 rResult[8][0]( 1, 0 ) = hx3 * gy3;
1048 rResult[8][0]( 1, 1 ) = gx3 * hy3;
1049 rResult[8][1]( 0, 0 ) = hx3 * gy3;
1050 rResult[8][1]( 0, 1 ) = gx3 * hy3;
1051 rResult[8][1]( 1, 0 ) = gx3 * hy3;
1052 rResult[8][1]( 1, 1 ) = 0.0;
1082 void save(
Serializer& rSerializer )
const override
1112 static Matrix CalculateShapeFunctionsIntegrationPointsValues(
1118 const int integration_points_number = integration_points.size();
1120 const int points_number = 9;
1122 Matrix shape_function_values( integration_points_number, points_number );
1125 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1127 double fx1 = 0.5 * ( integration_points[pnt].X() - 1 ) * integration_points[pnt].
X();
1128 double fx2 = 0.5 * ( integration_points[pnt].X() + 1 ) * integration_points[pnt].
X();
1129 double fx3 = 1 - integration_points[pnt].X() * integration_points[pnt].X();
1130 double fy1 = 0.5 * ( integration_points[pnt].Y() - 1 ) * integration_points[pnt].
Y();
1131 double fy2 = 0.5 * ( integration_points[pnt].Y() + 1 ) * integration_points[pnt].
Y();
1132 double fy3 = 1 - integration_points[pnt].Y() * integration_points[pnt].Y();
1134 shape_function_values( pnt, 0 ) = ( fx1 * fy1 );
1135 shape_function_values( pnt, 1 ) = ( fx2 * fy1 );
1136 shape_function_values( pnt, 2 ) = ( fx2 * fy2 );
1137 shape_function_values( pnt, 3 ) = ( fx1 * fy2 );
1138 shape_function_values( pnt, 4 ) = ( fx3 * fy1 );
1139 shape_function_values( pnt, 5 ) = ( fx2 * fy3 );
1140 shape_function_values( pnt, 6 ) = ( fx3 * fy2 );
1141 shape_function_values( pnt, 7 ) = ( fx1 * fy3 );
1142 shape_function_values( pnt, 8 ) = ( fx3 * fy3 );
1145 return shape_function_values;
1166 const int integration_points_number = integration_points.size();
1172 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1174 double fx1 = 0.5 * ( integration_points[pnt].X() - 1 ) * integration_points[pnt].
X();
1175 double fx2 = 0.5 * ( integration_points[pnt].X() + 1 ) * integration_points[pnt].
X();
1176 double fx3 = 1 - integration_points[pnt].X() * integration_points[pnt].X();
1177 double fy1 = 0.5 * ( integration_points[pnt].Y() - 1 ) * integration_points[pnt].
Y();
1178 double fy2 = 0.5 * ( integration_points[pnt].Y() + 1 ) * integration_points[pnt].
Y();
1179 double fy3 = 1 - integration_points[pnt].Y() * integration_points[pnt].Y();
1181 double gx1 = 0.5 * ( 2 * integration_points[pnt].X() - 1 );
1182 double gx2 = 0.5 * ( 2 * integration_points[pnt].X() + 1 );
1183 double gx3 = -2.0 * integration_points[pnt].X();
1184 double gy1 = 0.5 * ( 2 * integration_points[pnt].Y() - 1 );
1185 double gy2 = 0.5 * ( 2 * integration_points[pnt].Y() + 1 );
1186 double gy3 = -2.0 * integration_points[pnt].Y();
1189 result( 0, 0 ) = gx1 * fy1;
1190 result( 0, 1 ) = fx1 * gy1;
1191 result( 1, 0 ) = gx2 * fy1;
1192 result( 1, 1 ) = fx2 * gy1;
1193 result( 2, 0 ) = gx2 * fy2;
1194 result( 2, 1 ) = fx2 * gy2;
1195 result( 3, 0 ) = gx1 * fy2;
1196 result( 3, 1 ) = fx1 * gy2;
1197 result( 4, 0 ) = gx3 * fy1;
1198 result( 4, 1 ) = fx3 * gy1;
1199 result( 5, 0 ) = gx2 * fy3;
1200 result( 5, 1 ) = fx2 * gy3;
1201 result( 6, 0 ) = gx3 * fy2;
1202 result( 6, 1 ) = fx3 * gy2;
1203 result( 7, 0 ) = gx1 * fy3;
1204 result( 7, 1 ) = fx1 * gy3;
1205 result( 8, 0 ) = gx3 * fy3;
1206 result( 8, 1 ) = fx3 * gy3;
1208 d_shape_f_values[pnt] = result;
1211 return d_shape_f_values;
1222 Quadrature < QuadrilateralGaussLegendreIntegrationPoints1,
1223 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1224 Quadrature < QuadrilateralGaussLegendreIntegrationPoints2,
1225 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1226 Quadrature < QuadrilateralGaussLegendreIntegrationPoints3,
1227 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1228 Quadrature < QuadrilateralGaussLegendreIntegrationPoints4,
1229 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1230 Quadrature < QuadrilateralGaussLegendreIntegrationPoints5,
1231 2, IntegrationPoint<3> >::GenerateIntegrationPoints()
1234 return integration_points;
1245 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1247 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1249 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1251 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1253 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1257 return shape_functions_values;
1268 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1270 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1272 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1274 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1276 Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1280 return shape_functions_local_gradients;
1301 std::istream& rIStream,
1306 std::ostream& rOStream,
1310 rOStream << std::endl;
1315 template<
class TPo
intType>
1316 const GeometryData Quadrilateral2D9<TPointType>::msGeometryData(
1319 Quadrilateral2D9<TPointType>::AllIntegrationPoints(),
1320 Quadrilateral2D9<TPointType>::AllShapeFunctionsValues(),
1321 AllShapeFunctionsLocalGradients()
1324 template<
class TPo
intType>
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
@ Kratos_Quadrilateral2D9
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
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
GeometryData::IntegrationMethod IntegrationMethod
Definition: geometry.h:122
Short class definition.
Definition: integration_point.h:52
static double ComputeArea2DGeometry(const Geometry< TPointType > &rGeometry)
This method calculates and returns the volume of the geometry from a 3D geometry.
Definition: integration_utilities.h:107
Definition: amatrix_interface.h:41
void swap(Matrix &Other)
Definition: amatrix_interface.h:289
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
An three node 2D line geometry with quadratic shape functions.
Definition: line_2d_3.h:63
This class defines the node.
Definition: node.h:65
double Y() const
Definition: point.h:187
double X() const
Definition: point.h:181
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
void push_back(const TPointerType &x)
Definition: pointer_vector.h:270
A nine node 2D quadrilateral geometry with quadratic shape functions.
Definition: quadrilateral_2d_9.h:49
Quadrilateral2D9 & operator=(Quadrilateral2D9< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_9.h:335
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: quadrilateral_2d_9.h:89
BaseType::JacobiansType JacobiansType
Definition: quadrilateral_2d_9.h:154
SizeType PointsNumberInDirection(IndexType LocalDirectionIndex) const override
Returns number of points per direction.
Definition: quadrilateral_2d_9.h:381
Quadrilateral2D9(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)
Definition: quadrilateral_2d_9.h:207
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: quadrilateral_2d_9.h:291
void PrintData(std::ostream &rOStream) const override
Definition: quadrilateral_2d_9.h:644
Quadrilateral2D9 & operator=(const Quadrilateral2D9 &rOther)
Definition: quadrilateral_2d_9.h:316
BaseType::Pointer Create(const IndexType rNewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_9.h:352
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: quadrilateral_2d_9.h:583
bool IsInside(const CoordinatesArrayType &rPoint, CoordinatesArrayType &rResult, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Definition: quadrilateral_2d_9.h:468
double Volume() const override
This method calculates and returns the volume of this geometry.
Definition: quadrilateral_2d_9.h:435
ShapeFunctionsThirdDerivativesType & ShapeFunctionsThirdDerivatives(ShapeFunctionsThirdDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:924
BaseType::PointsArrayType PointsArrayType
Definition: quadrilateral_2d_9.h:111
std::string Info() const override
Definition: quadrilateral_2d_9.h:619
BaseType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: quadrilateral_2d_9.h:169
BaseType::SizeType SizeType
Definition: quadrilateral_2d_9.h:105
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: quadrilateral_2d_9.h:140
BaseType::NormalType NormalType
Definition: quadrilateral_2d_9.h:183
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients()
Definition: quadrilateral_2d_9.h:683
ShapeFunctionsSecondDerivativesType & ShapeFunctionsSecondDerivatives(ShapeFunctionsSecondDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:833
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: quadrilateral_2d_9.h:496
double DomainSize() const override
Definition: quadrilateral_2d_9.h:455
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: quadrilateral_2d_9.h:509
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients(IntegrationMethod ThisMethod)
Definition: quadrilateral_2d_9.h:662
KRATOS_CLASS_POINTER_DEFINITION(Quadrilateral2D9)
BaseType::IndexType IndexType
Definition: quadrilateral_2d_9.h:98
Quadrilateral2D9(Quadrilateral2D9< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_9.h:280
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_9.h:366
Quadrilateral2D9(const std::string &rGeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: quadrilateral_2d_9.h:246
double Length() const override
Definition: quadrilateral_2d_9.h:406
Quadrilateral2D9(IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: quadrilateral_2d_9.h:237
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:708
Quadrilateral2D9(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)
Definition: quadrilateral_2d_9.h:189
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: quadrilateral_2d_9.h:161
virtual Matrix & ShapeFunctionsGradients(Matrix &rResult, PointType &rPoint)
Definition: quadrilateral_2d_9.h:787
GeometryData::IntegrationMethod IntegrationMethod
Definition: quadrilateral_2d_9.h:73
Quadrilateral2D9(const PointsArrayType &ThisPoints)
Definition: quadrilateral_2d_9.h:229
Geometry< TPointType > BaseType
Definition: quadrilateral_2d_9.h:58
Quadrilateral2D9(Quadrilateral2D9 const &rOther)
Definition: quadrilateral_2d_9.h:263
TPointType PointType
Definition: quadrilateral_2d_9.h:84
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: quadrilateral_2d_9.h:754
double Area() const override
Definition: quadrilateral_2d_9.h:423
BaseType::GeometriesArrayType GeometriesArrayType
Definition: quadrilateral_2d_9.h:79
~Quadrilateral2D9() override
Definition: quadrilateral_2d_9.h:289
void PrintInfo(std::ostream &rOStream) const override
Definition: quadrilateral_2d_9.h:630
friend class Quadrilateral2D9
Definition: quadrilateral_2d_9.h:1287
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: quadrilateral_2d_9.h:296
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:535
BaseType::ShapeFunctionsThirdDerivativesType ShapeFunctionsThirdDerivativesType
Definition: quadrilateral_2d_9.h:177
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: quadrilateral_2d_9.h:133
Line2D3< TPointType > EdgeType
Definition: quadrilateral_2d_9.h:63
BaseType::IntegrationPointType IntegrationPointType
Definition: quadrilateral_2d_9.h:118
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: quadrilateral_2d_9.h:126
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: quadrilateral_2d_9.h:147
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
#define KRATOS_WARNING(label)
Definition: logger.h:265
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
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
const GeometryData Quadrilateral2D9< TPointType >::msGeometryData & msGeometryDimension(), Quadrilateral2D9< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
def load(f)
Definition: ode_solve.py:307
int j
Definition: quadrature.py:648
float temp
Definition: rotating_cone.py:85
integer i
Definition: TensorModule.f:17