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

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

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

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

Public Member Functions

 OpCalContactTractionOverLambdaMasterSlave (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 Lagrange multipliers virtual work, \( \delta W_{\text c}\) derivative with respect to Lagrange multipliers with respect to Lagrange multipliers on master side and assembles components to LHS global matrix. More...
 

Private Attributes

boost::shared_ptr< CommonDataSimpleContactcommonDataSimpleContact
 
MatrixDouble NN
 

Detailed Description

LHS-operator for the simple contact element.

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

Definition at line 1277 of file SimpleContact.hpp.

Constructor & Destructor Documentation

◆ OpCalContactTractionOverLambdaMasterSlave()

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

Definition at line 1279 of file SimpleContact.hpp.

1282 : ContactOp(field_name, lagrange_field_name, UserDataOperator::OPROWCOL,
1283 ContactOp::FACEMASTERSLAVE),
1284 commonDataSimpleContact(common_data_contact) {
1285 sYmm = false; // This will make sure to loop over all intities (e.g.
1286 // for order=2 it will make doWork to loop 16 time)
1287 }
constexpr auto field_name
boost::shared_ptr< CommonDataSimpleContact > commonDataSimpleContact
ContactPrismElementForcesAndSourcesCore::UserDataOperator ContactOp

Member Function Documentation

◆ doWork()

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

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

Computes linearisation of integrated on slave side complementarity function and assembles derivative of Lagrange multipliers virtual work \( \delta W_{\text c}\) with respect to Lagrange multipliers and assembles components to LHS global matrix

\[ {\text D} {\delta W^{(2)}_{\text c}(\lambda, \delta \mathbf{x}^{(2)}})[\Delta \lambda] \,\, = \int_{{\gamma}^{(1)}_{\text c}} \Delta \lambda \delta{\mathbf{x}^{(2)}} \,\,{ {\text d} {\gamma}} \]

where \({\gamma}^{(1)}_{\text c}\) is the surface integration domain of the slave surface, \( \lambda\) is the Lagrange multiplier, \(\mathbf{x}^{(2)}\) are the coordinates of the overlapping gauss points at master triangles.

Definition at line 1738 of file SimpleContact.cpp.

1740 {
1742
1743 // Both sides are needed since both sides contribute their shape
1744 // function to the stiffness matrix
1745 const int nb_row = row_data.getIndices().size();
1746 const int nb_col = col_data.getIndices().size();
1747
1748 if (nb_row && nb_col) {
1749 const int nb_gauss_pts = row_data.getN().size1();
1750
1751 int nb_base_fun_row = row_data.getFieldData().size() / 3;
1752 int nb_base_fun_col = col_data.getFieldData().size();
1753
1754 const double area_master =
1755 commonDataSimpleContact->areaMaster; // same area in master and slave
1756
1757 auto get_tensor_from_mat = [](MatrixDouble &m, const int r, const int c) {
1758 return FTensor::Tensor1<double *, 3>(&m(r + 0, c), &m(r + 1, c),
1759 &m(r + 2, c));
1760 };
1761
1762 auto get_tensor_vec = [](VectorDouble &n) {
1763 return FTensor::Tensor1<double *, 3>(&n(0), &n(1), &n(2));
1764 };
1765
1767
1768 NN.resize(3 * nb_base_fun_row, nb_base_fun_col, false);
1769 NN.clear();
1770
1771 auto const_unit_n =
1772 get_tensor_vec(*(commonDataSimpleContact->normalVectorSlavePtr.get()));
1773
1774 auto t_w = getFTensor0IntegrationWeightMaster();
1775
1776 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
1777
1778 double val_m = t_w * area_master;
1779 auto t_base_master = row_data.getFTensor0N(gg, 0);
1780
1781 for (int bbr = 0; bbr != nb_base_fun_row; ++bbr) {
1782 auto t_assemble_m = get_tensor_from_mat(NN, 3 * bbr, 0);
1783 auto t_base_lambda = col_data.getFTensor0N(gg, 0);
1784 const double m = val_m * t_base_master;
1785 for (int bbc = 0; bbc != nb_base_fun_col; ++bbc) {
1786 const double n = m * t_base_lambda;
1787 t_assemble_m(i) -= n * const_unit_n(i);
1788 ++t_assemble_m;
1789 ++t_base_lambda; // update cols slave
1790 }
1791 ++t_base_master; // update rows master
1792 }
1793 ++t_w;
1794 }
1795
1796 CHKERR MatSetValues(getSNESB(), row_data, col_data, &*NN.data().begin(),
1797 ADD_VALUES);
1798 }
1799
1801}
#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.
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.

Member Data Documentation

◆ commonDataSimpleContact

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

Definition at line 1319 of file SimpleContact.hpp.

◆ NN

MatrixDouble SimpleContactProblem::OpCalContactTractionOverLambdaMasterSlave::NN
private

Definition at line 1320 of file SimpleContact.hpp.


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