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

Calculate tangent operator for contact force for change of integration point positions, as result of change os spatial positions. More...

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

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

Public Member Functions

 OpLhsConvectIntegrationPtsContactTraction (const string row_field_name, const string col_field_name, boost::shared_ptr< CommonDataSimpleContact > common_data_contact, const ContactOp::FaceType face_type, boost::shared_ptr< MatrixDouble > diff_convect)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
 

Public Attributes

boost::shared_ptr< CommonDataSimpleContactcommonDataSimpleContact
 

Private Attributes

MatrixDouble matLhs
 
boost::shared_ptr< MatrixDouble > diffConvect
 

Detailed Description

Calculate tangent operator for contact force for change of integration point positions, as result of change os spatial positions.

Definition at line 2361 of file SimpleContact.hpp.

Constructor & Destructor Documentation

◆ OpLhsConvectIntegrationPtsContactTraction()

SimpleContactProblem::OpLhsConvectIntegrationPtsContactTraction::OpLhsConvectIntegrationPtsContactTraction ( const string  row_field_name,
const string  col_field_name,
boost::shared_ptr< CommonDataSimpleContact common_data_contact,
const ContactOp::FaceType  face_type,
boost::shared_ptr< MatrixDouble >  diff_convect 
)
inline

Definition at line 2365 of file SimpleContact.hpp.

2370 : ContactOp(row_field_name, col_field_name, UserDataOperator::OPROWCOL,
2371 face_type),
2372 commonDataSimpleContact(common_data_contact),
2373 diffConvect(diff_convect) {
2374 sYmm = false;
2375 }
boost::shared_ptr< CommonDataSimpleContact > commonDataSimpleContact
ContactPrismElementForcesAndSourcesCore::UserDataOperator ContactOp

Member Function Documentation

◆ doWork()

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

Definition at line 2873 of file SimpleContact.cpp.

2875 {
2877
2878 const int nb_row_dofs = row_data.getIndices().size();
2879 const int nb_col_dofs = col_data.getIndices().size();
2880 if (nb_row_dofs && (nb_col_dofs && col_type == MBVERTEX)) {
2881
2882 const int nb_gauss_pts = getGaussPtsSlave().size2();
2883 int nb_base_fun_row = nb_row_dofs / 3;
2884 int nb_base_fun_col = nb_col_dofs / 3;
2885
2886 matLhs.resize(nb_row_dofs, nb_col_dofs, false);
2887 matLhs.clear();
2888
2889 const double area_s =
2890 commonDataSimpleContact->areaSlave; // same area in master and slave
2891
2895
2896 auto get_tensor_vec = [](VectorDouble &n) {
2897 return FTensor::Tensor1<double *, 3>(&n[0], &n[1], &n[2]);
2898 };
2899
2900 auto t_const_unit_n =
2901 get_tensor_vec(*(commonDataSimpleContact->normalVectorSlavePtr));
2902 auto t_diff_lambda_xi = getFTensor1FromMat<2>(
2903 *(commonDataSimpleContact->gradKsiLambdaAtGaussPtsPtr));
2904 auto t_w = getFTensor0IntegrationWeightSlave();
2905
2906 auto get_diff_ksi = [](auto &m, auto gg) {
2908 &m(0, gg), &m(1, gg), &m(2, gg), &m(3, gg), &m(4, gg), &m(5, gg));
2909 };
2910
2911 auto t_base_row = row_data.getFTensor0N();
2912
2913 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
2914
2915 double val_s = t_w * area_s;
2916 auto t_base_row = row_data.getFTensor0N(gg, 0);
2917
2918 for (int rr = 0; rr != nb_base_fun_row; ++rr) {
2919
2921 &matLhs(3 * rr + 0, 0), &matLhs(3 * rr + 0, 1),
2922 &matLhs(3 * rr + 0, 2),
2923
2924 &matLhs(3 * rr + 1, 0), &matLhs(3 * rr + 1, 1),
2925 &matLhs(3 * rr + 1, 2),
2926
2927 &matLhs(3 * rr + 2, 0), &matLhs(3 * rr + 2, 1),
2928 &matLhs(3 * rr + 2, 2)};
2929
2930 auto t_diff_convect = get_diff_ksi(*diffConvect, 3 * gg);
2931
2932 for (int cc = 0; cc != nb_base_fun_col; ++cc) {
2933 t_mat(i, j) -= val_s * t_base_row * t_const_unit_n(i) *
2934 (t_diff_lambda_xi(I) * t_diff_convect(I, j));
2935
2936 ++t_diff_convect;
2937 ++t_mat;
2938 }
2939
2940 ++t_base_row;
2941 }
2942
2943 ++t_diff_lambda_xi;
2944 ++t_w;
2945 } // for gauss points
2946
2947 CHKERR MatSetValues(getSNESB(), row_data, col_data, &*matLhs.data().begin(),
2948 ADD_VALUES);
2949 }
2950
2952}
#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
FTensor::Index< 'j', 3 > j
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
constexpr IntegrationType I
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
const VectorInt & getIndices() const
Get global indices of dofs on entity.

Member Data Documentation

◆ commonDataSimpleContact

boost::shared_ptr<CommonDataSimpleContact> SimpleContactProblem::OpLhsConvectIntegrationPtsContactTraction::commonDataSimpleContact

Definition at line 2363 of file SimpleContact.hpp.

◆ diffConvect

boost::shared_ptr<MatrixDouble> SimpleContactProblem::OpLhsConvectIntegrationPtsContactTraction::diffConvect
private

Definition at line 2383 of file SimpleContact.hpp.

◆ matLhs

MatrixDouble SimpleContactProblem::OpLhsConvectIntegrationPtsContactTraction::matLhs
private

Definition at line 2382 of file SimpleContact.hpp.


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