v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Private Attributes | List of all members
SimpleContactProblem::OpCalContactAugmentedTractionOverLambdaSlaveSlave 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::OpCalContactAugmentedTractionOverLambdaSlaveSlave:
[legend]
Collaboration diagram for SimpleContactProblem::OpCalContactAugmentedTractionOverLambdaSlaveSlave:
[legend]

Public Member Functions

 OpCalContactAugmentedTractionOverLambdaSlaveSlave (const string field_name, const string lagrange_field_name, boost::shared_ptr< CommonDataSimpleContact > common_data_contact)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
 Integrates virtual work on slave side , \( \delta W_{\text c}\), derivative with respect to Lagrange multipliers on slave side and assembles its components to LHS global matrix. More...
 

Private Attributes

boost::shared_ptr< CommonDataSimpleContactcommonDataSimpleContact
 
MatrixDouble NN
 

Detailed Description

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

Integrates Lagrange multipliers virtual work, \( \delta W_{\text c}\) derivative with respect to Lagrange multipliers on slave side and assembles components of the LHS matrix.

Definition at line 1642 of file SimpleContact.hpp.

Constructor & Destructor Documentation

◆ OpCalContactAugmentedTractionOverLambdaSlaveSlave()

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

Definition at line 1645 of file SimpleContact.hpp.

1648 : ContactOp(field_name, lagrange_field_name, UserDataOperator::OPROWCOL,
1649 ContactOp::FACESLAVESLAVE),
1650 commonDataSimpleContact(common_data_contact) {
1651 sYmm = false; // This will make sure to loop over all intities (e.g.
1652 // for order=2 it will make doWork to loop 16 time)
1653 }
constexpr auto field_name
boost::shared_ptr< CommonDataSimpleContact > commonDataSimpleContact
ContactPrismElementForcesAndSourcesCore::UserDataOperator ContactOp

Member Function Documentation

◆ doWork()

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

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

Computes linearisation of virtual work on slave side integrated on the slave side and assembles the components of its derivative over Lagrange multipliers.

\[ {\text D} {\delta W^{(1)}_{\text c}(\lambda, \delta \mathbf{x}^{(1)}})[\Delta \lambda] \,\, = \left\{ \begin{array}{ll} \int_{{\gamma}^{(1)}_{\text c}} - \Delta \lambda {\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}^{(1)}\) are the coordinates of the overlapping gauss points at master triangles, \( 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 878 of file SimpleContact.cpp.

880 {
882
883 // Both sides are needed since both sides contribute their shape
884 // function to the stiffness matrix
885 const int nb_row = row_data.getIndices().size();
886 const int nb_col = col_data.getIndices().size();
887
888 if (nb_row && nb_col) {
889 const int nb_gauss_pts = row_data.getN().size1();
890
891 int nb_base_fun_row = row_data.getFieldData().size() / 3;
892 int nb_base_fun_col = col_data.getFieldData().size();
893
894 const double area_slave =
895 commonDataSimpleContact->areaSlave; // same area in master and slave
896
897 auto get_tensor_from_mat = [](MatrixDouble &m, const int r, const int c) {
898 return FTensor::Tensor1<double *, 3>(&m(r + 0, c), &m(r + 1, c),
899 &m(r + 2, c));
900 };
901
902 auto get_tensor_vec = [](VectorDouble &n) {
903 return FTensor::Tensor1<double *, 3>(&n(0), &n(1), &n(2));
904 };
905
907
908 NN.resize(3 * nb_base_fun_row, nb_base_fun_col, false);
909 NN.clear();
910
911 auto const_unit_n =
912 get_tensor_vec(*(commonDataSimpleContact->normalVectorSlavePtr.get()));
913
914 auto t_w = getFTensor0IntegrationWeightSlave();
915
916 auto t_aug_lambda_ptr =
917 getFTensor0FromVec(*commonDataSimpleContact->augmentedLambdasPtr);
918
919 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
920
921 if (t_aug_lambda_ptr > 0 && std::abs(t_aug_lambda_ptr) > ALM_TOL) {
922 ++t_w;
923 ++t_aug_lambda_ptr;
924 continue;
925 }
926
927 double val_m = t_w * area_slave;
928 auto t_base_master = row_data.getFTensor0N(gg, 0);
929
930 for (int bbr = 0; bbr != nb_base_fun_row; ++bbr) {
931 auto t_assemble_m = get_tensor_from_mat(NN, 3 * bbr, 0);
932 auto t_base_lambda = col_data.getFTensor0N(gg, 0);
933 const double m = val_m * t_base_master;
934 for (int bbc = 0; bbc != nb_base_fun_col; ++bbc) {
935 const double n = m * t_base_lambda;
936 t_assemble_m(i) -= n * const_unit_n(i);
937
938 ++t_assemble_m;
939 ++t_base_lambda; // update cols slave
940 }
941 ++t_base_master; // update rows master
942 }
943 ++t_w;
944 ++t_aug_lambda_ptr;
945 }
946
947 CHKERR MatSetValues(getSNESB(), row_data, col_data, &*NN.data().begin(),
948 ADD_VALUES);
949 }
951}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
int r
Definition: sdf.py:8
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorDouble & getFieldData() const
get dofs values
const VectorInt & getIndices() const
Get global indices of dofs on entity.
static constexpr double ALM_TOL

Member Data Documentation

◆ commonDataSimpleContact

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

Definition at line 1693 of file SimpleContact.hpp.

◆ NN

MatrixDouble SimpleContactProblem::OpCalContactAugmentedTractionOverLambdaSlaveSlave::NN
private

Definition at line 1694 of file SimpleContact.hpp.


The documentation for this struct was generated from the following files: