KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
Functions
Kratos::MortarUtilities Namespace Reference

Typedefs

Type Definitions
typedef Node NodeType
 
typedef Point PointType
 
typedef PointType::CoordinatesArrayType CoordinatesArrayType
 
typedef Geometry< NodeTypeGeometryType
 Definition of geometries. More...
 
typedef Geometry< PointTypeGeometryPointType
 
typedef std::size_t IndexType
 Index type definition. More...
 
typedef std::size_t SizeType
 Size type definition. More...
 
typedef std::unordered_map< IndexType, IndexTypeIntMap
 A map for integers. More...
 

Functions

array_1d< double, 3 > GaussPointUnitNormal (const Vector &rN, const GeometryType &rGeometry)
 This function calculates the r_normal in a specific GP with a given shape function. More...
 
template<>
void ResetAuxiliarValue< Variable< double > > (ModelPart &rThisModelPart)
 
template<>
void ResetAuxiliarValue< Variable< array_1d< double, 3 > > > (ModelPart &rThisModelPart)
 
template<>
const std::string GetAuxiliarVariable< Variable< double > > ()
 
template<>
const std::string GetAuxiliarVariable< Variable< array_1d< double, 3 > > > ()
 
template<>
double GetAuxiliarValue< Variable< double > > (NodeType &rThisNode, const std::size_t iSize)
 
template<>
double GetAuxiliarValue< Variable< array_1d< double, 3 > > > (NodeType &rThisNode, const std::size_t iSize)
 
Functions
bool LengthCheck (const GeometryPointType &rGeometryLine, const double Tolerance=1.0e-6)
 This functions checks if the length of the line is to short, with the potential of provoque ill condition in the dual LM formulation. More...
 
bool HeronCheck (const GeometryPointType &rGeometryTriangle)
 This functions checks if the semiperimeter is smaller than any of the sides of the triangle. More...
 
bool HeronCheck (const PointType &rPointOrig1, const PointType &rPointOrig2, const PointType &rPointOrig3)
 This functions checks if the semiperimeter is smaller than any of the sides of the triangle. More...
 
void RotatePoint (PointType &rPointToRotate, const PointType &rPointReferenceRotation, const array_1d< double, 3 > &rSlaveTangentXi, const array_1d< double, 3 > &rSlaveTangentEta, const bool Inversed)
 This function rotates to align the projected points to a parallel plane to XY. More...
 
void ComputeNodesMeanNormalModelPart (ModelPart &rModelPart, const bool ComputeConditions=true)
 It computes the mean of the normal in the condition in all the nodes. More...
 
void ComputeNodesTangentModelPart (ModelPart &rModelPart, const Variable< array_1d< double, 3 >> *pSlipVariable=NULL, const double SlipCoefficient=1.0, const bool SlipAlways=false)
 It computes the tangent in all the nodes of the model part. More...
 
void ComputeNodesTangentFromNormalModelPart (ModelPart &rModelPart)
 It computes the tangent in all the nodes of the model part from its normal. More...
 
void ComputeTangentsFromNormal (NodeType &rNode, const array_1d< double, 3 > &rNormal, const std::size_t Dimension=3)
 It computes the tangent on the given node using the normal provided. More...
 
void ComputeTangentNodeWithLMAndSlip (NodeType &rNode, const std::size_t StepLM=0, const Variable< array_1d< double, 3 >> *pSlipVariable=NULL, const double SlipCoefficient=1.0, const std::size_t Dimension=3)
 It computes the tangent on the given node using the LM direction and Slip direction. More...
 
void ComputeTangentNodeWithSlip (NodeType &rNode, const std::size_t StepLM=0, const Variable< array_1d< double, 3 >> *pSlipVariable=NULL, const double SlipCoefficient=1.0, const std::size_t Dimension=3)
 It computes the tangent on the given node using the Slip direction. More...
 
template<typename TType >
std::vector< std::size_t > SortIndexes (const std::vector< TType > &rThisVector)
 This function gives you the indexes needed to order a vector. More...
 
template<class TContainerType >
void InvertNormalForFlag (TContainerType &rContainer, const Flags Flag)
 It inverts the order of the nodes in the conditions of a model part in order to invert the normal when certain flag is active. More...
 
template<class TContainerType >
void InvertNormal (TContainerType &rContainer)
 It inverts the order of the nodes in the conditions of a model part in order to invert the normal. More...
 
template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix< double, TNumNodes, TDim > GetCoordinates (const GeometryType &rGeometry, const bool Current=true, const IndexType Step=0)
 It calculates the matrix of coordinates of a geometry. More...
 
template<SizeType TNumNodes, SizeType TDim>
BoundedMatrix< double, TNumNodes, TDim > ComputeTangentMatrix (const GeometryType &rGeometry)
 It calculates the matrix containing the tangent vector TANGENT_XI. More...
 
template<SizeType TNumNodes, class TVarType = Variable<double>>
array_1d< double, TNumNodes > GetVariableVector (const GeometryType &rGeometry, const TVarType &rVariable, const IndexType Step)
 It calculates the vector of an historical variable of a geometry. More...
 
template<SizeType TNumNodes, class TVarType = Variable<double>>
BoundedMatrix< double, TNumNodes, 1 > GetVariableVectorMatrix (const GeometryType &rGeometry, const TVarType &rVariable, const unsigned int Step)
 It calculates the vector of an historical variable of a geometry. More...
 
template<SizeType TNumNodes, class TVarType = Variable<double>>
array_1d< double, TNumNodes > GetVariableVector (const GeometryType &rGeometry, const TVarType &rVariable)
 It calculates the vector of a non-historical variable of a geometry. More...
 
template<SizeType TNumNodes, class TVarType = Variable<double>>
BoundedMatrix< double, TNumNodes, 1 > GetVariableVectorMatrix (const GeometryType &rGeometry, const TVarType &rVariable)
 It calculates the vector of a non-historical variable of a geometry. More...
 
template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix< double, TNumNodes, TDim > GetVariableMatrix (const GeometryType &rGeometry, const Variable< array_1d< double, 3 > > &rVariable, const unsigned int Step)
 It calculates the matrix of a variable of a geometry. More...
 
template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix< double, TNumNodes, TDim > GetVariableMatrix (const GeometryType &rGeometry, const Variable< array_1d< double, 3 > > &rVariable)
 It calculates the matrix of a non-historical variable of a geometry. More...
 
template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix< double, TNumNodes, TDim > GetAbsMatrix (const BoundedMatrix< double, TNumNodes, TDim > &rInputMatrix)
 It calculates the matrix containing the absolute value of another matrix. More...
 
template<SizeType TDim, class TVarType >
unsigned int SizeToCompute ()
 This method gives the size to be computed. More...
 
template<class TVarType , bool THistorical>
void ResetValue (ModelPart &rThisModelPart, const TVarType &rThisVariable)
 This method resets the value. More...
 
template<class TVarType >
void ResetAuxiliarValue (ModelPart &rThisModelPart)
 This method resets the auxiliar value. More...
 
template<class TVarType >
const std::string GetAuxiliarVariable ()
 This method returns the auxiliar variable. More...
 
template<class TVarType >
double GetAuxiliarValue (NodeType &rThisNode, const std::size_t iSize)
 This method returns the auxiliar variable. More...
 
template<class TVarType , bool THistorical>
void MatrixValue (const GeometryType &rThisGeometry, const TVarType &rThisVariable, Matrix &rThisValue)
 This method adds the value. More...
 
template<class TVarType , bool THistorical>
void AddValue (GeometryType &rThisGeometry, const TVarType &rThisVariable, const Matrix &rThisValue)
 This method adds the value. More...
 
template<class TVarType , bool THistorical>
void AddAreaWeightedNodalValue (NodeType &rThisNode, const TVarType &rThisVariable, const double RefArea=1.0, const double Tolerance=1.0e-4)
 This method adds the value. More...
 
template<class TVarType , bool THistorical>
void UpdateDatabase (ModelPart &rThisModelPart, const TVarType &rThisVariable, Vector &rDx, const std::size_t Index, IntMap &rConectivityDatabase)
 This method updates the database in the amster side. More...
 

Typedef Documentation

◆ CoordinatesArrayType

◆ GeometryPointType

◆ GeometryType

Definition of geometries.

◆ IndexType

Index type definition.

◆ IntMap

typedef std::unordered_map<IndexType, IndexType> Kratos::MortarUtilities::IntMap

A map for integers.

◆ NodeType

◆ PointType

◆ SizeType

typedef std::size_t Kratos::MortarUtilities::SizeType

Size type definition.

Function Documentation

◆ AddAreaWeightedNodalValue()

template<class TVarType , bool THistorical>
void Kratos::MortarUtilities::AddAreaWeightedNodalValue ( NodeType rThisNode,
const TVarType &  rThisVariable,
const double  RefArea = 1.0,
const double  Tolerance = 1.0e-4 
)

This method adds the value.

Parameters
rThisNodeThe node to update
rThisVariableThe variable to set

◆ AddValue()

template<class TVarType , bool THistorical>
void Kratos::MortarUtilities::AddValue ( GeometryType rThisGeometry,
const TVarType &  rThisVariable,
const Matrix rThisValue 
)

This method adds the value.

Warning
This operation is not threadsafe
Parameters
rThisGeometryThe geometrty to update
rThisVariableThe variable to set
rThisValueThe matrix to be updated

◆ ComputeNodesMeanNormalModelPart()

void Kratos::MortarUtilities::ComputeNodesMeanNormalModelPart ( ModelPart rModelPart,
const bool  ComputeConditions = true 
)

It computes the mean of the normal in the condition in all the nodes.

Parameters
rModelPartThe model part to compute
ComputeConditionsIf computed over conditions or elements

◆ ComputeNodesTangentFromNormalModelPart()

void Kratos::MortarUtilities::ComputeNodesTangentFromNormalModelPart ( ModelPart rModelPart)

It computes the tangent in all the nodes of the model part from its normal.

Parameters
rModelPartThe model part to compute

◆ ComputeNodesTangentModelPart()

void Kratos::MortarUtilities::ComputeNodesTangentModelPart ( ModelPart rModelPart,
const Variable< array_1d< double, 3 >> *  pSlipVariable = NULL,
const double  SlipCoefficient = 1.0,
const bool  SlipAlways = false 
)

It computes the tangent in all the nodes of the model part.

Parameters
rModelPartThe model part to compute
pSlipVariableThe pointer to the slip variable
SlipCoefficientThe slip contribution
SlipAlwaysUses the slip even in case that LM are available

◆ ComputeTangentMatrix()

template<SizeType TNumNodes, SizeType TDim>
BoundedMatrix<double, TNumNodes, TDim> Kratos::MortarUtilities::ComputeTangentMatrix ( const GeometryType rGeometry)

It calculates the matrix containing the tangent vector TANGENT_XI.

Parameters
rGeometryThe geometry to calculate
Returns
tangent_matrix The matrix containing the tangent vectors of the LM

◆ ComputeTangentNodeWithLMAndSlip()

void Kratos::MortarUtilities::ComputeTangentNodeWithLMAndSlip ( NodeType rNode,
const std::size_t  StepLM = 0,
const Variable< array_1d< double, 3 >> *  pSlipVariable = NULL,
const double  SlipCoefficient = 1.0,
const std::size_t  Dimension = 3 
)

It computes the tangent on the given node using the LM direction and Slip direction.

Parameters
rNodeThe node where to compute the tangent
StepLMThe considered step slip
pSlipVariableThe pointer to the slip variable
SlipCoefficientThe slip contribution
DimensionThe current working dimension

◆ ComputeTangentNodeWithSlip()

void Kratos::MortarUtilities::ComputeTangentNodeWithSlip ( NodeType rNode,
const std::size_t  StepLM = 0,
const Variable< array_1d< double, 3 >> *  pSlipVariable = NULL,
const double  SlipCoefficient = 1.0,
const std::size_t  Dimension = 3 
)

It computes the tangent on the given node using the Slip direction.

Parameters
rNodeThe node where to compute the tangent
StepLMThe considered step slip
pSlipVariableThe pointer to the slip variable
SlipCoefficientThe slip contribution
DimensionThe current working dimension

◆ ComputeTangentsFromNormal()

void Kratos::MortarUtilities::ComputeTangentsFromNormal ( NodeType rNode,
const array_1d< double, 3 > &  rNormal,
const std::size_t  Dimension = 3 
)

It computes the tangent on the given node using the normal provided.

Parameters
rNodeThe node where to compute the tangent
rNormalThe normal vector
DimensionThe current working dimension

◆ GaussPointUnitNormal()

array_1d<double,3> Kratos::MortarUtilities::GaussPointUnitNormal ( const Vector rN,
const GeometryType rGeometry 
)

This function calculates the r_normal in a specific GP with a given shape function.

Parameters
rNThe shape function considered
rGeometryThe geometry of condition of interest
Returns
The r_normal in the GP

◆ GetAbsMatrix()

template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix<double, TNumNodes, TDim> Kratos::MortarUtilities::GetAbsMatrix ( const BoundedMatrix< double, TNumNodes, TDim > &  rInputMatrix)

It calculates the matrix containing the absolute value of another matrix.

Parameters
rInputMatrixThe original matrix
Returns
AbsMatrix The matrix containing the absolute value of another matrix

◆ GetAuxiliarValue()

template<class TVarType >
double Kratos::MortarUtilities::GetAuxiliarValue ( NodeType rThisNode,
const std::size_t  iSize 
)

This method returns the auxiliar variable.

Parameters
rThisNodeReference to the node of interest
iSizeThe Index of the component
Returns
The value of the auxiliar variable

◆ GetAuxiliarValue< Variable< array_1d< double, 3 > > >()

template<>
double Kratos::MortarUtilities::GetAuxiliarValue< Variable< array_1d< double, 3 > > > ( NodeType rThisNode,
const std::size_t  iSize 
)

◆ GetAuxiliarValue< Variable< double > >()

template<>
double Kratos::MortarUtilities::GetAuxiliarValue< Variable< double > > ( NodeType rThisNode,
const std::size_t  iSize 
)

◆ GetAuxiliarVariable()

template<class TVarType >
const std::string Kratos::MortarUtilities::GetAuxiliarVariable ( )

This method returns the auxiliar variable.

Returns
The auxiliar variable

◆ GetAuxiliarVariable< Variable< array_1d< double, 3 > > >()

template<>
const std::string Kratos::MortarUtilities::GetAuxiliarVariable< Variable< array_1d< double, 3 > > > ( )

◆ GetAuxiliarVariable< Variable< double > >()

template<>
const std::string Kratos::MortarUtilities::GetAuxiliarVariable< Variable< double > > ( )

◆ GetCoordinates()

template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix<double, TNumNodes, TDim> Kratos::MortarUtilities::GetCoordinates ( const GeometryType rGeometry,
const bool  Current = true,
const IndexType  Step = 0 
)

It calculates the matrix of coordinates of a geometry.

Parameters
rGeometryThe geometry to calculate
CurrentIf we calculate the Current coordinates or the initial ones
StepThe time step where it is computed
Returns
coordinates The matrix containing the coordinates of the geometry

◆ GetVariableMatrix() [1/2]

template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix<double, TNumNodes, TDim> Kratos::MortarUtilities::GetVariableMatrix ( const GeometryType rGeometry,
const Variable< array_1d< double, 3 > > &  rVariable 
)

It calculates the matrix of a non-historical variable of a geometry.

Parameters
rGeometryThe geometry to calculate
rVariableThe name of the variable to calculate
Returns
var_matrix The matrix containing the variables of the geometry

◆ GetVariableMatrix() [2/2]

template<SizeType TDim, SizeType TNumNodes>
BoundedMatrix<double, TNumNodes, TDim> Kratos::MortarUtilities::GetVariableMatrix ( const GeometryType rGeometry,
const Variable< array_1d< double, 3 > > &  rVariable,
const unsigned int  Step 
)

It calculates the matrix of a variable of a geometry.

Parameters
rGeometryThe geometry to calculate
rVariableThe name of the variable to calculate
StepThe step where it is computed
Returns
var_matrix The matrix containing the variables of the geometry

◆ GetVariableVector() [1/2]

template<SizeType TNumNodes, class TVarType = Variable<double>>
array_1d<double, TNumNodes> Kratos::MortarUtilities::GetVariableVector ( const GeometryType rGeometry,
const TVarType &  rVariable 
)

It calculates the vector of a non-historical variable of a geometry.

Parameters
rGeometryThe geometry to calculate
rVariableThe name of the variable to calculate
Returns
var_vector The vector containing the variables of the geometry

◆ GetVariableVector() [2/2]

template<SizeType TNumNodes, class TVarType = Variable<double>>
array_1d<double, TNumNodes> Kratos::MortarUtilities::GetVariableVector ( const GeometryType rGeometry,
const TVarType &  rVariable,
const IndexType  Step 
)

It calculates the vector of an historical variable of a geometry.

Parameters
rGeometryThe geometry to calculate
rVariableThe name of the variable to calculate
StepThe step where it is computed
Returns
var_vector The vector containing the variables of the geometry

◆ GetVariableVectorMatrix() [1/2]

template<SizeType TNumNodes, class TVarType = Variable<double>>
BoundedMatrix<double, TNumNodes, 1> Kratos::MortarUtilities::GetVariableVectorMatrix ( const GeometryType rGeometry,
const TVarType &  rVariable 
)

It calculates the vector of a non-historical variable of a geometry.

Parameters
rGeometryThe geometry to calculate
rVariableThe name of the variable to calculate
Returns
var_vector The vector containing the variables of the geometry

◆ GetVariableVectorMatrix() [2/2]

template<SizeType TNumNodes, class TVarType = Variable<double>>
BoundedMatrix<double, TNumNodes, 1> Kratos::MortarUtilities::GetVariableVectorMatrix ( const GeometryType rGeometry,
const TVarType &  rVariable,
const unsigned int  Step 
)

It calculates the vector of an historical variable of a geometry.

Parameters
rGeometryThe geometry to calculate
rVariableThe name of the variable to calculate
StepThe step where it is computed
Returns
var_vector The vector containing the variables of the geometry

◆ HeronCheck() [1/2]

bool Kratos::MortarUtilities::HeronCheck ( const GeometryPointType rGeometryTriangle)

This functions checks if the semiperimeter is smaller than any of the sides of the triangle.

Parameters
rGeometryTriangleThe triangle to be checked
Returns
True if the triangle is in bad shape, false otherwise

◆ HeronCheck() [2/2]

bool Kratos::MortarUtilities::HeronCheck ( const PointType rPointOrig1,
const PointType rPointOrig2,
const PointType rPointOrig3 
)

This functions checks if the semiperimeter is smaller than any of the sides of the triangle.

Parameters
rPointOrig1The triangle first point
rPointOrig2The triangle second point
rPointOrig3The triangle third point
Returns
True if the triangle is in bad shape, false otherwise

◆ InvertNormal()

template<class TContainerType >
void Kratos::MortarUtilities::InvertNormal ( TContainerType &  rContainer)

It inverts the order of the nodes in the conditions of a model part in order to invert the normal.

Parameters
rContainerReference to the objective container

◆ InvertNormalForFlag()

template<class TContainerType >
void Kratos::MortarUtilities::InvertNormalForFlag ( TContainerType &  rContainer,
const Flags  Flag 
)

It inverts the order of the nodes in the conditions of a model part in order to invert the normal when certain flag is active.

Parameters
rContainerReference to the objective container
FlagThe flag of the entities inverted

◆ LengthCheck()

bool Kratos::MortarUtilities::LengthCheck ( const GeometryPointType rGeometryLine,
const double  Tolerance = 1.0e-6 
)

This functions checks if the length of the line is to short, with the potential of provoque ill condition in the dual LM formulation.

Parameters
rGeometryLineThe line to be checked
ToleranceThe threshold length
Returns
True if the line is too short, false otherwise

◆ MatrixValue()

template<class TVarType , bool THistorical>
void Kratos::MortarUtilities::MatrixValue ( const GeometryType rThisGeometry,
const TVarType &  rThisVariable,
Matrix rThisValue 
)

This method adds the value.

Parameters
rThisGeometryThe geometrty to update
rThisVariableThe variable to set
rThisValueThe matrix to be updated

◆ ResetAuxiliarValue()

template<class TVarType >
void Kratos::MortarUtilities::ResetAuxiliarValue ( ModelPart rThisModelPart)

This method resets the auxiliar value.

Parameters
rThisModelPartThe model part to update

◆ ResetAuxiliarValue< Variable< array_1d< double, 3 > > >()

template<>
void Kratos::MortarUtilities::ResetAuxiliarValue< Variable< array_1d< double, 3 > > > ( ModelPart rThisModelPart)

◆ ResetAuxiliarValue< Variable< double > >()

template<>
void Kratos::MortarUtilities::ResetAuxiliarValue< Variable< double > > ( ModelPart rThisModelPart)

◆ ResetValue()

template<class TVarType , bool THistorical>
void Kratos::MortarUtilities::ResetValue ( ModelPart rThisModelPart,
const TVarType &  rThisVariable 
)

This method resets the value.

Parameters
rThisModelPartThe model part to update
rThisVariableThe variable to set

◆ RotatePoint()

void Kratos::MortarUtilities::RotatePoint ( PointType rPointToRotate,
const PointType rPointReferenceRotation,
const array_1d< double, 3 > &  rSlaveTangentXi,
const array_1d< double, 3 > &  rSlaveTangentEta,
const bool  Inversed 
)

This function rotates to align the projected points to a parallel plane to XY.

Parameters
rPointToRotateThe points from the origin geometry and the the point rotated
rPointReferenceRotationThe center point used as reference to rotate
rSlaveTangentXiThe first tangent vector of the slave condition
rSlaveTangentEtaThe second tangent vector of the slave condition
InversedIf we rotate to the XY or we recover from XY

◆ SizeToCompute()

template<SizeType TDim, class TVarType >
unsigned int Kratos::MortarUtilities::SizeToCompute ( )

This method gives the size to be computed.

◆ SortIndexes()

template<typename TType >
std::vector<std::size_t> Kratos::MortarUtilities::SortIndexes ( const std::vector< TType > &  rThisVector)

This function gives you the indexes needed to order a vector.

Parameters
rThisVectorThe vector to order
Returns
idx The vector of indexes

◆ UpdateDatabase()

template<class TVarType , bool THistorical>
void Kratos::MortarUtilities::UpdateDatabase ( ModelPart rThisModelPart,
const TVarType &  rThisVariable,
Vector rDx,
const std::size_t  Index,
IntMap rConectivityDatabase 
)

This method updates the database in the amster side.

Parameters
rThisModelPartThe model part
rThisVariableThe variable to set
rDxThe vector with the increment of the value
IndexThe index used in the case of a vector variable
rConectivityDatabaseThe database that will be used to assemble the system