210 Quadrilateral2D8(
typename PointType::Pointer pPoint1,
typename PointType::Pointer pPoint2,
211 typename PointType::Pointer pPoint3,
typename PointType::Pointer pPoint4,
212 typename PointType::Pointer pPoint5,
typename PointType::Pointer pPoint6,
213 typename PointType::Pointer pPoint7,
typename PointType::Pointer pPoint8 )
227 :
BaseType( ThisPoints, &msGeometryData )
237 ) :
BaseType(GeometryId, rThisPoints, &msGeometryData)
244 const std::string& rGeometryName,
246 ) :
BaseType(rGeometryName, rThisPoints, &msGeometryData)
330 template<
class TOtherPo
intType>
353 return typename BaseType::Pointer(
new Quadrilateral2D8( NewGeometryId, rThisPoints ) );
367 auto p_geometry =
typename BaseType::Pointer(
new Quadrilateral2D8( NewGeometryId, rGeometry.
Points() ) );
368 p_geometry->SetData(rGeometry.
GetData());
379 if ((LocalDirectionIndex == 0) || (LocalDirectionIndex == 1)) {
382 KRATOS_ERROR <<
"Possible direction index reaches from 0-1. Given direction index: "
383 << LocalDirectionIndex << std::endl;
406 point[0] = 1.0/3.0; point[1] = 1.0/3.0; point[2] = 1.0/3.0;
468 const double Tolerance = std::numeric_limits<double>::epsilon()
473 if ( std::abs( rResult[0] ) <= (1.0 + Tolerance) )
475 if ( std::abs( rResult[1] ) <= (1.0 + Tolerance) )
510 CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
512 Matrix shape_functions_values =
513 CalculateShapeFunctionsIntegrationPointsValues( ThisMethod );
516 if ( rResult.size() != this->IntegrationPointsNumber( ThisMethod ) )
533 jacobian( 0, 0 ) += ( this->
GetPoint(
i ).X() ) * ( shape_functions_gradients[pnt](
i, 0 ) );
534 jacobian( 0, 1 ) += ( this->
GetPoint(
i ).X() ) * ( shape_functions_gradients[pnt](
i, 1 ) );
535 jacobian( 1, 0 ) += ( this->
GetPoint(
i ).Y() ) * ( shape_functions_gradients[pnt](
i, 0 ) );
536 jacobian( 1, 1 ) += ( this->
GetPoint(
i ).Y() ) * ( shape_functions_gradients[pnt](
i, 1 ) );
539 rResult[pnt] = jacobian;
565 Matrix & DeltaPosition )
const override
569 CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
571 Matrix shape_functions_values =
572 CalculateShapeFunctionsIntegrationPointsValues( ThisMethod );
575 if ( rResult.size() != this->IntegrationPointsNumber( ThisMethod ) )
592 jacobian( 0, 0 ) += ( this->
GetPoint(
i ).X() - DeltaPosition(
i,0) ) * ( shape_functions_gradients[pnt](
i, 0 ) );
593 jacobian( 0, 1 ) += ( this->
GetPoint(
i ).X() - DeltaPosition(
i,0) ) * ( shape_functions_gradients[pnt](
i, 1 ) );
594 jacobian( 1, 0 ) += ( this->
GetPoint(
i ).Y() - DeltaPosition(
i,1) ) * ( shape_functions_gradients[pnt](
i, 0 ) );
595 jacobian( 1, 1 ) += ( this->
GetPoint(
i ).Y() - DeltaPosition(
i,1) ) * ( shape_functions_gradients[pnt](
i, 1 ) );
598 rResult[pnt] = jacobian;
627 if (rResult.size1() != 2 || rResult.size2() != 2)
628 rResult.
resize( 2, 2 ,
false);
632 Matrix ShapeFunctionsGradientInIntegrationPoint = shape_functions_gradients( IntegrationPointIndex );
635 ShapeFunctionsValuesInIntegrationPoint =
row( CalculateShapeFunctionsIntegrationPointsValues( ThisMethod ), IntegrationPointIndex );
642 ( this->
GetPoint(
i ).X() ) * ( ShapeFunctionsGradientInIntegrationPoint(
i, 0 ) );
644 ( this->
GetPoint(
i ).X() ) * ( ShapeFunctionsGradientInIntegrationPoint(
i, 1 ) );
646 ( this->
GetPoint(
i ).Y() ) * ( ShapeFunctionsGradientInIntegrationPoint(
i, 0 ) );
648 ( this->
GetPoint(
i ).Y() ) * ( ShapeFunctionsGradientInIntegrationPoint(
i, 1 ) );
673 if (rResult.size1() != 2 || rResult.size2() != 2)
674 rResult.
resize( 2, 2 ,
false);
677 Matrix shape_functions_gradients;
684 rResult( 0, 0 ) += ( this->
GetPoint(
i ).X() ) * ( shape_functions_gradients(
i, 0 ) );
685 rResult( 0, 1 ) += ( this->
GetPoint(
i ).X() ) * ( shape_functions_gradients(
i, 1 ) );
686 rResult( 1, 0 ) += ( this->
GetPoint(
i ).Y() ) * ( shape_functions_gradients(
i, 0 ) );
687 rResult( 1, 1 ) += ( this->
GetPoint(
i ).Y() ) * ( shape_functions_gradients(
i, 1 ) );
710 if( rResult.size() != this->IntegrationPointsNumber( ThisMethod ) )
715 this->
Jacobian( J, pnt, ThisMethod);
716 rResult[pnt] = ((
J( 0, 0 )*
J( 1, 1 ) ) - (
J( 0, 1 )*
J( 1, 0 ) ) );
746 jacobian =
Jacobian( jacobian, IntegrationPointIndex, ThisMethod );
747 return(( jacobian( 0, 0 )*jacobian( 1, 1 ) ) - ( jacobian( 0, 1 )*jacobian( 1, 0 ) ) );
770 jacobian =
Jacobian( jacobian, rPoint );
771 return(( jacobian( 0, 0 )*jacobian( 1, 1 ) ) - ( jacobian( 0, 1 )*jacobian( 1, 0 ) ) );
797 if ( rResult.size() != this->IntegrationPointsNumber( ThisMethod ) )
808 Matrix tempMatrix( 2, 2 );
843 tempMatrix =
Jacobian( tempMatrix, IntegrationPointIndex, ThisMethod );
851 KRATOS_ERROR <<
"Zero determinant of jacobian." << *
this << std::endl;
855 rResult.
resize( 2, 2 ,
false);
858 rResult( 0, 0 ) = ( tempMatrix( 1, 1 ) ) / ( det_j );
860 rResult( 1, 0 ) = -( tempMatrix( 1, 0 ) ) / ( det_j );
862 rResult( 0, 1 ) = -( tempMatrix( 0, 1 ) ) / ( det_j );
888 tempMatrix.
resize( 2, 2,
false );
889 tempMatrix = this->
Jacobian( tempMatrix, rPoint );
896 KRATOS_ERROR <<
"Zero determinant of jacobian." << *
this << std::endl;
900 rResult.
resize( 2, 2 ,
false);
903 rResult( 0, 0 ) = ( tempMatrix( 1, 1 ) ) / ( det_j );
905 rResult( 1, 0 ) = -( tempMatrix( 1, 0 ) ) / ( det_j );
907 rResult( 0, 1 ) = -( tempMatrix( 0, 1 ) ) / ( det_j );
909 rResult( 1, 1 ) = ( tempMatrix( 0, 0 ) ) / ( det_j );
971 switch ( ShapeFunctionIndex )
975 return -(( 1.0 - rPoint[0] )*( 1.0 - rPoint[1] )
976 *( 1.0 + rPoint[0] + rPoint[1] ) ) / 4.0;
978 return -(( 1.0 + rPoint[0] )*( 1.0 - rPoint[1] )
979 *( 1.0 - rPoint[0] + rPoint[1] ) ) / 4.0;
981 return -(( 1.0 + rPoint[0] )*( 1.0 + rPoint[1] )
982 *( 1.0 - rPoint[0] - rPoint[1] ) ) / 4.0;
984 return -(( 1.0 - rPoint[0] )*( 1.0 + rPoint[1] )
985 *( 1.0 + rPoint[0] - rPoint[1] ) ) / 4.0;
988 return (( 1.0 -rPoint[0]*rPoint[0] )
989 *( 1.0 - rPoint[1] ) ) / 2.0;
991 return (( 1.0 + rPoint[0] )
992 *( 1.0 - rPoint[1]*rPoint[1] ) ) / 2.0 ;
994 return (( 1.0 -rPoint[0]*rPoint[0] )
995 *( 1.0 + rPoint[1] ) ) / 2.0 ;
997 return (( 1.0 -rPoint[0] )
998 *( 1.0 - rPoint[1]*rPoint[1] ) ) / 2.0 ;
1000 KRATOS_ERROR <<
"Wrong index of shape function!" << *
this << std::endl;
1019 if(rResult.size() != 8)
1025 rResult[0] = -(( 1.0 - rCoordinates[0] )*( 1.0 - rCoordinates[1] )
1026 *( 1.0 + rCoordinates[0] + rCoordinates[1] ) ) / 4.0;
1027 rResult[1] = -(( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] )
1028 *( 1.0 - rCoordinates[0] + rCoordinates[1] ) ) / 4.0;
1029 rResult[2] = -(( 1.0 + rCoordinates[0] )*( 1.0 + rCoordinates[1] )
1030 *( 1.0 - rCoordinates[0] - rCoordinates[1] ) ) / 4.0;
1031 rResult[3] = -(( 1.0 - rCoordinates[0] )*( 1.0 + rCoordinates[1] )
1032 *( 1.0 + rCoordinates[0] - rCoordinates[1] ) ) / 4.0;
1035 rResult[4] = (( 1.0 -rCoordinates[0]*rCoordinates[0] )
1036 *( 1.0 - rCoordinates[1] ) ) / 2.0;
1037 rResult[5] = (( 1.0 + rCoordinates[0] )
1038 *( 1.0 - rCoordinates[1]*rCoordinates[1] ) ) / 2.0 ;
1039 rResult[6] = (( 1.0 -rCoordinates[0]*rCoordinates[0] )
1040 *( 1.0 + rCoordinates[1] ) ) / 2.0 ;
1041 rResult[7] = (( 1.0 -rCoordinates[0] )
1042 *( 1.0 - rCoordinates[1]*rCoordinates[1] ) ) / 2.0 ;
1068 const unsigned int integration_points_number =
1071 if ( integration_points_number == 0 )
1073 KRATOS_ERROR <<
"This integration method is not supported" << *
this << std::endl;
1077 if ( rResult.size() != integration_points_number )
1087 CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
1095 for (
unsigned int pnt = 0; pnt < integration_points_number; pnt++ )
1097 rResult[pnt].
resize( 4, 2,
false );
1099 for (
unsigned int i = 0;
i < 4;
i++ )
1101 for (
unsigned int j = 0;
j < 2;
j++ )
1103 rResult[pnt](
i,
j ) =
1104 ( locG[pnt](
i, 0 ) * invJ[pnt](
j, 0 ) )
1105 + ( locG[pnt](
i, 1 ) * invJ[pnt](
j, 1 ) );
1123 return "2 dimensional quadrilateral with eight nodes in 2D space";
1134 rOStream <<
"2 dimensional quadrilateral with eight nodes in 2D space";
1150 std::cout << std::endl;
1156 rOStream <<
" Jacobian in the origin\t : " << jacobian;
1168 = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
1169 const int integration_points_number
1173 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1175 Result[pnt] = localGradients[pnt];
1189 = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
1190 const int integration_points_number
1194 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1196 Result[pnt] = localGradients[pnt];
1215 rResult.
resize( 8, 2 ,
false);
1219 rResult( 0, 0 ) = - ((-1.0 + rPoint[1])*( 2.0 * rPoint[0] + rPoint[1])) / 4.0;
1220 rResult( 0, 1 ) = - ((-1.0 + rPoint[0])*( 2.0 * rPoint[1] + rPoint[0])) / 4.0;
1222 rResult( 1, 0 ) = ((-1.0 + rPoint[1])*(-2.0 * rPoint[0] + rPoint[1])) / 4.0;
1223 rResult( 1, 1 ) = (( 1.0 + rPoint[0])*( 2.0 * rPoint[1] - rPoint[0])) / 4.0;
1225 rResult( 2, 0 ) = (( 1.0 + rPoint[1])*( 2.0 * rPoint[0] + rPoint[1])) / 4.0;
1226 rResult( 2, 1 ) = (( 1.0 + rPoint[0])*( 2.0 * rPoint[1] + rPoint[0])) / 4.0;
1228 rResult( 3, 0 ) = - (( 1.0 + rPoint[1])*(-2.0 * rPoint[0] + rPoint[1])) / 4.0;
1229 rResult( 3, 1 ) = - ((-1.0 + rPoint[0])*( 2.0 * rPoint[1] - rPoint[0])) / 4.0;
1232 rResult( 4, 0 ) = rPoint[0] * (-1.0 + rPoint[1]);
1233 rResult( 4, 1 ) = ((1.0 + rPoint[0]) * (-1.0 + rPoint[0])) / 2.0;
1235 rResult( 5, 0 ) = - ((1.0 + rPoint[1]) * (-1.0 + rPoint[1])) / 2.0;
1236 rResult( 5, 1 ) = - rPoint[1] * ( 1.0 + rPoint[0]);
1238 rResult( 6, 0 ) = - rPoint[0] * ( 1.0 + rPoint[1]);
1239 rResult( 6, 1 ) = - ((1.0 + rPoint[0]) * (-1.0 + rPoint[0])) / 2.0;
1241 rResult( 7, 0 ) = ((1.0 + rPoint[1]) * (-1.0 + rPoint[1])) / 2.0;
1242 rResult( 7, 1 ) = rPoint[1] * (-1.0 + rPoint[0]);
1254 rResult.
resize( 8, 2 ,
false);
1256 rResult( 0, 0 ) = -1.0;
1257 rResult( 0, 1 ) = -1.0;
1258 rResult( 1, 0 ) = 1.0;
1259 rResult( 1, 1 ) = -1.0;
1260 rResult( 2, 0 ) = 1.0;
1261 rResult( 2, 1 ) = 1.0;
1262 rResult( 3, 0 ) = -1.0;
1263 rResult( 3, 1 ) = 1.0;
1264 rResult( 4, 0 ) = 0.0;
1265 rResult( 4, 1 ) = -1.0;
1266 rResult( 5, 0 ) = 1.0;
1267 rResult( 5, 1 ) = 0.0;
1268 rResult( 6, 0 ) = 0.0;
1269 rResult( 6, 1 ) = 1.0;
1270 rResult( 7, 0 ) = -1.0;
1271 rResult( 7, 1 ) = 0.0;
1285 rResult.
resize( 8, 2 ,
false);
1303 if ( rResult.size() != this->PointsNumber() )
1313 rResult[
i].
resize( 2, 2 ,
false);
1317 rResult[0]( 0, 0 ) = ( 4.0 - 4 * rPoint[1] ) / 8.0;
1319 rResult[0]( 0, 1 ) = ( -2.0 ) * ( 1.0 + 2.0 * rPoint[0] + rPoint[1] - 1.0 ) / 8.0
1320 + (( -1.0 + rPoint[1] ) * ( -2.0 ) ) / 8.0;
1321 rResult[0]( 1, 0 ) = ( -2.0 ) * ( 1.0 + rPoint[0] + 2.0 * rPoint[1] - 1.0 ) / 8.0
1322 + (( -1.0 + rPoint[0] ) * ( -2.0 ) ) / 8.0;
1323 rResult[0]( 1, 1 ) = (( -1.0 + rPoint[0] ) * ( -2.0 ) * ( 2.0 ) ) / 8.0;
1325 rResult[1]( 0, 0 ) = -(( -1.0 + rPoint[1] ) * ( -2.0 ) * ( -2.0 ) ) / 8.0;
1326 rResult[1]( 0, 1 ) = ( 2.0 ) * ( 1.0 - 2.0 * rPoint[0] + rPoint[1] - 1.0 ) / 8.0
1327 - (( -1.0 + rPoint[1] ) * ( -2.0 ) ) / 8.0;
1328 rResult[1]( 1, 0 ) = ( -1.0 + rPoint[0] - 2.0 * rPoint[1] + 1.0 ) * ( -2.0 ) / 8.0
1329 + (( 1.0 + rPoint[0] ) * ( -2.0 ) ) / 8.0;
1330 rResult[1]( 1, 1 ) = (( 1.0 + rPoint[0] ) * ( -2.0 ) * ( -2.0 ) ) / 8.0;
1332 rResult[2]( 0, 0 ) = -(( 1.0 + rPoint[1] ) * ( 2.0 ) * ( -2.0 ) ) / 8.0;
1333 rResult[2]( 0, 1 ) = -(( 1.0 ) * ( -1.0 + 2.0 * rPoint[0] + rPoint[1] + 1.0 ) * ( -2.0 ) ) / 8.0
1334 - (( 1.0 + rPoint[1] ) * ( -2.0 ) ) / 8.0;
1335 rResult[2]( 1, 0 ) = -(( 1.0 ) * ( -1.0 + rPoint[0] + 2.0 * rPoint[1] + 1.0 ) * ( -2.0 ) ) / 8.0
1336 - (( 1.0 + rPoint[0] ) * ( -2.0 ) ) / 8.0;
1337 rResult[2]( 1, 1 ) = -(( 1.0 + rPoint[0] ) * ( 2.0 ) * ( -2.0 ) ) / 8.0;
1339 rResult[3]( 0, 0 ) = (( 1.0 + rPoint[1] ) * ( -2.0 ) * ( -2.0 ) ) / 8.0;
1340 rResult[3]( 0, 1 ) = (( 1.0 ) * ( -1.0 - 2.0 * rPoint[0] + rPoint[1] + 1.0 ) * ( -2.0 ) ) / 8.0
1341 + (( 1.0 + rPoint[1] ) * ( -2.0 ) ) / 8.0;
1342 rResult[3]( 1, 0 ) = -(( -2.0 ) * ( 1.0 + rPoint[0] - 2.0 * rPoint[1] - 1.0 ) ) / 8.0
1343 - (( -1.0 + rPoint[0] ) * ( -2.0 ) ) / 8.0;
1344 rResult[3]( 1, 1 ) = -(( -1.0 + rPoint[0] ) * ( -2.0 ) * ( -2.0 ) ) / 8.0;
1346 rResult[4]( 0, 0 ) = -(( -1.0 + rPoint[1] ) * ( -2.0 ) ) / 2.0;
1347 rResult[4]( 0, 1 ) = -( rPoint[0] * ( -2.0 ) ) / 2.0;
1348 rResult[4]( 1, 0 ) = -(( 2.0 * rPoint[0] ) * ( -2.0 ) ) / 4.0;
1349 rResult[4]( 1, 1 ) = 0.0;
1351 rResult[5]( 0, 0 ) = 0.0;
1352 rResult[5]( 0, 1 ) = (( 2.0 * rPoint[1] ) * ( -2.0 ) ) / 4.0;
1353 rResult[5]( 1, 0 ) = ( rPoint[1] * ( -2.0 ) ) / 2.0;
1354 rResult[5]( 1, 1 ) = (( 1.0 + rPoint[0] ) * ( -2.0 ) ) / 2.0;
1356 rResult[6]( 0, 0 ) = (( 1.0 + rPoint[1] ) * ( -2.0 ) ) / 2.0;
1357 rResult[6]( 0, 1 ) = ( rPoint[0] * ( -2.0 ) ) / 2.0;
1358 rResult[6]( 1, 0 ) = (( 2.0 * rPoint[0] ) * ( -2.0 ) ) / 4.0;
1359 rResult[6]( 1, 1 ) = 0.0;
1361 rResult[7]( 0, 0 ) = 0.0;
1362 rResult[7]( 0, 1 ) = -(( 2.0 * rPoint[1] ) * ( -2.0 ) ) / 4.0;
1363 rResult[7]( 1, 0 ) = -( rPoint[1] * ( -2.0 ) ) / 2.0;
1364 rResult[7]( 1, 1 ) = -(( -1.0 + rPoint[0] ) * ( -2.0 ) ) / 2.0;
1382 if ( rResult.size() != this->PointsNumber() )
1399 for (
int j = 0;
j < 2;
j++ )
1401 rResult[
i][
j].
resize( 2, 2 ,
false);
1406 rResult[0][0]( 0, 0 ) = 0.0;
1408 rResult[0][0]( 0, 1 ) = -0.5;
1409 rResult[0][0]( 1, 0 ) = -0.5;
1410 rResult[0][0]( 1, 1 ) = -0.5;
1411 rResult[0][1]( 0, 0 ) = -0.5;
1412 rResult[0][1]( 0, 1 ) = -0.5;
1413 rResult[0][1]( 1, 0 ) = -0.5;
1414 rResult[0][1]( 1, 1 ) = 0.0;
1415 rResult[1][0]( 0, 0 ) = 0.0;
1416 rResult[1][0]( 0, 1 ) = -0.5;
1417 rResult[1][0]( 1, 0 ) = -0.5;
1418 rResult[1][0]( 1, 1 ) = 0.5;
1419 rResult[1][1]( 0, 0 ) = -0.5;
1420 rResult[1][1]( 0, 1 ) = 0.5;
1421 rResult[1][1]( 1, 0 ) = 0.5;
1422 rResult[1][1]( 1, 1 ) = 0.0;
1423 rResult[2][0]( 0, 0 ) = 0.0;
1424 rResult[2][0]( 0, 1 ) = 0.5;
1425 rResult[2][0]( 1, 0 ) = 0.5;
1426 rResult[2][0]( 1, 1 ) = 0.5;
1427 rResult[2][1]( 0, 0 ) = 0.5;
1428 rResult[2][1]( 0, 1 ) = 0.5;
1429 rResult[2][1]( 1, 0 ) = 0.5;
1430 rResult[2][1]( 1, 1 ) = 0.0;
1431 rResult[3][0]( 0, 0 ) = 0.0;
1432 rResult[3][0]( 0, 1 ) = 0.5;
1433 rResult[3][0]( 1, 0 ) = 0.5;
1434 rResult[3][0]( 1, 1 ) = -0.5;
1435 rResult[3][1]( 0, 0 ) = 0.5;
1436 rResult[3][1]( 0, 1 ) = -0.5;
1437 rResult[3][1]( 1, 0 ) = -0.5;
1438 rResult[3][1]( 1, 1 ) = 0.0;
1439 rResult[4][0]( 0, 0 ) = 0.0;
1440 rResult[4][0]( 0, 1 ) = 1.0;
1441 rResult[4][0]( 1, 0 ) = 1.0;
1442 rResult[4][0]( 1, 1 ) = 0.0;
1443 rResult[4][1]( 0, 0 ) = 1.0;
1444 rResult[4][1]( 0, 1 ) = 0.0;
1445 rResult[4][1]( 1, 0 ) = 0.0;
1446 rResult[4][1]( 1, 1 ) = 0.0;
1447 rResult[5][0]( 0, 0 ) = 0.0;
1448 rResult[5][0]( 0, 1 ) = 0.0;
1449 rResult[5][0]( 1, 0 ) = 0.0;
1450 rResult[5][0]( 1, 1 ) = -1.0;
1451 rResult[5][1]( 0, 0 ) = 0.0;
1452 rResult[5][1]( 0, 1 ) = -1.0;
1453 rResult[5][1]( 1, 0 ) = 1.0;
1454 rResult[5][1]( 1, 1 ) = 0.0;
1455 rResult[6][0]( 0, 0 ) = 0.0;
1456 rResult[6][0]( 0, 1 ) = -1.0;
1457 rResult[6][0]( 1, 0 ) = -1.0;
1458 rResult[6][0]( 1, 1 ) = 0.0;
1459 rResult[6][1]( 0, 0 ) = -1.0;
1460 rResult[6][1]( 0, 1 ) = 0.0;
1461 rResult[6][1]( 1, 0 ) = 0.0;
1462 rResult[6][1]( 1, 1 ) = 0.0;
1463 rResult[7][0]( 0, 0 ) = 0.0;
1464 rResult[7][0]( 0, 1 ) = 0.0;
1465 rResult[7][0]( 1, 0 ) = 0.0;
1466 rResult[7][0]( 1, 1 ) = 1.0;
1467 rResult[7][1]( 0, 0 ) = 0.0;
1468 rResult[7][1]( 0, 1 ) = 1.0;
1469 rResult[7][1]( 1, 0 ) = -1.0;
1470 rResult[7][1]( 1, 0 ) = 0.0;
1498 void save(
Serializer& rSerializer )
const override
1531 const unsigned int integration_points_number = integration_points.size();
1534 const unsigned int points_number = 8;
1537 Matrix shape_function_values( integration_points_number, points_number );
1540 for (
unsigned int pnt = 0; pnt < integration_points_number; pnt++ )
1543 row( shape_function_values, pnt )[0] =
1544 -(( 1.0 - integration_points[pnt].X() )
1545 * ( 1.0 - integration_points[pnt].
Y() )
1546 * ( 1.0 + integration_points[pnt].
X()
1547 + integration_points[pnt].Y() ) ) / 4.0;
1548 row( shape_function_values, pnt )[1] =
1549 -(( 1.0 + integration_points[pnt].X() )
1550 * ( 1.0 - integration_points[pnt].
Y() ) * ( 1.0
1551 - integration_points[pnt].
X()
1552 + integration_points[pnt].Y() ) ) / 4.0;
1553 row( shape_function_values, pnt )[2] =
1554 -(( 1.0 + integration_points[pnt].X() )
1555 * ( 1.0 + integration_points[pnt].
Y() ) * ( 1.0
1556 - integration_points[pnt].
X()
1557 - integration_points[pnt].Y() ) ) / 4.0;
1558 row( shape_function_values, pnt )[3] =
1559 -(( 1.0 - integration_points[pnt].X() ) * ( 1.0
1560 + integration_points[pnt].
Y() ) * ( 1.0 ) * ( 1.0
1561 + integration_points[pnt].X()
1562 - integration_points[pnt].Y() ) ) / 4.0;
1565 row( shape_function_values, pnt )[4] =
1566 (( 1.0 - integration_points[pnt].X()
1567 * integration_points[pnt].X() )
1568 * ( 1.0 - integration_points[pnt].
Y() ) ) / 2.0;
1569 row( shape_function_values, pnt )[5] =
1570 (( 1.0 + integration_points[pnt].X() )
1571 * ( 1.0 - integration_points[pnt].
Y()
1572 * integration_points[pnt].Y() ) ) / 2.0 ;
1573 row( shape_function_values, pnt )[6] =
1574 (( 1.0 - integration_points[pnt].X()
1575 * integration_points[pnt].X() )
1576 * ( 1.0 + integration_points[pnt].
Y() ) ) / 2.0 ;
1577 row( shape_function_values, pnt )[7] =
1578 (( 1.0 - integration_points[pnt].X() )
1579 * ( 1.0 - integration_points[pnt].
Y()
1580 * integration_points[pnt].Y() ) ) / 2.0 ;
1583 return shape_function_values;
1603 const unsigned int integration_points_number = integration_points.size();
1609 for (
unsigned int pnt = 0; pnt < integration_points_number; pnt++ )
1614 result( 0, 0 ) = - ((-1.0 + integration_points[pnt].Y())*( 2.0 * integration_points[pnt].
X() + integration_points[pnt].Y())) / 4.0;
1615 result( 0, 1 ) = - ((-1.0 + integration_points[pnt].X())*( 2.0 * integration_points[pnt].
Y() + integration_points[pnt].X())) / 4.0;
1617 result( 1, 0 ) = ((-1.0 + integration_points[pnt].Y())*(-2.0 * integration_points[pnt].
X() + integration_points[pnt].Y())) / 4.0;
1618 result( 1, 1 ) = (( 1.0 + integration_points[pnt].X())*( 2.0 * integration_points[pnt].
Y() - integration_points[pnt].X())) / 4.0;
1620 result( 2, 0 ) = (( 1.0 + integration_points[pnt].Y())*( 2.0 * integration_points[pnt].
X() + integration_points[pnt].Y())) / 4.0;
1621 result( 2, 1 ) = (( 1.0 + integration_points[pnt].X())*( 2.0 * integration_points[pnt].
Y() + integration_points[pnt].X())) / 4.0;
1623 result( 3, 0 ) = - (( 1.0 + integration_points[pnt].Y())*(-2.0 * integration_points[pnt].
X() + integration_points[pnt].Y())) / 4.0;
1624 result( 3, 1 ) = - ((-1.0 + integration_points[pnt].X())*( 2.0 * integration_points[pnt].
Y() - integration_points[pnt].X())) / 4.0;
1627 result( 4, 0 ) = integration_points[pnt].X() * (-1.0 + integration_points[pnt].Y());
1628 result( 4, 1 ) = ((1.0 + integration_points[pnt].X()) * (-1.0 + integration_points[pnt].
X())) / 2.0;
1630 result( 5, 0 ) = - ((1.0 + integration_points[pnt].Y()) * (-1.0 + integration_points[pnt].
Y())) / 2.0;
1631 result( 5, 1 ) = - integration_points[pnt].Y() * ( 1.0 + integration_points[pnt].X());
1633 result( 6, 0 ) = - integration_points[pnt].X() * ( 1.0 + integration_points[pnt].Y());
1634 result( 6, 1 ) = - ((1.0 + integration_points[pnt].X()) * (-1.0 + integration_points[pnt].
X())) / 2.0;
1636 result( 7, 0 ) = ((1.0 + integration_points[pnt].Y()) * (-1.0 + integration_points[pnt].
Y())) / 2.0;
1637 result( 7, 1 ) = integration_points[pnt].Y() * (-1.0 + integration_points[pnt].X());
1639 d_shape_f_values[pnt] = result;
1642 return d_shape_f_values;
1654 Quadrature < QuadrilateralGaussLegendreIntegrationPoints1,
1655 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1656 Quadrature < QuadrilateralGaussLegendreIntegrationPoints2,
1657 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1658 Quadrature < QuadrilateralGaussLegendreIntegrationPoints3,
1659 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1660 Quadrature < QuadrilateralGaussLegendreIntegrationPoints4,
1661 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1662 Quadrature < QuadrilateralGaussLegendreIntegrationPoints5,
1663 2, IntegrationPoint<3> >::GenerateIntegrationPoints()
1666 return integration_points;
1677 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1679 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1681 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1683 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1685 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1689 return shape_functions_values;
1696 AllShapeFunctionsLocalGradients()
1701 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1703 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1705 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1707 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1709 Quadrilateral2D8<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1713 return shape_functions_local_gradients;
1734 std::istream& rIStream,
1741 std::ostream& rOStream,
1745 rOStream << std::endl;
1750 template<
class TPo
intType>
1751 const GeometryData Quadrilateral2D8<TPointType>::msGeometryData(
1754 Quadrilateral2D8<TPointType>::AllIntegrationPoints(),
1755 Quadrilateral2D8<TPointType>::AllShapeFunctionsValues(),
1756 AllShapeFunctionsLocalGradients()
1759 template<
class TPo
intType>
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
@ Kratos_Quadrilateral2D8
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
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
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
SizeType WorkingSpaceDimension() const
Definition: geometry.h:1287
TPointType const & GetPoint(const int Index) const
Definition: geometry.h:1816
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
SizeType IntegrationPointsNumber() const
Definition: geometry.h:2257
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
void clear()
Definition: amatrix_interface.h:284
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
CoordinatesArrayType const & Coordinates() const
Definition: point.h:215
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 2D quadrilateral geometry with quadratic shape functions.
Definition: quadrilateral_2d_8.h:49
double Length() const override
Definition: quadrilateral_2d_8.h:403
double DomainSize() const override
Definition: quadrilateral_2d_8.h:452
virtual Matrix & ShapeFunctionsGradients(Matrix &rResult, PointType &rPoint)
Definition: quadrilateral_2d_8.h:1283
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: quadrilateral_2d_8.h:1252
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: quadrilateral_2d_8.h:287
BaseType::JacobiansType JacobiansType
Definition: quadrilateral_2d_8.h:153
void PrintInfo(std::ostream &rOStream) const override
Definition: quadrilateral_2d_8.h:1132
Quadrilateral2D8 & operator=(const Quadrilateral2D8 &rOther)
Definition: quadrilateral_2d_8.h:312
BaseType::Pointer Create(const IndexType NewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_8.h:348
~Quadrilateral2D8() override
Definition: quadrilateral_2d_8.h:285
std::string Info() const override
Definition: quadrilateral_2d_8.h:1121
ShapeFunctionsThirdDerivativesType & ShapeFunctionsThirdDerivatives(ShapeFunctionsThirdDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:1376
Quadrilateral2D8(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)
Definition: quadrilateral_2d_8.h:210
BaseType::IntegrationPointType IntegrationPointType
Definition: quadrilateral_2d_8.h:117
Line2D3< TPointType > EdgeType
Definition: quadrilateral_2d_8.h:62
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients(IntegrationMethod ThisMethod)
Definition: quadrilateral_2d_8.h:1164
Quadrilateral2D8(Quadrilateral2D8 const &rOther)
Definition: quadrilateral_2d_8.h:260
BaseType::PointsArrayType PointsArrayType
Definition: quadrilateral_2d_8.h:110
BaseType::SizeType SizeType
Definition: quadrilateral_2d_8.h:104
Quadrilateral2D8(Quadrilateral2D8< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_8.h:277
BaseType::IndexType IndexType
Definition: quadrilateral_2d_8.h:97
Matrix & Jacobian(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:670
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: quadrilateral_2d_8.h:146
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: quadrilateral_2d_8.h:1014
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: quadrilateral_2d_8.h:292
double DeterminantOfJacobian(IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:743
void PrintData(std::ostream &rOStream) const override
Definition: quadrilateral_2d_8.h:1146
SizeType PointsNumberInDirection(IndexType LocalDirectionIndex) const override
Returns number of points per direction.
Definition: quadrilateral_2d_8.h:377
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:1209
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:968
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: quadrilateral_2d_8.h:942
double Area() const override
Definition: quadrilateral_2d_8.h:422
double DeterminantOfJacobian(const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:767
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: quadrilateral_2d_8.h:465
KRATOS_CLASS_POINTER_DEFINITION(Quadrilateral2D8)
TPointType PointType
Definition: quadrilateral_2d_8.h:83
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: quadrilateral_2d_8.h:187
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: quadrilateral_2d_8.h:160
Geometry< TPointType > BaseType
Definition: quadrilateral_2d_8.h:57
JacobiansType & Jacobian(JacobiansType &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:505
JacobiansType & Jacobian(JacobiansType &rResult, IntegrationMethod ThisMethod, Matrix &DeltaPosition) const override
Definition: quadrilateral_2d_8.h:563
Quadrilateral2D8 & operator=(Quadrilateral2D8< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_8.h:331
Quadrilateral2D8(const PointType &Point1, const PointType &Point2, const PointType &Point3, const PointType &Point4, const PointType &Point5, const PointType &Point6, const PointType &Point7, const PointType &Point8)
Definition: quadrilateral_2d_8.h:193
void ShapeFunctionsIntegrationPointsGradients(ShapeFunctionsGradientsType &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:1064
BaseType::NormalType NormalType
Definition: quadrilateral_2d_8.h:182
ShapeFunctionsSecondDerivativesType & ShapeFunctionsSecondDerivatives(ShapeFunctionsSecondDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:1301
JacobiansType & InverseOfJacobian(JacobiansType &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:794
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients()
Definition: quadrilateral_2d_8.h:1185
Quadrilateral2D8(const PointsArrayType &ThisPoints)
Definition: quadrilateral_2d_8.h:226
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: quadrilateral_2d_8.h:929
GeometryData::IntegrationMethod IntegrationMethod
Definition: quadrilateral_2d_8.h:72
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: quadrilateral_2d_8.h:132
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: quadrilateral_2d_8.h:139
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_8.h:362
Matrix & InverseOfJacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:839
BaseType::GeometriesArrayType GeometriesArrayType
Definition: quadrilateral_2d_8.h:78
Quadrilateral2D8(const std::string &rGeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: quadrilateral_2d_8.h:243
BaseType::ShapeFunctionsThirdDerivativesType ShapeFunctionsThirdDerivativesType
Definition: quadrilateral_2d_8.h:177
Vector & DeterminantOfJacobian(Vector &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:705
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: quadrilateral_2d_8.h:125
BaseType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: quadrilateral_2d_8.h:168
Matrix & Jacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_2d_8.h:624
Quadrilateral2D8(const IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: quadrilateral_2d_8.h:234
Matrix & InverseOfJacobian(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_8.h:884
friend class Quadrilateral2D8
Definition: quadrilateral_2d_8.h:1720
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
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
const GeometryData Quadrilateral2D8< TPointType >::msGeometryData & msGeometryDimension(), Quadrilateral2D8< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
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
AMatrix::MatrixRow< const TExpressionType > row(AMatrix::MatrixExpression< TExpressionType, TCategory > const &TheExpression, std::size_t RowIndex)
Definition: amatrix_interface.h:649
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
int j
Definition: quadrature.py:648
float temp
Definition: rotating_cone.py:85
J
Definition: sensitivityMatrix.py:58
integer i
Definition: TensorModule.f:17