v0.9.2
lesson8_contact.cpp
Go to the documentation of this file.
1 /**
2  * \file lesson8_contact.cpp
3  * \example lesson8_contact.cpp
4  *
5  * Example of contact problem
6  *
7  */
8 
9 /* This file is part of MoFEM.
10  * MoFEM is free software: you can redistribute it and/or modify it under
11  * the terms of the GNU Lesser General Public License as published by the
12  * Free Software Foundation, either version 3 of the License, or (at your
13  * option) any later version.
14  *
15  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
16  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
18  * License for more details.
19  *
20  * You should have received a copy of the GNU Lesser General Public
21  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
22 
23 #include <MoFEM.hpp>
24 
25 using namespace MoFEM;
26 
32 
33 constexpr int order = 4;
34 constexpr double young_modulus = 1;
35 constexpr double poisson_ratio = 0.25;
36 constexpr double cn = 1;
37 constexpr double spring_stiffness = 1e-6;
38 
39 #include <ElasticOps.hpp>
40 #include <ContactOps.hpp>
41 using namespace OpElasticTools;
42 using namespace OpContactTools;
43 
44 struct Example {
45 
46  Example(MoFEM::Interface &m_field) : mField(m_field) {}
47 
48  MoFEMErrorCode runProblem();
49 
50 private:
51  MoFEM::Interface &mField;
52 
53  MoFEMErrorCode setUP();
54  MoFEMErrorCode createCommonData();
55  MoFEMErrorCode bC();
56  MoFEMErrorCode OPs();
57  MoFEMErrorCode tsSolve();
58  MoFEMErrorCode postProcess();
59  MoFEMErrorCode checkResults();
60 
62  boost::shared_ptr<OpContactTools::CommonData> commonDataPtr;
63  boost::shared_ptr<PostProcFaceOnRefinedMeshFor2D> postProcFe;
64  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uXScatter;
65  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uYScatter;
66 
67  Range getEntsOnMeshSkin();
68 };
69 
70 //! [Run problem]
73  CHKERR setUP();
74  CHKERR createCommonData();
75  CHKERR bC();
76  CHKERR OPs();
77  CHKERR tsSolve();
78  CHKERR postProcess();
79  CHKERR checkResults();
81 }
82 //! [Run problem]
83 
84 //! [Set up problem]
87  Simple *simple = mField.getInterface<Simple>();
88  // Add field
89  CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 2);
90  CHKERR simple->addBoundaryField("U", H1, AINSWORTH_LEGENDRE_BASE, 2);
91 
92  CHKERR simple->addDomainField("SIGMA", HCURL, DEMKOWICZ_JACOBI_BASE, 2);
93  CHKERR simple->addBoundaryField("SIGMA", HCURL, DEMKOWICZ_JACOBI_BASE, 2);
94 
95  CHKERR simple->setFieldOrder("U", order);
96  CHKERR simple->setFieldOrder("SIGMA", 0);
97 
98  auto skin_edges = getEntsOnMeshSkin();
99 
100  // filter not owned entities, those are not on boundary
101  Range boundary_ents;
102  ParallelComm *pcomm =
103  ParallelComm::get_pcomm(&mField.get_moab(), MYPCOMM_INDEX);
104  if (pcomm == NULL) {
105  pcomm = new ParallelComm(&mField.get_moab(), mField.get_comm());
106  }
107 
108  CHKERR pcomm->filter_pstatus(skin_edges, PSTATUS_SHARED | PSTATUS_MULTISHARED,
109  PSTATUS_NOT, -1, &boundary_ents);
110 
111  CHKERR simple->setFieldOrder("SIGMA", order - 1, &boundary_ents);
112  // CHKERR simple->setFieldOrder("U", order + 1, &skin_edges);
113 
114  CHKERR simple->setUp();
116 }
117 //! [Set up problem]
118 
119 //! [Create common data]
122 
123  auto get_matrial_stiffens = [&](FTensor::Ddg<double, 2, 2> &t_D) {
129  t_D(i, j, k, l) = 0;
130 
131  constexpr double c = young_modulus / (1 - poisson_ratio * poisson_ratio);
132  constexpr double o = poisson_ratio * c;
133 
134  t_D(0, 0, 0, 0) = c;
135  t_D(0, 0, 1, 1) = o;
136 
137  t_D(1, 1, 0, 0) = o;
138  t_D(1, 1, 1, 1) = c;
139 
140  t_D(0, 1, 0, 1) = (1 - poisson_ratio) * c;
142  };
143 
144  commonDataPtr = boost::make_shared<OpContactTools::CommonData>();
145  CHKERR get_matrial_stiffens(commonDataPtr->tD);
146 
147  commonDataPtr->mGradPtr = boost::make_shared<MatrixDouble>();
148  commonDataPtr->mStrainPtr = boost::make_shared<MatrixDouble>();
149  commonDataPtr->mStressPtr = boost::make_shared<MatrixDouble>();
150  commonDataPtr->contactStressPtr = boost::make_shared<MatrixDouble>();
151  commonDataPtr->contactStressDivergencePtr =
152  boost::make_shared<MatrixDouble>();
153  commonDataPtr->contactTractionPtr = boost::make_shared<MatrixDouble>();
154  commonDataPtr->contactDispPtr = boost::make_shared<MatrixDouble>();
155 
156  jAc.resize(2, 2, false);
157  invJac.resize(2, 2, false);
158 
160 }
161 //! [Create common data]
162 
163 //! [Boundary condition]
166 
167  auto fix_disp = [&](const std::string blockset_name) {
168  Range fix_ents;
170  if (it->getName().compare(0, blockset_name.length(), blockset_name) ==
171  0) {
172  CHKERR mField.get_moab().get_entities_by_handle(it->meshset, fix_ents,
173  true);
174  }
175  }
176  return fix_ents;
177  };
178 
179  auto remove_ents = [&](const Range &&ents, const bool fix_x,
180  const bool fix_y) {
181  auto prb_mng = mField.getInterface<ProblemsManager>();
182  auto simple = mField.getInterface<Simple>();
184  Range verts;
185  CHKERR mField.get_moab().get_connectivity(ents, verts, true);
186  verts.merge(ents);
187  const int lo_coeff = fix_x ? 0 : 1;
188  const int hi_coeff = fix_y ? 1 : 0;
189  CHKERR prb_mng->removeDofsOnEntities(simple->getProblemName(), "U", verts,
190  lo_coeff, hi_coeff);
192  };
193 
194  Range boundary_ents;
195  boundary_ents.merge(fix_disp("FIX_X"));
196  boundary_ents.merge(fix_disp("FIX_Y"));
197  boundary_ents.merge(fix_disp("FIX_ALL"));
198  CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
199  mField.getInterface<Simple>()->getProblemName(), "SIGMA", boundary_ents,
200  0, 2);
201 
202  CHKERR remove_ents(fix_disp("FIX_X"), true, false);
203  CHKERR remove_ents(fix_disp("FIX_Y"), false, true);
204  CHKERR remove_ents(fix_disp("FIX_ALL"), true, true);
205 
207 }
208 //! [Boundary condition]
209 
210 //! [Push operators to pipeline]
213  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
214 
215  auto add_domain_base_ops = [&](auto &pipeline) {
216  pipeline.push_back(new OpCalculateJacForFace(jAc));
217  pipeline.push_back(new OpCalculateInvJacForFace(invJac));
218  pipeline.push_back(new OpSetInvJacH1ForFace(invJac));
219  pipeline.push_back(new OpMakeHdivFromHcurl());
220  pipeline.push_back(new OpSetContravariantPiolaTransformFace(jAc));
221  pipeline.push_back(new OpSetInvJacHcurlFace(invJac));
222  };
223 
224  auto add_domain_ops_lhs = [&](auto &pipeline) {
225  pipeline.push_back(new OpStiffnessMatrixLhs("U", "U", commonDataPtr));
226  pipeline.push_back(
227  new OpConstrainDomainLhs_dU("SIGMA", "U", commonDataPtr));
228  };
229 
230  auto add_domain_ops_rhs = [&](auto &pipeline) {
231  auto gravity = [](double x, double y) {
232  return FTensor::Tensor1<double, 2>{0., 1.};
233  };
234  pipeline.push_back(new OpForceRhs("U", commonDataPtr, gravity));
235 
236  pipeline.push_back(
237  new OpCalculateVectorFieldGradient<2, 2>("U", commonDataPtr->mGradPtr));
238  pipeline.push_back(new OpStrain("U", commonDataPtr));
239  pipeline.push_back(new OpStress("U", commonDataPtr));
240  pipeline.push_back(new OpInternalForceRhs("U", commonDataPtr));
241 
242  pipeline.push_back(new OpCalculateVectorFieldValues<2>(
243  "U", commonDataPtr->contactDispPtr));
244 
245  pipeline.push_back(new OpCalculateHVecTensorField<2, 2>(
246  "SIGMA", commonDataPtr->contactStressPtr));
247  pipeline.push_back(new OpCalculateHVecTensorDivergence<2, 2>(
248  "SIGMA", commonDataPtr->contactStressDivergencePtr));
249  pipeline.push_back(new OpConstrainDomainRhs("SIGMA", commonDataPtr));
250  pipeline.push_back(new OpInternalDomainContactRhs("U", commonDataPtr));
251  };
252 
253  auto add_boundary_base_ops = [&](auto &pipeline) {
254  pipeline.push_back(new OpSetContrariantPiolaTransformOnEdge());
255  pipeline.push_back(new OpCalculateVectorFieldValues<2>(
256  "U", commonDataPtr->contactDispPtr));
257  pipeline.push_back(new OpConstrainBoundaryTraction("SIGMA", commonDataPtr));
258  };
259 
260  auto add_boundary_ops_lhs = [&](auto &pipeline) {
261  pipeline.push_back(
262  new OpConstrainBoundaryLhs_dU("SIGMA", "U", commonDataPtr));
263  pipeline.push_back(
264  new OpConstrainBoundaryLhs_dTraction("SIGMA", "SIGMA", commonDataPtr));
265  pipeline.push_back(new OpSpringLhs("U", "U"));
266  };
267 
268  auto add_boundary_ops_rhs = [&](auto &pipeline) {
269  pipeline.push_back(new OpConstrainBoundaryRhs("SIGMA", commonDataPtr));
270  pipeline.push_back(new OpSpringRhs("U", commonDataPtr));
271  };
272 
273  add_domain_base_ops(pipeline_mng->getOpDomainLhsPipeline());
274  add_domain_base_ops(pipeline_mng->getOpDomainRhsPipeline());
275  add_domain_ops_lhs(pipeline_mng->getOpDomainLhsPipeline());
276  add_domain_ops_rhs(pipeline_mng->getOpDomainRhsPipeline());
277 
278  add_boundary_base_ops(pipeline_mng->getOpBoundaryLhsPipeline());
279  add_boundary_base_ops(pipeline_mng->getOpBoundaryRhsPipeline());
280  add_boundary_ops_lhs(pipeline_mng->getOpBoundaryLhsPipeline());
281  add_boundary_ops_rhs(pipeline_mng->getOpBoundaryRhsPipeline());
282 
283  auto integration_rule = [](int, int, int approx_order) { return 2 * order; };
284  CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
285  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
286  CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule);
287  CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(integration_rule);
288 
290 }
291 //! [Push operators to pipeline]
292 
293 //! [Solve]
296 
297  Simple *simple = mField.getInterface<Simple>();
298  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
299  ISManager *is_manager = mField.getInterface<ISManager>();
300 
301  auto solver = pipeline_mng->createTS();
302 
303  auto dm = simple->getDM();
304  auto D = smartCreateDMVector(dm);
305 
306  CHKERR TSSetSolution(solver, D);
307  CHKERR TSSetFromOptions(solver);
308  CHKERR TSSetUp(solver);
309 
310  auto set_section_monitor = [&]() {
312  SNES snes;
313  CHKERR TSGetSNES(solver, &snes);
314  PetscViewerAndFormat *vf;
315  CHKERR PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,
316  PETSC_VIEWER_DEFAULT, &vf);
317  CHKERR SNESMonitorSet(
318  snes,
319  (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal, void *))SNESMonitorFields,
320  vf, (MoFEMErrorCode(*)(void **))PetscViewerAndFormatDestroy);
322  };
323 
324  auto create_post_process_element = [&]() {
326  postProcFe = boost::make_shared<PostProcFaceOnRefinedMeshFor2D>(mField);
327  postProcFe->generateReferenceElementMesh();
328 
329  postProcFe->getOpPtrVector().push_back(new OpCalculateJacForFace(jAc));
330  postProcFe->getOpPtrVector().push_back(
332  postProcFe->getOpPtrVector().push_back(new OpSetInvJacH1ForFace(invJac));
333  postProcFe->getOpPtrVector().push_back(new OpMakeHdivFromHcurl());
334  postProcFe->getOpPtrVector().push_back(
336  postProcFe->getOpPtrVector().push_back(new OpSetInvJacHcurlFace(invJac));
337 
338  postProcFe->getOpPtrVector().push_back(
339  new OpCalculateVectorFieldGradient<2, 2>("U", commonDataPtr->mGradPtr));
340  postProcFe->getOpPtrVector().push_back(new OpStrain("U", commonDataPtr));
341  postProcFe->getOpPtrVector().push_back(new OpStress("U", commonDataPtr));
342  postProcFe->getOpPtrVector().push_back(
344  "SIGMA", commonDataPtr->contactStressDivergencePtr));
345  postProcFe->getOpPtrVector().push_back(new OpCalculateHVecTensorField<2, 2>(
346  "SIGMA", commonDataPtr->contactStressPtr));
347 
348  postProcFe->getOpPtrVector().push_back(new OpPostProcElastic(
349  "U", postProcFe->postProcMesh, postProcFe->mapGaussPts, commonDataPtr));
350  postProcFe->getOpPtrVector().push_back(
351  new OpPostProcContact("SIGMA", postProcFe->postProcMesh,
352  postProcFe->mapGaussPts, commonDataPtr));
353  postProcFe->addFieldValuesPostProc("U");
355  };
356 
357  auto scatter_create = [&](auto coeff) {
359  CHKERR is_manager->isCreateProblemFieldAndRank(simple->getProblemName(),
360  ROW, "U", coeff, coeff, is);
361  int loc_size;
362  CHKERR ISGetLocalSize(is, &loc_size);
363  Vec v;
364  CHKERR VecCreateMPI(mField.get_comm(), loc_size, PETSC_DETERMINE, &v);
365  VecScatter scatter;
366  CHKERR VecScatterCreate(D, is, v, PETSC_NULL, &scatter);
367  return std::make_tuple(SmartPetscObj<Vec>(v),
368  SmartPetscObj<VecScatter>(scatter));
369  };
370 
371  auto set_time_monitor = [&]() {
373  boost::shared_ptr<Monitor> monitor_ptr(
374  new Monitor(dm, postProcFe, uXScatter, uYScatter));
375  boost::shared_ptr<ForcesAndSourcesCore> null;
376  CHKERR DMMoFEMTSSetMonitor(dm, solver, simple->getDomainFEName(),
377  monitor_ptr, null, null);
379  };
380 
381  CHKERR set_section_monitor();
382  CHKERR create_post_process_element();
383  uXScatter = scatter_create(0);
384  uYScatter = scatter_create(1);
385  CHKERR set_time_monitor();
386 
387  CHKERR TSSolve(solver, D);
388 
389  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
390  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
391  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
392 
394 }
395 //! [Solve]
396 
397 //! [Postprocess results]
398 MoFEMErrorCode Example::postProcess() { return 0; }
399 //! [Postprocess results]
400 
401 //! [Check]
403 //! [Check]
404 
406  Range faces;
407  CHKERR mField.get_moab().get_entities_by_type(0, MBTRI, faces);
408  Skinner skin(&mField.get_moab());
409  Range skin_edges;
410  CHKERR skin.find_skin(0, faces, false, skin_edges);
411  return skin_edges;
412 };
413 
414 static char help[] = "...\n\n";
415 
416 int main(int argc, char *argv[]) {
417 
418  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
419 
420  try {
421 
422  //! [Register MoFEM discrete manager in PETSc]
423  DMType dm_name = "DMMOFEM";
424  CHKERR DMRegister_MoFEM(dm_name);
425  //! [Register MoFEM discrete manager in PETSc
426 
427  //! [Create MoAB]
428  moab::Core mb_instance; ///< mesh database
429  moab::Interface &moab = mb_instance; ///< mesh database interface
430  //! [Create MoAB]
431 
432  //! [Create MoFEM]
433  MoFEM::Core core(moab); ///< finite element database
434  MoFEM::Interface &m_field = core; ///< finite element database insterface
435  //! [Create MoFEM]
436 
437  //! [Load mesh]
438  Simple *simple = m_field.getInterface<Simple>();
439  CHKERR simple->getOptions();
440  CHKERR simple->loadFile("");
441  //! [Load mesh]
442 
443  //! [Example]
444  Example ex(m_field);
445  CHKERR ex.runProblem();
446  //! [Example]
447  }
448  CATCH_ERRORS;
449 
451 }
constexpr double poisson_ratio
constexpr double spring_stiffness
boost::ptr_vector< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
Deprecated interface functions.
Problem manager is used to build and partition problems.
Range getEntsOnMeshSkin()
[Check]
boost::shared_ptr< PostProcFaceOnRefinedMeshFor2D > postProcFe
MoFEMErrorCode createCommonData()
[Set up problem]
MoFEM::EdgeElementForcesAndSourcesCoreSwitch< EdgeElementForcesAndSourcesCore::NO_COVARIANT_TRANSFORM_HCURL > EdgeEle2D
Definition: ContactOps.hpp:856
Section manager is used to create indexes and sectionsFIXME: ISManager is not properly testsed by ato...
Definition: ISManager.hpp:36
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:445
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
static constexpr int approx_order
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:75
MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule)
MatrixDouble invJac
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
Example(MoFEM::Interface &m_field)
boost::shared_ptr< OpContactTools::CommonData > commonDataPtr
Core (interface) class.
Definition: Core.hpp:50
Get values at integration pts for tensor filed rank 1, i.e. vector field.
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
int main(int argc, char *argv[])
Simple interface for fast problem set-up.
Definition: Simple.hpp:36
boost::ptr_vector< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
constexpr double young_modulus
MoFEMErrorCode tsSolve()
[Push operators to pipeline]
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:51
Edge finite elementUser is implementing own operator at Gauss points level, by own object derived fro...
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:939
MoFEMErrorCode runProblem()
[Operator]
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
MoFEMErrorCode OPs()
[Boundary condition]
MoFEMErrorCode checkResults()
[Print results]
Calculate inverse of jacobian for face element.
MoFEMErrorCode bC()
[Create common data]
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:48
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
MatrixDouble jAc
field with continuous tangents
Definition: definitions.h:178
Calculate divergence of tonsorial field using vectorial base.
static char help[]
PipelineManager interface.
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:302
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMMoFEM.cpp:915
MoFEMErrorCode setUP()
[Set up problem]
Transform local reference derivatives of shape functions to global derivatives.
FTensor::Index< 'k', 2 > k
Definition: ContactOps.hpp:28
constexpr int order
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:291
MoFEMErrorCode postProcess()
[Integrate]
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
brief Transform local reference derivatives of shape function to global derivatives for face
boost::ptr_vector< UserDataOperator > & getOpBoundaryLhsPipeline()
Get the Op Boundary Lhs Pipeline object.
Definition: ElasticOps.hpp:86
MoFEM::FaceElementForcesAndSourcesCoreSwitch< FaceElementForcesAndSourcesCore::NO_CONTRAVARIANT_TRANSFORM_HDIV|FaceElementForcesAndSourcesCore::NO_COVARIANT_TRANSFORM_HCURL > FaceEle2D
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
Data on single entity (This is passed as argument to DataOperator::doWork)
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
Apply contravariant (Piola) transfer to Hdiv space on face.
FTensor::Index< 'j', 2 > j
Definition: ContactOps.hpp:27
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:439
Calculate jacobian for face element.
continuous field
Definition: definitions.h:177
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
constexpr double cn
intrusive_ptr for managing petsc objects
Definition: AuxPETSc.hpp:128
boost::ptr_vector< UserDataOperator > & getOpBoundaryRhsPipeline()
Get the Op Boundary Rhs Pipeline object.
DataForcesAndSourcesCore::EntData EntData
FTensor::Index< 'l', 2 > l
Definition: ContactOps.hpp:29
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:69
Face finite elementUser is implementing own operator at Gauss point level, by own object derived from...
Make Hdiv space from Hcurl space in 2d.
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:26