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

LHS-operator for the simple contact element with Augmented Lagrangian Method. More...

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

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

Public Member Functions

 OpCalContactAugmentedTractionOverSpatialSlaveMaster (const string field_name, const string field_name_2, const double cn_value, boost::shared_ptr< CommonDataSimpleContact > common_data_contact)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
 Integrates virtual work on slave side, \( \delta W_{\text c}\), derivative with respect to spatial positions on master side and assembles its components to the global LHS matrix. More...
 

Private Attributes

boost::shared_ptr< CommonDataSimpleContactcommonDataSimpleContact
 
const double cN
 
MatrixDouble NN
 

Detailed Description

LHS-operator for the simple contact element with Augmented Lagrangian Method.

Integrates virtual work of spatial position on slave side, \( \delta W_{\text c}\), derivative with respect to spatial positions on master side and assembles its components to the global LHS matrix.

Definition at line 1917 of file SimpleContact.hpp.

Constructor & Destructor Documentation

◆ OpCalContactAugmentedTractionOverSpatialSlaveMaster()

SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::OpCalContactAugmentedTractionOverSpatialSlaveMaster ( const string  field_name,
const string  field_name_2,
const double  cn_value,
boost::shared_ptr< CommonDataSimpleContact common_data_contact 
)
inline

Definition at line 1920 of file SimpleContact.hpp.

1925  field_name, field_name_2, UserDataOperator::OPROWCOL,
1927  FACESLAVEMASTER),
1928  cN(cn_value), commonDataSimpleContact(common_data_contact) {
1929  sYmm = false; // This will make sure to loop over all intities (e.g.
1930  // for order=2 it will make doWork to loop 16 time)
1931  }

Member Function Documentation

◆ doWork()

MoFEMErrorCode SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
EntitiesFieldData::EntData row_data,
EntitiesFieldData::EntData col_data 
)

Integrates virtual work on slave side, \( \delta W_{\text c}\), derivative with respect to spatial positions on master side and assembles its components to the global LHS matrix.

Computes linearisation of virtual work of spatial position on slave side , \( \delta W_{\text c}\), over master side spatial positions and assembles its components to LHS global matrix

\[ {\text D} {\delta W^{(1)}_{\text c}(\lambda, \Delta \mathbf{x}^{(2)}})[\delta {\mathbf{x}^{(1)}}] \,\, = \left\{ \begin{array}{ll} \int_{{\gamma}^{(1)}_{\text c}} -c_{\textrm n}\Delta {\mathbf{x}^{(2)}}\cdot [{\mathbf{n}}_{\rm c} \otimes {\mathbf{n}}_{\rm c}] \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, \(\mathbf{x}^{(2)}\) and \(\mathbf{x}^{(1)}\) are the coordinates of the overlapping gauss points at master and slave triangles, respectively. Also, \( c_{\textrm n}\) is the regularisation/augmentation parameter of stress dimensions and \( g_{\textrm{n}}\) is the gap evaluated on the corresponding slave side.

Definition at line 1421 of file SimpleContact.cpp.

1423  {
1425 
1426  // Both sides are needed since both sides contribute their shape
1427  // function to the stiffness matrix
1428  const int nb_row = row_data.getIndices().size();
1429  if (!nb_row)
1431  const int nb_col = col_data.getIndices().size();
1432  if (!nb_col)
1434  const int nb_gauss_pts = row_data.getN().size1();
1435 
1436  int nb_base_fun_row = row_data.getFieldData().size() / 3;
1437  int nb_base_fun_col = col_data.getFieldData().size() / 3;
1438 
1439  auto get_tensor_from_mat = [](MatrixDouble &m, const int r, const int c) {
1441  &m(r + 0, c + 0), &m(r + 0, c + 1), &m(r + 0, c + 2), &m(r + 1, c + 0),
1442  &m(r + 1, c + 1), &m(r + 1, c + 2), &m(r + 2, c + 0), &m(r + 2, c + 1),
1443  &m(r + 2, c + 2));
1444  };
1445 
1446  auto get_tensor_vec = [](VectorDouble &n) {
1447  return FTensor::Tensor1<double *, 3>(&n(0), &n(1), &n(2));
1448  };
1449 
1452 
1453  NN.resize(3 * nb_base_fun_row, 3 * nb_base_fun_col, false);
1454  NN.clear();
1455 
1456  auto t_w = getFTensor0IntegrationWeightSlave();
1457 
1458  const double area_slave = commonDataSimpleContact->areaSlave;
1459 
1460  auto normal =
1461  get_tensor_vec(commonDataSimpleContact->normalVectorSlavePtr.get()[0]);
1462 
1463  auto t_aug_lambda_ptr =
1464  getFTensor0FromVec(*commonDataSimpleContact->augmentedLambdasPtr);
1465 
1466  for (int gg = 0; gg != nb_gauss_pts; gg++) {
1467 
1468  if (t_aug_lambda_ptr > 0 && std::abs(t_aug_lambda_ptr) > ALM_TOL) {
1469  ++t_w;
1470  ++t_aug_lambda_ptr;
1471  continue;
1472  }
1473 
1474  const double val_s = t_w * area_slave * cN;
1475 
1476  FTensor::Tensor0<double *> t_base_master_col(&col_data.getN()(gg, 0));
1477 
1478  for (int bbc = 0; bbc != nb_base_fun_col; bbc++) {
1479 
1480  FTensor::Tensor0<double *> t_base_slave_row(&row_data.getN()(gg, 0));
1481  const double s = val_s * t_base_master_col;
1482 
1483  for (int bbr = 0; bbr != nb_base_fun_row; bbr++) {
1484 
1485  auto t_assemble_s = get_tensor_from_mat(NN, 3 * bbr, 3 * bbc);
1486 
1487  t_assemble_s(i, j) -= s * normal(i) * normal(j) * t_base_slave_row;
1488 
1489  ++t_base_slave_row; // update rows
1490  }
1491  ++t_base_master_col; // update cols slave
1492  }
1493  ++t_w;
1494  ++t_aug_lambda_ptr;
1495  }
1496 
1497  CHKERR MatSetValues(getSNESB(), row_data, col_data, &*NN.data().begin(),
1498  ADD_VALUES);
1499 
1501 }

Member Data Documentation

◆ cN

const double SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::cN
private

Definition at line 1974 of file SimpleContact.hpp.

◆ commonDataSimpleContact

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

Definition at line 1973 of file SimpleContact.hpp.

◆ NN

MatrixDouble SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::NN
private

Definition at line 1975 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
MoFEM::MatSetValues
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
Definition: EntitiesFieldData.hpp:1644
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::commonDataSimpleContact
boost::shared_ptr< CommonDataSimpleContact > commonDataSimpleContact
Definition: SimpleContact.hpp:1973
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
MoFEM::EntitiesFieldData::EntData::getFieldData
const VectorDouble & getFieldData() const
get dofs values
Definition: EntitiesFieldData.hpp:1254
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROWCOL
@ OPROWCOL
operator doWork is executed on FE rows &columns
Definition: ForcesAndSourcesCore.hpp:569
sdf.r
int r
Definition: sdf.py:8
FTensor::Tensor2< double *, 3, 3 >
c
const double c
speed of light (cm/ns)
Definition: initial_diffusion.cpp:39
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
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
SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::cN
const double cN
Definition: SimpleContact.hpp:1974
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
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
SimpleContactProblem::OpCalContactAugmentedTractionOverSpatialSlaveMaster::NN
MatrixDouble NN
Definition: SimpleContact.hpp:1975
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
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