v0.14.0
Public Member Functions | Private Attributes | List of all members
SimpleContactProblem::OpCalAugmentedTractionRhsSlave Struct Reference

RHS-operator for the simple contact element. More...

#include <users_modules/mortar_contact/src/SimpleContact.hpp>

Inheritance diagram for SimpleContactProblem::OpCalAugmentedTractionRhsSlave:
[legend]
Collaboration diagram for SimpleContactProblem::OpCalAugmentedTractionRhsSlave:
[legend]

Public Member Functions

 OpCalAugmentedTractionRhsSlave (const string field_name, boost::shared_ptr< CommonDataSimpleContact > common_data_contact)
 
MoFEMErrorCode doWork (int side, EntityType type, EntitiesFieldData::EntData &data)
 Integrates Lagrange multipliers virtual work on slave surface and assembles components to the RHS global vector. More...
 

Private Attributes

boost::shared_ptr< CommonDataSimpleContactcommonDataSimpleContact
 
VectorDouble vecF
 

Detailed Description

RHS-operator for the simple contact element.

Integrates Augmented Lagrange multipliers virtual work on slave surface and assembles components to the RHS global vector.

Definition at line 1104 of file SimpleContact.hpp.

Constructor & Destructor Documentation

◆ OpCalAugmentedTractionRhsSlave()

SimpleContactProblem::OpCalAugmentedTractionRhsSlave::OpCalAugmentedTractionRhsSlave ( const string  field_name,
boost::shared_ptr< CommonDataSimpleContact common_data_contact 
)
inline

Member Function Documentation

◆ doWork()

MoFEMErrorCode SimpleContactProblem::OpCalAugmentedTractionRhsSlave::doWork ( int  side,
EntityType  type,
EntitiesFieldData::EntData data 
)

Integrates Lagrange multipliers virtual work on slave surface and assembles components to the RHS global vector.

Integrates Lagrange multipliers virtual work \( \delta W_{\text c}\) on slave surface and assembles components to the RHS global vector

\[ {\delta W^{(1)}_{\text c}(\lambda, \delta \mathbf{x}^{(1)}}) \,\, = - \left\{ \begin{array}{ll} \int_{{\gamma}^{(1)}_{\text c}} (\lambda + c_{\textrm n} g_{\textrm{n}}) \mathbf{n}(\mathbf{x}^{(1)}) \cdot \delta{\mathbf{x}^{(1)}} \,\,{ {\text d} {\gamma}} & \lambda + c_{\text n} g_{\textrm{n}}\leq 0 \\ 0 & \lambda + c_{\text n} g_{\textrm{n}} > 0 \\ \end{array} \right. \]

where \({\gamma}^{(1)}_{\text c}\) is the surface integration domain of the slave surface, \( \lambda\) is the Lagrange multiplier, \( c_{\textrm n}\) is the regularisation/augmentation parameter of stress dimensions, \( g_{\textrm{n}}\) is the value of gap at the associated gauss point, \(\mathbf{x}^{(1)}\) are the coordinates of the overlapping gauss points at slave triangles.

Definition at line 1803 of file SimpleContact.cpp.

1804  {
1806 
1807  if (data.getIndices().size() == 0)
1809 
1810  const int nb_dofs = data.getIndices().size();
1811  if (nb_dofs) {
1812  const int nb_gauss_pts = data.getN().size1();
1813  int nb_base_fun_col = data.getFieldData().size() / 3;
1814 
1815  vecF.resize(nb_dofs, false);
1816  vecF.clear();
1817 
1818  const double area_s =
1819  commonDataSimpleContact->areaSlave; // same area in master and slave
1820 
1821  auto get_tensor_vec = [](VectorDouble &n, const int r) {
1822  return FTensor::Tensor1<double *, 3>(&n(r + 0), &n(r + 1), &n(r + 2));
1823  };
1824 
1826 
1827  auto const_unit_n = get_tensor_vec(
1828  commonDataSimpleContact->normalVectorSlavePtr.get()[0], 0);
1829 
1830  auto t_w = getFTensor0IntegrationWeightSlave();
1831  auto t_aug_lambda_ptr =
1832  getFTensor0FromVec(*commonDataSimpleContact->augmentedLambdasPtr);
1833 
1834  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
1835  if (t_aug_lambda_ptr > 0 && std::abs(t_aug_lambda_ptr) > ALM_TOL) {
1836  ++t_aug_lambda_ptr;
1837  ++t_w;
1838  continue;
1839  }
1840  const double val_s = t_w * area_s;
1841 
1842  FTensor::Tensor0<double *> t_base_slave(&data.getN()(gg, 0));
1843 
1844  for (int bbc = 0; bbc != nb_base_fun_col; ++bbc) {
1845 
1846  auto t_assemble_s = get_tensor_vec(vecF, 3 * bbc);
1847 
1848  t_assemble_s(i) -=
1849  val_s * const_unit_n(i) * t_aug_lambda_ptr * t_base_slave;
1850 
1851  ++t_base_slave;
1852  }
1853  ++t_aug_lambda_ptr;
1854  ++t_w;
1855  } // for gauss points
1856 
1857  CHKERR VecSetValues(getSNESf(), data, &*vecF.begin(), ADD_VALUES);
1858  }
1860 }

Member Data Documentation

◆ commonDataSimpleContact

boost::shared_ptr<CommonDataSimpleContact> SimpleContactProblem::OpCalAugmentedTractionRhsSlave::commonDataSimpleContact
private

Definition at line 1148 of file SimpleContact.hpp.

◆ vecF

VectorDouble SimpleContactProblem::OpCalAugmentedTractionRhsSlave::vecF
private

Definition at line 1149 of file SimpleContact.hpp.


The documentation for this struct was generated from the following files:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
SimpleContactProblem::OpCalAugmentedTractionRhsSlave::vecF
VectorDouble vecF
Definition: SimpleContact.hpp:1149
MoFEM::EntitiesFieldData::EntData::getFieldData
const VectorDouble & getFieldData() const
get dofs values
Definition: EntitiesFieldData.hpp:1254
SimpleContactProblem::OpCalAugmentedTractionRhsSlave::commonDataSimpleContact
boost::shared_ptr< CommonDataSimpleContact > commonDataSimpleContact
Definition: SimpleContact.hpp:1148
MoFEM::VecSetValues
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
Definition: EntitiesFieldData.hpp:1589
sdf.r
int r
Definition: sdf.py:8
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPCOL
@ OPCOL
operator doWork function is executed on FE columns
Definition: ForcesAndSourcesCore.hpp:568
MoFEM::getFTensor0FromVec
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
MoFEM::EntitiesFieldData::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Definition: EntitiesFieldData.hpp:1214
SimpleContactProblem::ALM_TOL
static constexpr double ALM_TOL
Definition: SimpleContact.hpp:65
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
FTensor::Index< 'i', 3 >
convert.n
n
Definition: convert.py:82
FTensor::Tensor0
Definition: Tensor0.hpp:16
MoFEM::EntitiesFieldData::EntData::getN
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Definition: EntitiesFieldData.hpp:1318
UserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Definition: HookeElement.hpp:75
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359