v0.14.0
Loading...
Searching...
No Matches
MortarContactInterface.hpp
Go to the documentation of this file.
1/** \file MortarContactInterface.hpp
2 */
3
4/* This file is part of MoFEM.
5 * MoFEM is free software: you can redistribute it and/or modify it under
6 * the terms of the GNU Lesser General Public License as published by the
7 * Free Software Foundation, either version 3 of the License, or (at your
8 * option) any later version.
9 *
10 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13 * License for more details.
14 *
15 * You should have received a copy of the GNU Lesser General Public
16 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
17#pragma once
18
19#ifndef __MORTARCONTACTINTERFACE_HPP__
20#define __MORTARCONTACTINTERFACE_HPP__
21
23
25 SmartPetscObj<DM> dM;
26
27 // params
28 double cnValue;
29 PetscInt oRder;
30 PetscInt orderContact;
31 PetscInt orderLambda;
32 PetscInt nbHoLevels;
33 PetscBool isNewtonCotes;
34 PetscBool convectPts;
36 PetscBool almFlag;
37 PetscBool dEbug;
38 PetscBool isPartitioned;
39
40 // data structs
41 boost::shared_ptr<MortarContactProblem> contactProblemPtr;
44 std::vector<std::pair<Range, Range>> contactSurfacePairs; // <Slave, Master>
45
48
51
54
55 std::vector<BitRefLevel> bitLevels;
56
58 boost::shared_ptr<contactMIndex> contactCommondataMultiIndex;
59 boost::shared_ptr<PostProcFaceOnRefinedMesh> postProcContactPtr;
60
61 boost::shared_ptr<MortarContactProblem::MortarContactElement>
63 boost::shared_ptr<MortarContactProblem::CommonDataMortarContact>
65
66 moab::Core mb_post;
67 moab::Interface &moabPostProcMesh;
68 std::array<double, 2> nbGaussPts;
69 std::array<double, 2> contactArea;
70
71 MortarContactInterface(MoFEM::Interface &m_field, string postion_field,
72 string mesh_posi_field_name = "MESH_NODE_POSITIONS",
73 bool is_displacement_field = false,
74 bool is_quasi_static = true)
75 : mField(m_field), positionField(postion_field),
76 meshNodeField(mesh_posi_field_name),
77 isDisplacementField(is_displacement_field),
79
80 oRder = -1;
81 orderContact = 1;
82 nbHoLevels = 0;
83 orderLambda = 1;
84 cnValue = -1;
85 postProcGapTol = 0.;
86 almFlag = PETSC_FALSE;
87 isNewtonCotes = PETSC_FALSE;
88 convectPts = PETSC_FALSE;
89 printContactState = PETSC_FALSE;
90 dEbug = PETSC_FALSE;
91 isPartitioned = PETSC_FALSE;
92 }
94
95 // interface functions
96 MoFEMErrorCode getCommandLineParameters() override;
97 MoFEMErrorCode addElementFields() override;
98
99 MoFEMErrorCode createElements() override;
100 MoFEMErrorCode setOperators() override;
101 BitRefLevel getBitRefLevel() override;
102 MoFEMErrorCode addElementsToDM(SmartPetscObj<DM> dm) override;
103
104 // MoFEMErrorCode setupSolverKSP() override;
105 MoFEMErrorCode setupSolverJacobianSNES() override;
106 MoFEMErrorCode setupSolverFunctionSNES() override;
107 MoFEMErrorCode setupSolverJacobianTS(const TSType type) override;
108 MoFEMErrorCode setupSolverFunctionTS(const TSType type) override;
109
110 MoFEMErrorCode updateElementVariables() override;
111 MoFEMErrorCode postProcessElement(int step) override;
112
113 // contact interface
114 std::vector<BitRefLevel> &getBitRefLevelVector();
115 MoFEMErrorCode setupElementEntities();
116 MoFEMErrorCode setPositionFieldOrder(int order);
117
118 boost::ptr_deque<MoFEM::ForcesAndSourcesCore::UserDataOperator> &
120
121 MoFEMErrorCode findContactSurfacePairs();
122 MoFEMErrorCode findPostProcSurface();
123
124 template <typename T, bool RHS>
125 MoFEMErrorCode setupSolverImpl(const TSType type = IM);
126 template <typename T>
127 MoFEMErrorCode setupSolverFunction(const TSType type = IM);
128 template <typename T>
129 MoFEMErrorCode setupSolverJacobian(const TSType type = IM);
130
131 MoFEMErrorCode postProcessContactState(int step);
132};
133
134template <>
135MoFEMErrorCode
136MortarContactInterface::setupSolverJacobian<SNES>(const TSType type);
137template <>
138MoFEMErrorCode
139MortarContactInterface::setupSolverJacobian<SNES>(const TSType type);
140template <>
141MoFEMErrorCode
142MortarContactInterface::setupSolverFunction<TS>(const TSType type);
143template <>
144MoFEMErrorCode
145MortarContactInterface::setupSolverFunction<TS>(const TSType type);
146
147#endif // __MORTARCONTACTINTERFACE_HPP__
constexpr int order
PetscBool is_quasi_static
Definition: plastic.cpp:190
multi_index_container< boost::shared_ptr< ContactCommonData >, indexed_by< ordered_unique< tag< Prism_tag >, member< ContactCommonData, EntityHandle, &ContactCommonData::pRism > > > > ContactCommonData_multiIndex
Set of functions declaring elements and setting operators for generic element interface.
Deprecated interface functions.
MoFEMErrorCode setupSolverFunctionTS(const TSType type) override
std::vector< BitRefLevel > & getBitRefLevelVector()
MoFEMErrorCode getCommandLineParameters() override
MoFEMErrorCode findContactSurfacePairs()
MoFEMErrorCode postProcessContactState(int step)
boost::shared_ptr< PostProcFaceOnRefinedMesh > postProcContactPtr
MoFEMErrorCode setupSolverFunction(const TSType type=IM)
MoFEMErrorCode setupElementEntities()
boost::shared_ptr< MortarContactProblem::MortarContactElement > fePostProcSimpleContact
MoFEMErrorCode createElements() override
MoFEMErrorCode setupSolverJacobianSNES() override
MoFEMErrorCode setupSolverJacobianTS(const TSType type) override
MoFEMErrorCode setPositionFieldOrder(int order)
std::vector< BitRefLevel > bitLevels
std::array< double, 2 > nbGaussPts
MoFEMErrorCode setupSolverImpl(const TSType type=IM)
MoFEMErrorCode postProcessElement(int step) override
MoFEMErrorCode setOperators() override
boost::shared_ptr< contactMIndex > contactCommondataMultiIndex
MoFEMErrorCode addElementsToDM(SmartPetscObj< DM > dm) override
MoFEMErrorCode setupSolverJacobian(const TSType type=IM)
MoFEMErrorCode setupSolverFunctionSNES() override
BitRefLevel getBitRefLevel() override
boost::ptr_deque< MoFEM::ForcesAndSourcesCore::UserDataOperator > & getPostProcSimpleContactPipeline()
boost::shared_ptr< MortarContactProblem::CommonDataMortarContact > commonDataPostProc
MortarContactInterface(MoFEM::Interface &m_field, string postion_field, string mesh_posi_field_name="MESH_NODE_POSITIONS", bool is_displacement_field=false, bool is_quasi_static=true)
moab::Interface & moabPostProcMesh
std::array< double, 2 > contactArea
ContactSearchKdTree::ContactCommonData_multiIndex contactMIndex
MoFEMErrorCode updateElementVariables() override
MoFEMErrorCode addElementFields() override
MortarContactInterface()=delete
boost::shared_ptr< MortarContactProblem > contactProblemPtr
std::vector< std::pair< Range, Range > > contactSurfacePairs