13 #if !defined(KRATOS_QUADRILATERAL_INTERFACE_2D_4_H_INCLUDED )
14 #define KRATOS_QUADRILATERAL_INTERFACE_2D_4_H_INCLUDED
220 typename PointType::Pointer pSecondPoint,
221 typename PointType::Pointer pThirdPoint,
222 typename PointType::Pointer pFourthPoint )
225 double p2_p1_X = (pThirdPoint->X() + pFourthPoint->X() - pFirstPoint->X() - pSecondPoint->X())*0.5;
226 double p2_p1_Y = (pThirdPoint->Y() + pFourthPoint->Y() - pFirstPoint->Y() - pSecondPoint->Y())*0.5;
227 double p4_p3_X = (pThirdPoint->X() + pSecondPoint->X() - pFirstPoint->X() - pFourthPoint->X())*0.5;
228 double p4_p3_Y = (pThirdPoint->Y() + pSecondPoint->Y() - pFirstPoint->Y() - pFourthPoint->Y())*0.5;
230 double lx = std::sqrt(p4_p3_X*p4_p3_X + p4_p3_Y*p4_p3_Y);
231 double ly = std::sqrt(p2_p1_X*p2_p1_X + p2_p1_Y*p2_p1_Y);
250 :
BaseType( ThisPoints, &msGeometryData )
334 template<
class TOtherPo
intType>
358 rResult( 0, 0 ) = -1.0;
359 rResult( 0, 1 ) = -1.0;
360 rResult( 1, 0 ) = 1.0;
361 rResult( 1, 1 ) = -1.0;
362 rResult( 2, 0 ) = 1.0;
363 rResult( 2, 1 ) = 1.0;
364 rResult( 3, 0 ) = -1.0;
365 rResult( 3, 1 ) = 1.0;
398 return std::sqrt(
v[0]*
v[0] +
v[1]*
v[1] );
452 const double Tolerance = std::numeric_limits<double>::epsilon()
457 if ( std::abs(rResult[0]) <= (1.0+Tolerance) )
459 if ( std::abs(rResult[1]) <= (1.0+Tolerance) )
492 Normal[0] = SecondPoint[1] - FirstPoint[1];
493 Normal[1] = FirstPoint[0] - SecondPoint[0];
499 VectorPoint[0] = rPoint[0] - FirstPoint[0];
500 VectorPoint[1] = rPoint[1] - FirstPoint[1];
501 double dist_proy = VectorPoint[0] *
Normal[0] + VectorPoint[1] *
Normal[1];
507 double l1 = (rPoint[0] - FirstPoint[0]) * (rPoint[0] - FirstPoint[0])
508 + (rPoint[1] - FirstPoint[1]) * (rPoint[1] - FirstPoint[1]);
511 double l2 = (rPoint[0] - SecondPoint[0]) * (rPoint[0] - SecondPoint[0])
512 + (rPoint[1] - SecondPoint[1]) * (rPoint[1] - SecondPoint[1]);
517 if (l1 <= (
L +
tol) && l2 <= (
L +
tol))
519 rResult[0] = 2.0 * l1/(
L +
tol) - 1.0;
561 if(rResult.size1() != 2 || rResult.size2() != 1)
562 rResult.
resize(2, 1,
false);
563 rResult( 0, 0 ) = ( p1[0] - p0[0] ) * 0.5;
564 rResult( 1, 0 ) = ( p1[1] - p0[1] ) * 0.5;
591 const Matrix& DeltaPosition )
const override
598 dp0[0] = (DeltaPosition(0,0) + DeltaPosition(3,0))*0.5;
599 dp0[1] = (DeltaPosition(0,1) + DeltaPosition(3,1))*0.5;
600 dp1[0] = (DeltaPosition(1,0) + DeltaPosition(2,0))*0.5;
601 dp1[1] = (DeltaPosition(1,1) + DeltaPosition(2,1))*0.5;
603 if(rResult.size1() != 2 || rResult.size2() != 1)
604 rResult.
resize(2, 1,
false);
605 rResult( 0, 0 ) = ( (p1[0]-dp1[0]) - (p0[0]-dp0[0]) ) * 0.5;
606 rResult( 1, 0 ) = ( (p1[1]-dp1[1]) - (p0[1]-dp0[1]) ) * 0.5;
630 if(rResult.size1() != 2 || rResult.size2() != 1)
631 rResult.
resize(2, 1,
false);
632 rResult( 0, 0 ) = ( p1[0] - p0[0] ) * 0.5;
633 rResult( 1, 0 ) = ( p1[1] - p0[1] ) * 0.5;
656 if ( rResult.size() != this->IntegrationPointsNumber( ThisMethod ) )
660 double detJ =
Length()*0.5;
864 switch ( ShapeFunctionIndex )
867 return( 0.25*( 1.0 - rPoint[0] )*( 1.0 - rPoint[1] ) );
869 return( 0.25*( 1.0 + rPoint[0] )*( 1.0 - rPoint[1] ) );
871 return( 0.25*( 1.0 + rPoint[0] )*( 1.0 + rPoint[1] ) );
873 return( 0.25*( 1.0 - rPoint[0] )*( 1.0 + rPoint[1] ) );
875 KRATOS_ERROR <<
"Wrong index of shape function!" << *
this << std::endl;
894 if(rResult.size() != 4) rResult.
resize(4,
false);
895 rResult[0] = 0.25*( 1.0 - rCoordinates[0] )*( 1.0 - rCoordinates[1] );
896 rResult[1] = 0.25*( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] );
897 rResult[2] = 0.25*( 1.0 + rCoordinates[0] )*( 1.0 + rCoordinates[1] );
898 rResult[3] = 0.25*( 1.0 - rCoordinates[0] )*( 1.0 + rCoordinates[1] );
923 const unsigned int integration_points_number =
926 if ( integration_points_number == 0 )
927 KRATOS_ERROR <<
"This integration method is not supported" << *
this << std::endl;
930 if ( rResult.size() != integration_points_number )
940 CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
948 for (
unsigned int pnt = 0; pnt < integration_points_number; pnt++ )
950 rResult[pnt].
resize( 4, 2 ,
false);
952 for (
int i = 0;
i < 4;
i++ )
954 for (
int j = 0;
j < 2;
j++ )
956 rResult[pnt](
i,
j ) =
957 ( locG[pnt](
i, 0 ) * invJ[pnt](
j, 0 ) )
958 + ( locG[pnt](
i, 1 ) * invJ[pnt](
j, 1 ) );
975 std::string
Info()
const override
977 return "2 dimensional quadrilateral with four nodes in 2D space";
988 rOStream <<
"2 dimensional quadrilateral with four nodes in 2D space";
1009 std::cout << std::endl;
1015 rOStream <<
" Jacobian in the origin\t : " << jacobian;
1030 rResult.
resize( 4, 2 ,
false);
1032 rResult( 0, 0 ) = -0.25 * ( 1.0 - rPoint[1] );
1033 rResult( 0, 1 ) = -0.25 * ( 1.0 - rPoint[0] );
1034 rResult( 1, 0 ) = 0.25 * ( 1.0 - rPoint[1] );
1035 rResult( 1, 1 ) = -0.25 * ( 1.0 + rPoint[0] );
1036 rResult( 2, 0 ) = 0.25 * ( 1.0 + rPoint[1] );
1037 rResult( 2, 1 ) = 0.25 * ( 1.0 + rPoint[0] );
1038 rResult( 3, 0 ) = -0.25 * ( 1.0 + rPoint[1] );
1039 rResult( 3, 1 ) = 0.25 * ( 1.0 - rPoint[0] );
1052 if ( rResult.size() != this->PointsNumber() )
1060 rResult[0].
resize( 2, 2 ,
false);
1062 rResult[1].
resize( 2, 2 ,
false);
1063 rResult[2].
resize( 2, 2 ,
false);
1064 rResult[3].
resize( 2, 2 ,
false);
1065 rResult[0]( 0, 0 ) = 0.0;
1066 rResult[0]( 0, 1 ) = 0.25;
1067 rResult[0]( 1, 0 ) = 0.25;
1068 rResult[0]( 1, 1 ) = 0.0;
1069 rResult[1]( 0, 0 ) = 0.0;
1070 rResult[1]( 0, 1 ) = -0.25;
1071 rResult[1]( 1, 0 ) = -0.25;
1072 rResult[1]( 1, 1 ) = 0.0;
1073 rResult[2]( 0, 0 ) = 0.0;
1074 rResult[2]( 0, 1 ) = 0.25;
1075 rResult[2]( 1, 0 ) = 0.25;
1076 rResult[2]( 1, 1 ) = 0.0;
1077 rResult[3]( 0, 0 ) = 0.0;
1078 rResult[3]( 0, 1 ) = -0.25;
1079 rResult[3]( 1, 0 ) = -0.25;
1080 rResult[3]( 1, 1 ) = 0.0;
1093 if ( rResult.size() != this->PointsNumber() )
1110 for (
unsigned int j = 0;
j < 2;
j++ )
1112 rResult[
i][
j].
resize( 2, 2,
false );
1126 for (
int i = 0;
i < 4;
i++ )
1128 rResult[
i][0]( 0, 0 ) = 0.0;
1129 rResult[
i][0]( 0, 1 ) = 0.0;
1130 rResult[
i][0]( 1, 0 ) = 0.0;
1131 rResult[
i][0]( 1, 1 ) = 0.0;
1132 rResult[
i][1]( 0, 0 ) = 0.0;
1133 rResult[
i][1]( 0, 1 ) = 0.0;
1134 rResult[
i][1]( 1, 0 ) = 0.0;
1135 rResult[
i][1]( 1, 1 ) = 0.0;
1171 void save(
Serializer& rSerializer )
const override
1203 static Matrix CalculateShapeFunctionsIntegrationPointsValues(
1207 AllIntegrationPoints();
1210 const int integration_points_number = integration_points.size();
1212 const int points_number = 4;
1214 Matrix shape_function_values( integration_points_number, points_number );
1217 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1219 shape_function_values( pnt, 0 ) =
1220 0.25 * ( 1.0 - integration_points[pnt].X() )
1221 * ( 1.0 - integration_points[pnt].
Y() );
1222 shape_function_values( pnt, 1 ) =
1223 0.25 * ( 1.0 + integration_points[pnt].X() )
1224 * ( 1.0 - integration_points[pnt].
Y() );
1225 shape_function_values( pnt, 2 ) =
1226 0.25 * ( 1.0 + integration_points[pnt].X() )
1227 * ( 1.0 + integration_points[pnt].
Y() );
1228 shape_function_values( pnt, 3 ) =
1229 0.25 * ( 1.0 - integration_points[pnt].X() )
1230 * ( 1.0 + integration_points[pnt].
Y() );
1232 return shape_function_values;
1251 AllIntegrationPoints();
1254 const int integration_points_number = integration_points.size();
1260 for (
int pnt = 0; pnt < integration_points_number; pnt++ )
1263 result( 0, 0 ) = -0.25 * ( 1.0 - integration_points[pnt].Y() );
1264 result( 0, 1 ) = -0.25 * ( 1.0 - integration_points[pnt].X() );
1265 result( 1, 0 ) = 0.25 * ( 1.0 - integration_points[pnt].Y() );
1266 result( 1, 1 ) = -0.25 * ( 1.0 + integration_points[pnt].X() );
1267 result( 2, 0 ) = 0.25 * ( 1.0 + integration_points[pnt].Y() );
1268 result( 2, 1 ) = 0.25 * ( 1.0 + integration_points[pnt].X() );
1269 result( 3, 0 ) = -0.25 * ( 1.0 + integration_points[pnt].Y() );
1270 result( 3, 1 ) = 0.25 * ( 1.0 - integration_points[pnt].X() );
1271 d_shape_f_values[pnt] = result;
1273 return d_shape_f_values;
1284 Quadrature < QuadrilateralGaussLobattoIntegrationPoints1,
1285 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1286 Quadrature < QuadrilateralGaussLobattoIntegrationPoints2,
1287 2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1292 return integration_points;
1303 QuadrilateralInterface2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1305 QuadrilateralInterface2D4<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1311 return shape_functions_values;
1318 AllShapeFunctionsLocalGradients()
1329 return shape_functions_local_gradients;
1368 std::istream& rIStream,
1374 std::ostream& rOStream,
1378 rOStream << std::endl;
1385 template<
class TPo
intType>
const
1386 GeometryData QuadrilateralInterface2D4<TPointType>::msGeometryData(
1389 QuadrilateralInterface2D4<TPointType>::AllIntegrationPoints(),
1390 QuadrilateralInterface2D4<TPointType>::AllShapeFunctionsValues(),
1391 AllShapeFunctionsLocalGradients()
1394 template<
class TPo
intType>
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
@ Kratos_Quadrilateral2D4
SizeType IntegrationPointsNumber() const
Definition: geometry_data.h:430
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
virtual array_1d< double, 3 > Normal(const CoordinatesArrayType &rPointLocalCoordinates) const
It returns a vector that is normal to its corresponding geometry in the given local point.
Definition: geometry.h:1543
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
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
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
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 two node 2D line geometry with linear shape functions.
Definition: line_2d_2.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
Definition: quadrilateral_interface_2d_4.h:54
double DeterminantOfJacobian(const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:720
BaseType::IndexType IndexType
Definition: quadrilateral_interface_2d_4.h:97
BaseType::ShapeFunctionsThirdDerivativesType ShapeFunctionsThirdDerivativesType
Definition: quadrilateral_interface_2d_4.h:179
void PrintInfo(std::ostream &rOStream) const override
Definition: quadrilateral_interface_2d_4.h:986
double Area() const override
Definition: quadrilateral_interface_2d_4.h:416
BaseType::GeometriesArrayType GeometriesArrayType
Definition: quadrilateral_interface_2d_4.h:84
Matrix & InverseOfJacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_interface_2d_4.h:776
GeometryData::IntegrationMethod IntegrationMethod
Definition: quadrilateral_interface_2d_4.h:78
CoordinatesArrayType & PointLocalCoordinates(CoordinatesArrayType &rResult, const CoordinatesArrayType &rPoint) const override
Returns the local coordinates of a given arbitrary point.
Definition: quadrilateral_interface_2d_4.h:475
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: quadrilateral_interface_2d_4.h:137
QuadrilateralInterface2D4(typename PointType::Pointer pFirstPoint, typename PointType::Pointer pSecondPoint, typename PointType::Pointer pThirdPoint, typename PointType::Pointer pFourthPoint)
Definition: quadrilateral_interface_2d_4.h:219
Vector & DeterminantOfJacobian(Vector &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_interface_2d_4.h:652
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: quadrilateral_interface_2d_4.h:143
ShapeFunctionsSecondDerivativesType & ShapeFunctionsSecondDerivatives(ShapeFunctionsSecondDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:1049
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: quadrilateral_interface_2d_4.h:355
double DomainSize() const override
Definition: quadrilateral_interface_2d_4.h:436
QuadrilateralInterface2D4 & operator=(const QuadrilateralInterface2D4 &rOther)
Definition: quadrilateral_interface_2d_4.h:317
JacobiansType & InverseOfJacobian(JacobiansType &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_interface_2d_4.h:746
BaseType::JacobiansType JacobiansType
Definition: quadrilateral_interface_2d_4.h:156
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:1028
BaseType::PointsArrayType PointsArrayType
Definition: quadrilateral_interface_2d_4.h:110
QuadrilateralInterface2D4(QuadrilateralInterface2D4 const &rOther)
Definition: quadrilateral_interface_2d_4.h:265
void ShapeFunctionsIntegrationPointsGradients(ShapeFunctionsGradientsType &rResult, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_interface_2d_4.h:919
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: quadrilateral_interface_2d_4.h:835
Matrix & InverseOfJacobian(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:800
friend class QuadrilateralInterface2D4
Definition: quadrilateral_interface_2d_4.h:1346
BaseType::IntegrationPointType IntegrationPointType
Definition: quadrilateral_interface_2d_4.h:121
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: quadrilateral_interface_2d_4.h:892
Line2D2< TPointType > EdgeType
Definition: quadrilateral_interface_2d_4.h:68
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: quadrilateral_interface_2d_4.h:822
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:861
TPointType PointType
Definition: quadrilateral_interface_2d_4.h:89
std::string Info() const override
Definition: quadrilateral_interface_2d_4.h:975
virtual ~QuadrilateralInterface2D4()
Definition: quadrilateral_interface_2d_4.h:290
void PrintData(std::ostream &rOStream) const override
Definition: quadrilateral_interface_2d_4.h:1005
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: quadrilateral_interface_2d_4.h:115
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: quadrilateral_interface_2d_4.h:163
Geometry< TPointType > BaseType
Definition: quadrilateral_interface_2d_4.h:63
double Length() const override
Definition: quadrilateral_interface_2d_4.h:392
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: quadrilateral_interface_2d_4.h:297
BaseType::NormalType NormalType
Definition: quadrilateral_interface_2d_4.h:184
Matrix & Jacobian(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:626
BaseType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: quadrilateral_interface_2d_4.h:171
ShapeFunctionsThirdDerivativesType & ShapeFunctionsThirdDerivatives(ShapeFunctionsThirdDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_interface_2d_4.h:1090
Matrix & Jacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_interface_2d_4.h:555
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: quadrilateral_interface_2d_4.h:292
BaseType::Pointer Create(PointsArrayType const &ThisPoints) const override
Creates a new geometry pointer.
Definition: quadrilateral_interface_2d_4.h:345
BaseType::SizeType SizeType
Definition: quadrilateral_interface_2d_4.h:104
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: quadrilateral_interface_2d_4.h:130
Matrix & Jacobian(Matrix &rResult, IndexType IntegrationPointIndex, IntegrationMethod ThisMethod, const Matrix &DeltaPosition) const override
Definition: quadrilateral_interface_2d_4.h:588
QuadrilateralInterface2D4(const PointsArrayType &ThisPoints)
Definition: quadrilateral_interface_2d_4.h:249
KRATOS_CLASS_POINTER_DEFINITION(QuadrilateralInterface2D4)
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: quadrilateral_interface_2d_4.h:149
QuadrilateralInterface2D4 & operator=(QuadrilateralInterface2D4< TOtherPointType > const &rOther)
Definition: quadrilateral_interface_2d_4.h:335
QuadrilateralInterface2D4(QuadrilateralInterface2D4< TOtherPointType > const &rOther)
Definition: quadrilateral_interface_2d_4.h:282
double DeterminantOfJacobian(IndexType IntegrationPointIndex, IntegrationMethod ThisMethod) const override
Definition: quadrilateral_interface_2d_4.h:689
bool IsInside(const CoordinatesArrayType &rPoint, CoordinatesArrayType &rResult, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Definition: quadrilateral_interface_2d_4.h:449
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
BOOST_UBLAS_INLINE void clear()
Definition: array_1d.h:325
#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
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
Internals::Matrix< double, AMatrix::dynamic, AMatrix::dynamic > Matrix
Definition: amatrix_interface.h:470
const GeometryData QuadrilateralInterface2D4< TPointType >::msGeometryData & msGeometryDimension(), QuadrilateralInterface2D4< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
const GeometryData BrepCurve< TContainerPointType, TContainerPointEmbeddedType >::msGeometryData & msGeometryDimension
Definition: brep_curve.h:483
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
v
Definition: generate_convection_diffusion_explicit_element.py:114
int tol
Definition: hinsberg_optimization.py:138
def load(f)
Definition: ode_solve.py:307
int L
Definition: ode_solve.py:390
int j
Definition: quadrature.py:648
float temp
Definition: rotating_cone.py:85
integer i
Definition: TensorModule.f:17
e
Definition: run_cpp_mpi_tests.py:31