v0.9.1
basic_contact.cpp
Go to the documentation of this file.
1 /**
2  * \file basic_contact.cpp
3  * \example basic_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 
61  MatrixDouble invJac, jAc;
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  CHKERR simple->setFieldOrder("SIGMA", order - 1, &skin_edges);
100  // CHKERR simple->setFieldOrder("U", order + 1, &skin_edges);
101 
102  CHKERR simple->setUp();
104 }
105 //! [Set up problem]
106 
107 //! [Create common data]
110 
111  auto get_matrial_stiffens = [&](FTensor::Ddg<double, 2, 2> &t_D) {
117  t_D(i, j, k, l) = 0;
118 
119  constexpr double c = young_modulus / (1 - poisson_ratio * poisson_ratio);
120  constexpr double o = poisson_ratio * c;
121 
122  t_D(0, 0, 0, 0) = c;
123  t_D(0, 0, 1, 1) = o;
124 
125  t_D(1, 1, 0, 0) = o;
126  t_D(1, 1, 1, 1) = c;
127 
128  t_D(0, 1, 0, 1) = (1 - poisson_ratio) * c;
130  };
131 
132  commonDataPtr = boost::make_shared<OpContactTools::CommonData>();
133  CHKERR get_matrial_stiffens(commonDataPtr->tD);
134 
135  commonDataPtr->mGradPtr = boost::make_shared<MatrixDouble>();
136  commonDataPtr->mStrainPtr = boost::make_shared<MatrixDouble>();
137  commonDataPtr->mStressPtr = boost::make_shared<MatrixDouble>();
138  commonDataPtr->contactStressPtr = boost::make_shared<MatrixDouble>();
139  commonDataPtr->contactStressDivergencePtr =
140  boost::make_shared<MatrixDouble>();
141  commonDataPtr->contactTractionPtr = boost::make_shared<MatrixDouble>();
142  commonDataPtr->contactDispPtr = boost::make_shared<MatrixDouble>();
143 
144  jAc.resize(2, 2, false);
145  invJac.resize(2, 2, false);
146 
148 }
149 //! [Create common data]
150 
151 //! [Boundary condition]
154 
155  auto fix_disp = [&](const std::string blockset_name) {
156  Range fix_ents;
158  if (it->getName().compare(0, blockset_name.length(), blockset_name) ==
159  0) {
160  CHKERR mField.get_moab().get_entities_by_handle(it->meshset, fix_ents,
161  true);
162  }
163  }
164  return fix_ents;
165  };
166 
167  auto remove_ents = [&](const Range &&ents, const bool fix_x,
168  const bool fix_y) {
169  auto prb_mng = mField.getInterface<ProblemsManager>();
170  auto simple = mField.getInterface<Simple>();
172  Range verts;
173  CHKERR mField.get_moab().get_connectivity(ents, verts, true);
174  verts.merge(ents);
175  const int lo_coeff = fix_x ? 0 : 1;
176  const int hi_coeff = fix_y ? 1 : 0;
177  CHKERR prb_mng->removeDofsOnEntities(simple->getProblemName(), "U", verts,
178  lo_coeff, hi_coeff);
180  };
181 
182  Range boundary_ents;
183  boundary_ents.merge(fix_disp("FIX_X"));
184  boundary_ents.merge(fix_disp("FIX_Y"));
185  boundary_ents.merge(fix_disp("FIX_ALL"));
186  CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
187  mField.getInterface<Simple>()->getProblemName(), "SIGMA", boundary_ents,
188  0, 2);
189 
190  CHKERR remove_ents(fix_disp("FIX_X"), true, false);
191  CHKERR remove_ents(fix_disp("FIX_Y"), false, true);
192  CHKERR remove_ents(fix_disp("FIX_ALL"), true, true);
193 
195 }
196 //! [Boundary condition]
197 
198 //! [Push operators to pipeline]
201  Basic *basic = mField.getInterface<Basic>();
202 
203  auto add_domain_base_ops = [&](auto &pipeline) {
204  pipeline.push_back(new OpCalculateJacForFace(jAc));
205  pipeline.push_back(new OpCalculateInvJacForFace(invJac));
206  pipeline.push_back(new OpSetInvJacH1ForFace(invJac));
207  pipeline.push_back(new OpMakeHdivFromHcurl());
208  pipeline.push_back(new OpSetContravariantPiolaTransformFace(jAc));
209  pipeline.push_back(new OpSetInvJacHcurlFace(invJac));
210  };
211 
212  auto add_domain_ops_lhs = [&](auto &pipeline) {
213  pipeline.push_back(new OpStiffnessMatrixLhs("U", "U", commonDataPtr));
214  pipeline.push_back(
215  new OpConstrainDomainLhs_dU("SIGMA", "U", commonDataPtr));
216  };
217 
218  auto add_domain_ops_rhs = [&](auto &pipeline) {
219  auto gravity = [](double x, double y) {
220  return FTensor::Tensor1<double, 2>{0., 1.};
221  };
222  pipeline.push_back(new OpForceRhs("U", commonDataPtr, gravity));
223 
224  pipeline.push_back(
225  new OpCalculateVectorFieldGradient<2, 2>("U", commonDataPtr->mGradPtr));
226  pipeline.push_back(new OpStrain("U", commonDataPtr));
227  pipeline.push_back(new OpStress("U", commonDataPtr));
228  pipeline.push_back(new OpInternalForceRhs("U", commonDataPtr));
229 
230  pipeline.push_back(new OpCalculateVectorFieldValues<2>(
231  "U", commonDataPtr->contactDispPtr));
232 
233  pipeline.push_back(new OpCalculateHVecTensorField<2, 2>(
234  "SIGMA", commonDataPtr->contactStressPtr));
235  pipeline.push_back(new OpCalculateHVecTensorDivergence<2, 2>(
236  "SIGMA", commonDataPtr->contactStressDivergencePtr));
237  pipeline.push_back(new OpConstrainDomainRhs("SIGMA", commonDataPtr));
238  pipeline.push_back(new OpInternalDomainContactRhs("U", commonDataPtr));
239  };
240 
241  auto add_boundary_base_ops = [&](auto &pipeline) {
242  pipeline.push_back(new OpSetContrariantPiolaTransformOnEdge());
243  pipeline.push_back(new OpCalculateVectorFieldValues<2>(
244  "U", commonDataPtr->contactDispPtr));
245  pipeline.push_back(new OpConstrainBoundaryTraction("SIGMA", commonDataPtr));
246  };
247 
248  auto add_boundary_ops_lhs = [&](auto &pipeline) {
249  pipeline.push_back(
250  new OpConstrainBoundaryLhs_dU("SIGMA", "U", commonDataPtr));
251  pipeline.push_back(
252  new OpConstrainBoundaryLhs_dTraction("SIGMA", "SIGMA", commonDataPtr));
253  pipeline.push_back(new OpSpringLhs("U", "U"));
254  };
255 
256  auto add_boundary_ops_rhs = [&](auto &pipeline) {
257  pipeline.push_back(new OpConstrainBoundaryRhs("SIGMA", commonDataPtr));
258  pipeline.push_back(new OpSpringRhs("U", commonDataPtr));
259  };
260 
261  add_domain_base_ops(basic->getOpDomainLhsPipeline());
262  add_domain_base_ops(basic->getOpDomainRhsPipeline());
263  add_domain_ops_lhs(basic->getOpDomainLhsPipeline());
264  add_domain_ops_rhs(basic->getOpDomainRhsPipeline());
265 
266  add_boundary_base_ops(basic->getOpBoundaryLhsPipeline());
267  add_boundary_base_ops(basic->getOpBoundaryRhsPipeline());
268  add_boundary_ops_lhs(basic->getOpBoundaryLhsPipeline());
269  add_boundary_ops_rhs(basic->getOpBoundaryRhsPipeline());
270 
271  auto integration_rule = [](int, int, int approx_order) { return 2 * order; };
272  CHKERR basic->setDomainRhsIntegrationRule(integration_rule);
273  CHKERR basic->setDomainLhsIntegrationRule(integration_rule);
274  CHKERR basic->setBoundaryRhsIntegrationRule(integration_rule);
275  CHKERR basic->setBoundaryLhsIntegrationRule(integration_rule);
276 
278 }
279 //! [Push operators to pipeline]
280 
281 //! [Solve]
284 
285  Simple *simple = mField.getInterface<Simple>();
286  Basic *basic = mField.getInterface<Basic>();
287  ISManager *is_manager = mField.getInterface<ISManager>();
288 
289  auto solver = basic->createTS();
290 
291  auto dm = simple->getDM();
292  auto D = smartCreateDMVector(dm);
293 
294  CHKERR TSSetSolution(solver, D);
295  CHKERR TSSetFromOptions(solver);
296  CHKERR TSSetUp(solver);
297 
298  auto set_section_monitor = [&]() {
300  SNES snes;
301  CHKERR TSGetSNES(solver, &snes);
302  PetscViewerAndFormat *vf;
303  CHKERR PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,
304  PETSC_VIEWER_DEFAULT, &vf);
305  CHKERR SNESMonitorSet(
306  snes,
307  (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal, void *))SNESMonitorFields,
308  vf, (MoFEMErrorCode(*)(void **))PetscViewerAndFormatDestroy);
310  };
311 
312  auto create_post_process_element = [&]() {
314  postProcFe = boost::make_shared<PostProcFaceOnRefinedMeshFor2D>(mField);
315  postProcFe->generateReferenceElementMesh();
316 
317  postProcFe->getOpPtrVector().push_back(new OpCalculateJacForFace(jAc));
318  postProcFe->getOpPtrVector().push_back(
319  new OpCalculateInvJacForFace(invJac));
320  postProcFe->getOpPtrVector().push_back(new OpSetInvJacH1ForFace(invJac));
321  postProcFe->getOpPtrVector().push_back(new OpMakeHdivFromHcurl());
322  postProcFe->getOpPtrVector().push_back(
324  postProcFe->getOpPtrVector().push_back(new OpSetInvJacHcurlFace(invJac));
325 
326  postProcFe->getOpPtrVector().push_back(
327  new OpCalculateVectorFieldGradient<2, 2>("U", commonDataPtr->mGradPtr));
328  postProcFe->getOpPtrVector().push_back(new OpStrain("U", commonDataPtr));
329  postProcFe->getOpPtrVector().push_back(new OpStress("U", commonDataPtr));
330  postProcFe->getOpPtrVector().push_back(
332  "SIGMA", commonDataPtr->contactStressDivergencePtr));
333  postProcFe->getOpPtrVector().push_back(new OpCalculateHVecTensorField<2, 2>(
334  "SIGMA", commonDataPtr->contactStressPtr));
335 
336  postProcFe->getOpPtrVector().push_back(new OpPostProcElastic(
337  "U", postProcFe->postProcMesh, postProcFe->mapGaussPts, commonDataPtr));
338  postProcFe->getOpPtrVector().push_back(
339  new OpPostProcContact("SIGMA", postProcFe->postProcMesh,
340  postProcFe->mapGaussPts, commonDataPtr));
341  postProcFe->addFieldValuesPostProc("U");
343  };
344 
345  auto scatter_create = [&](auto coeff) {
347  CHKERR is_manager->isCreateProblemFieldAndRank(simple->getProblemName(),
348  ROW, "U", coeff, coeff, is);
349  int loc_size;
350  CHKERR ISGetLocalSize(is, &loc_size);
351  Vec v;
352  CHKERR VecCreateMPI(mField.get_comm(), loc_size, PETSC_DETERMINE, &v);
353  VecScatter scatter;
354  CHKERR VecScatterCreate(D, is, v, PETSC_NULL, &scatter);
355  return std::make_tuple(SmartPetscObj<Vec>(v),
356  SmartPetscObj<VecScatter>(scatter));
357  };
358 
359  auto set_time_monitor = [&]() {
361  boost::shared_ptr<Monitor> monitor_ptr(
362  new Monitor(dm, postProcFe, uXScatter, uYScatter));
363  boost::shared_ptr<ForcesAndSourcesCore> null;
364  CHKERR DMMoFEMTSSetMonitor(dm, solver, simple->getDomainFEName(),
365  monitor_ptr, null, null);
367  };
368 
369  CHKERR set_section_monitor();
370  CHKERR create_post_process_element();
371  uXScatter = scatter_create(0);
372  uYScatter = scatter_create(1);
373  CHKERR set_time_monitor();
374 
375  CHKERR TSSolve(solver, D);
376 
377  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
378  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
379  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
380 
382 }
383 //! [Solve]
384 
385 //! [Postprocess results]
386 MoFEMErrorCode Example::postProcess() { return 0; }
387 //! [Postprocess results]
388 
389 //! [Check]
391 //! [Check]
392 
394  Range faces;
395  CHKERR mField.get_moab().get_entities_by_type(0, MBTRI, faces);
396  Skinner skin(&mField.get_moab());
397  Range skin_edges;
398  CHKERR skin.find_skin(0, faces, false, skin_edges);
399  return skin_edges;
400 };
401 
402 static char help[] = "...\n\n";
403 
404 int main(int argc, char *argv[]) {
405 
406  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
407 
408  try {
409 
410  //! [Register MoFEM discrete manager in PETSc]
411  DMType dm_name = "DMMOFEM";
412  CHKERR DMRegister_MoFEM(dm_name);
413  //! [Register MoFEM discrete manager in PETSc
414 
415  //! [Create MoAB]
416  moab::Core mb_instance; ///< mesh database
417  moab::Interface &moab = mb_instance; ///< mesh database interface
418  //! [Create MoAB]
419 
420  //! [Create MoFEM]
421  MoFEM::Core core(moab); ///< finite element database
422  MoFEM::Interface &m_field = core; ///< finite element database insterface
423  //! [Create MoFEM]
424 
425  //! [Load mesh]
426  Simple *simple = m_field.getInterface<Simple>();
427  CHKERR simple->getOptions();
428  CHKERR simple->loadFile("");
429  //! [Load mesh]
430 
431  //! [Example]
432  Example ex(m_field);
433  CHKERR ex.runProblem();
434  //! [Example]
435  }
436  CATCH_ERRORS;
437 
439 }
FTensor::Index< 'i', 3 > i
Deprecated interface functions.
FTensor::Index< 'j', 3 > j
Problem manager is used to build and partition problems.
Range getEntsOnMeshSkin()
[Check]
boost::shared_ptr< PostProcFaceOnRefinedMeshFor2D > postProcFe
MoFEMErrorCode createCommonData()
[Set up problem]
constexpr double poisson_ratio
Definition: ContactOps.hpp:856
MoFEM::FaceElementForcesAndSourcesCoreSwitch< FaceElementForcesAndSourcesCore::NO_CONTRAVARIANT_TRANSFORM_HDIV|FaceElementForcesAndSourcesCore::NO_COVARIANT_TRANSFORM_HCURL > FaceEle2D
Definition: Basic.hpp:48
boost::ptr_vector< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
Definition: Basic.hpp:466
Section manager is used to create indexes and sectionsFIXME: ISManager is not properly testsed by ato...
Definition: ISManager.hpp:36
DataForcesAndSourcesCore::EntData EntData
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:431
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:74
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
Example(MoFEM::Interface &m_field)
boost::shared_ptr< OpContactTools::CommonData > commonDataPtr
boost::ptr_vector< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
Definition: Basic.hpp:442
Core (interface) class.
Definition: Core.hpp:50
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
Definition: Basic.hpp:363
Get values at integration pts for tensor filed rank 1, i.e. vector field.
Edge finite elementUser is implementing own operator at Gauss points level, by own object derived fro...
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
Simple interface for fast problem set-up.
Definition: Simple.hpp:36
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
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:907
constexpr double young_modulus
MoFEMErrorCode runProblem()
[Run problem]
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]
boost::ptr_vector< UserDataOperator > & getOpBoundaryLhsPipeline()
Get the Op Boundary Lhs Pipeline object.
Definition: Basic.hpp:490
Basic interface.
Definition: Basic.hpp:36
MoFEMErrorCode checkResults()
[Postprocess results]
Calculate inverse of jacobian for face element.
MoFEMErrorCode bC()
[Create common data]
MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule)
Definition: Basic.hpp:337
constexpr double cn
constexpr int order
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
static char help[]
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:48
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:151
MatrixDouble jAc
field with continuous tangents
Definition: definitions.h:177
Calculate divergence of tonsorial field using vectorial base.
boost::ptr_vector< UserDataOperator > & getOpBoundaryRhsPipeline()
Get the Op Boundary Rhs Pipeline object.
Definition: Basic.hpp:514
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:299
#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:601
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:901
MoFEMErrorCode setUP()
[Set up problem]
Transform local reference derivatives of shape functions to global derivatives.
constexpr double spring_stiffness
FTensor::Index< 'k', 2 > k
Definition: ContactOps.hpp:28
MoFEMErrorCode postProcess()
[Solve]
brief Transform local reference derivatives of shape function to global derivatives for face
Definition: ElasticOps.hpp:86
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
Data on single entity (This is passed as argument to DataOperator::doWork)
Apply contravariant (Piola) transfer to Hdiv space on face.
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
MoFEM::EdgeElementForcesAndSourcesCoreSwitch< EdgeElementForcesAndSourcesCore::NO_COVARIANT_TRANSFORM_HCURL > EdgeEle2D
Definition: Basic.hpp:50
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Definition: Basic.hpp:285
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
Definition: Basic.hpp:311
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:438
Calculate jacobian for face element.
continuous field
Definition: definitions.h:176
intrusive_ptr for managing petsc objects
Definition: AuxPETSc.hpp:128
int main(int argc, char *argv[])
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:61
Face finite elementUser is implementing own operator at Gauss point level, by own object derived from...
Make Hdiv space from Hcurl space in 2d.