v0.13.2
Loading...
Searching...
No Matches
Public Member Functions | Private Attributes | List of all members
SimpleContactProblem::OpCalDerIntCompFunOverSpatPosSlaveSlave Struct Reference

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

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

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

Public Member Functions

 OpCalDerIntCompFunOverSpatPosSlaveSlave (const string lagrange_field_name, const string field_name, boost::shared_ptr< CommonDataSimpleContact > common_data_contact, boost::shared_ptr< double > cn)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
 Integrates linearisation of the complementarity function at slave face gauss points and assembles components to LHS global matrix. More...
 

Private Attributes

boost::shared_ptr< CommonDataSimpleContactcommonDataSimpleContact
 
boost::shared_ptr< doublecNPtr
 
MatrixDouble NN
 

Detailed Description

LHS-operator for the simple contact element.

Integrates the variation with respect to slave spatial positions of the complementarity function to fulfill KKT conditions in the integral sense and assembles components to LHS global matrix.

Definition at line 1512 of file SimpleContact.hpp.

Constructor & Destructor Documentation

◆ OpCalDerIntCompFunOverSpatPosSlaveSlave()

SimpleContactProblem::OpCalDerIntCompFunOverSpatPosSlaveSlave::OpCalDerIntCompFunOverSpatPosSlaveSlave ( const string  lagrange_field_name,
const string  field_name,
boost::shared_ptr< CommonDataSimpleContact common_data_contact,
boost::shared_ptr< double cn 
)
inline

Definition at line 1514 of file SimpleContact.hpp.

1518 : ContactOp(lagrange_field_name, field_name, UserDataOperator::OPROWCOL,
1519 ContactOp::FACESLAVESLAVE),
1520 cNPtr(cn), commonDataSimpleContact(common_data_contact) {
1521 sYmm = false; // This will make sure to loop over all entities (e.g.
1522 // for order=2 it will make doWork to loop 16 time)
1523 }
double cn
Definition: contact.cpp:124
constexpr auto field_name
boost::shared_ptr< CommonDataSimpleContact > commonDataSimpleContact
ContactPrismElementForcesAndSourcesCore::UserDataOperator ContactOp

Member Function Documentation

◆ doWork()

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

Integrates linearisation of the complementarity function at slave face gauss points and assembles components to LHS global matrix.

Integrates and assembles the variation with respect to slave spatial positions of the complementarity function to fulfill KKT conditions in the integral sense and assembles components to LHS global matrix.

\[ {\text D}{\overline C(\lambda, \mathbf{x}^{(i)}, \delta \lambda)}[\Delta \mathbf{x}^{(1)}] = \int_{{\gamma}^{(1)}_{\text c}} -\Delta \mathbf{x}^{(1)} \cdot \mathbf{n}(\mathbf{x}^{(1)}) c_{\text n} \left( 1 + {\text {sign}}\left( \lambda - c_{\text n} g_{\textrm{n}} \right) {\left| \lambda - c_{\text n} g_{\textrm{n}}\right|}^{r-1}\right) \delta{{\lambda}} \,\,{ {\text d} {\gamma}} \]

where \({\gamma}^{(1)}_{\text c}\) is the surface integration domain of the slave surface, \( \lambda\) is the Lagrange multiplier, \(\mathbf{x}^{(i)}\) are the coordinates of the overlapping gauss points at slave and master triangles for \(i = 1\) and \(i = 2\), respectively. Furthermore, \( c_{\text n}\) works as an augmentation parameter and affects convergence, \(r\) is regularisation parameter that can be chosen in \([1, 1.1]\) ( \(r = 1\) is the default value) and \( g_{\textrm{n}}\) is the gap function evaluated at the slave triangle gauss points as:

\[ g_{\textrm{n}} = - \mathbf{n}(\mathbf{x}^{(1)}) \cdot \left( \mathbf{x}^{(1)} - \mathbf{x}^{(2)} \right) \]

Definition at line 2127 of file SimpleContact.cpp.

2129 {
2131
2132 const int nb_row = row_data.getIndices().size();
2133 const int nb_col = col_data.getIndices().size();
2134
2135 if (nb_row && nb_col) {
2136
2137 const int nb_gauss_pts = row_data.getN().size1();
2138 int nb_base_fun_row = row_data.getFieldData().size();
2139 int nb_base_fun_col = col_data.getFieldData().size() / 3;
2140
2141 const double area_slave =
2142 commonDataSimpleContact->areaSlave; // same area in master and slave
2143
2144 NN.resize(nb_base_fun_row, 3 * nb_base_fun_col, false);
2145 NN.clear();
2146
2147 auto get_tensor_from_vec = [](VectorDouble &n) {
2148 return FTensor::Tensor1<double *, 3>(&n(0), &n(1), &n(2));
2149 };
2150
2152
2153 auto t_lagrange_slave =
2154 getFTensor0FromVec(*commonDataSimpleContact->lagMultAtGaussPtsPtr);
2155 auto t_gap_gp = getFTensor0FromVec(*commonDataSimpleContact->gapPtr);
2156
2157 auto t_const_unit_n =
2158 get_tensor_from_vec(*(commonDataSimpleContact->normalVectorSlavePtr));
2159
2160 auto t_w = getFTensor0IntegrationWeightSlave();
2161 const double cn_value = *cNPtr.get();
2162 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
2163 const double val_m = SimpleContactProblem::ConstrainFunction_dg(
2164 cn_value, t_gap_gp, t_lagrange_slave) *
2165 t_w * area_slave;
2166
2167 auto t_base_lambda = row_data.getFTensor0N(gg, 0);
2168 for (int bbr = 0; bbr != nb_base_fun_row; ++bbr) {
2169
2170 auto t_base_slave = col_data.getFTensor0N(gg, 0);
2172 &NN(bbr, 0), &NN(bbr, 1), &NN(bbr, 2)};
2173 const double m = val_m * t_base_lambda;
2174
2175 for (int bbc = 0; bbc != nb_base_fun_col; ++bbc) {
2176
2177 t_mat(i) -= t_const_unit_n(i) * m * t_base_slave;
2178
2179 ++t_base_slave; // update rows
2180 ++t_mat;
2181 }
2182 ++t_base_lambda; // update cols master
2183 }
2184
2185 ++t_gap_gp;
2186 ++t_lagrange_slave;
2187 ++t_w;
2188 }
2189
2190 CHKERR MatSetValues(getSNESB(), row_data, col_data, &*NN.data().begin(),
2191 ADD_VALUES);
2192 }
2193
2195}
#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
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
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 double ConstrainFunction_dg(const double cn, const double g, const double l)

Member Data Documentation

◆ cNPtr

boost::shared_ptr<double> SimpleContactProblem::OpCalDerIntCompFunOverSpatPosSlaveSlave::cNPtr
private

Definition at line 1563 of file SimpleContact.hpp.

◆ commonDataSimpleContact

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

Definition at line 1562 of file SimpleContact.hpp.

◆ NN

MatrixDouble SimpleContactProblem::OpCalDerIntCompFunOverSpatPosSlaveSlave::NN
private

Definition at line 1564 of file SimpleContact.hpp.


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