v0.13.0
nonlinear_elastic.cpp
Go to the documentation of this file.
1 /**
2  * \file nonlinear_elastic.cpp
3  * \example nonlinear_elastic.cpp
4  *
5  * Plane stress elastic dynamic 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 #include <MatrixFunction.hpp>
25 
26 using namespace MoFEM;
27 
28 template <int DIM> struct ElementsAndOps {};
29 
30 template <> struct ElementsAndOps<2> {
36 };
37 
38 template <> struct ElementsAndOps<3> {
44 };
45 
46 constexpr int SPACE_DIM = 3; //< Space dimension of problem, mesh
47 
54 
56  GAUSS>::OpGradTensorGrad<1, SPACE_DIM, SPACE_DIM, 1>;
59 
61  GAUSS>::OpSource<1, SPACE_DIM>;
62 
65 
66 constexpr double young_modulus = 100;
67 constexpr double poisson_ratio = 0.3;
68 constexpr double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
69 constexpr double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
70 
71 #include <HenckyOps.hpp>
72 using namespace HenckyOps;
73 
74 struct Example {
75 
76  Example(MoFEM::Interface &m_field) : mField(m_field) {}
77 
79 
80 private:
81  MoFEM::Interface &mField;
82 
91 
92  boost::shared_ptr<MatrixDouble> matGradPtr;
93  boost::shared_ptr<MatrixDouble> matStrainPtr;
94  boost::shared_ptr<MatrixDouble> matStressPtr;
95  boost::shared_ptr<MatrixDouble> matDPtr;
96  boost::shared_ptr<HenckyOps::CommonData> commonHenckyDataPtr;
97 };
98 
99 //! [Create common data]
102 
103  auto set_material_stiffness = [&]() {
110  constexpr double A =
111  (SPACE_DIM == 2) ? 2 * shear_modulus_G /
112  (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
113  : 1;
114  auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*matDPtr);
115  t_D(i, j, k, l) = 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
116  A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) *
117  t_kd(i, j) * t_kd(k, l);
119  };
120 
121  matGradPtr = boost::make_shared<MatrixDouble>();
122  matStrainPtr = boost::make_shared<MatrixDouble>();
123  matStressPtr = boost::make_shared<MatrixDouble>();
124  matDPtr = boost::make_shared<MatrixDouble>();
125 
126  commonHenckyDataPtr = boost::make_shared<HenckyOps::CommonData>();
127  commonHenckyDataPtr->matGradPtr = matGradPtr;
128  commonHenckyDataPtr->matDPtr = matDPtr;
129 
130  constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
131  matDPtr->resize(size_symm * size_symm, 1);
132 
133  CHKERR set_material_stiffness();
134 
136 }
137 //! [Create common data]
138 
139 //! [Run problem]
142  CHKERR readMesh();
143  CHKERR setupProblem();
144  CHKERR createCommonData();
145  CHKERR boundaryCondition();
146  CHKERR assembleSystem();
147  CHKERR solveSystem();
148  CHKERR outputResults();
149  CHKERR checkResults();
151 }
152 //! [Run problem]
153 
154 //! [Read mesh]
157  auto simple = mField.getInterface<Simple>();
158  CHKERR simple->getOptions();
159  CHKERR simple->loadFile();
161 }
162 //! [Read mesh]
163 
164 //! [Set up problem]
167  Simple *simple = mField.getInterface<Simple>();
168  // Add field
169  CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, SPACE_DIM);
170  CHKERR simple->addBoundaryField("U", H1, AINSWORTH_LEGENDRE_BASE, SPACE_DIM);
171  int order = 2;
172  CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
173  CHKERR simple->setFieldOrder("U", order);
174  CHKERR simple->setUp();
176 }
177 //! [Set up problem]
178 
179 //! [Boundary condition]
182  auto *pipeline_mng = mField.getInterface<PipelineManager>();
183  auto simple = mField.getInterface<Simple>();
184  auto bc_mng = mField.getInterface<BcManager>();
185 
186  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_X",
187  "U", 0, 0);
188  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_Y",
189  "U", 1, 1);
190  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_Z",
191  "U", 2, 2);
192  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_ALL",
193  "U", 0, 3);
194 
195  auto get_time = [&](double, double, double) {
196  auto *pipeline_mng = mField.getInterface<PipelineManager>();
197  auto &fe_domain_rhs = pipeline_mng->getDomainRhsFE();
198  return fe_domain_rhs->ts_t;
199  };
200 
202  if (it->getName().compare(0, 5, "FORCE") == 0) {
203  Range force_edges;
204  std::vector<double> attr_vec;
205  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 1,
206  force_edges, true);
207  it->getAttributes(attr_vec);
208  if (attr_vec.size() < SPACE_DIM)
209  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
210  "Wrong size of boundary attributes vector. Set right block "
211  "size attributes.");
212  auto force_vec_ptr = boost::make_shared<MatrixDouble>(SPACE_DIM, 1);
213  std::copy(&attr_vec[0], &attr_vec[SPACE_DIM],
214  force_vec_ptr->data().begin());
215  pipeline_mng->getOpBoundaryRhsPipeline().push_back(
216  new OpBoundaryVec("U", force_vec_ptr, get_time,
217  boost::make_shared<Range>(force_edges)));
218  }
219  }
220 
221  //! [Define gravity vector]
222  auto get_body_force = [this](const double, const double, const double) {
223  auto *pipeline_mng = mField.getInterface<PipelineManager>();
226  auto fe_domain_rhs = pipeline_mng->getDomainRhsFE();
227  const auto time = fe_domain_rhs->ts_t;
228  // hardcoded gravity load in y direction
229  t_source(i) = 0;
230  t_source(1) = 1 * time;
231  return t_source;
232  };
233  //! [Define gravity vector]
234 
235  pipeline_mng->getOpDomainRhsPipeline().push_back(
236  new OpBodyForce("U", get_body_force));
237 
239 }
240 //! [Boundary condition]
241 
242 //! [Push operators to pipeline]
245  auto *simple = mField.getInterface<Simple>();
246  auto *pipeline_mng = mField.getInterface<PipelineManager>();
247 
248  auto add_domain_base_ops = [&](auto &pipeline) {
250  if (SPACE_DIM == 2) {
251  auto det_ptr = boost::make_shared<VectorDouble>();
252  auto jac_ptr = boost::make_shared<MatrixDouble>();
253  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
254  pipeline.push_back(new OpCalculateHOJacForFace(jac_ptr));
255  pipeline.push_back(new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
256  pipeline.push_back(new OpSetInvJacH1ForFace(inv_jac_ptr));
257  }
258 
260  "U", matGradPtr));
261  pipeline.push_back(
262  new OpCalculateEigenVals<SPACE_DIM>("U", commonHenckyDataPtr));
263  pipeline.push_back(
264  new OpCalculateLogC<SPACE_DIM>("U", commonHenckyDataPtr));
265  pipeline.push_back(
266  new OpCalculateLogC_dC<SPACE_DIM>("U", commonHenckyDataPtr));
267  pipeline.push_back(
268  new OpCalculateHenckyStress<SPACE_DIM>("U", commonHenckyDataPtr));
269  pipeline.push_back(
270  new OpCalculatePiolaStress<SPACE_DIM>("U", commonHenckyDataPtr));
271 
273  };
274 
275  auto add_domain_ops_lhs = [&](auto &pipeline) {
277 
278  pipeline.push_back(
279  new OpHenckyTangent<SPACE_DIM>("U", commonHenckyDataPtr));
280  pipeline.push_back(new OpK("U", "U", commonHenckyDataPtr->getMatTangent()));
281 
283  };
284 
285  auto add_domain_ops_rhs = [&](auto &pipeline) {
287  // Calculate internal forece
288  pipeline.push_back(new OpInternalForce(
289  "U", commonHenckyDataPtr->getMatFirstPiolaStress()));
290 
292  };
293 
294  CHKERR add_domain_base_ops(pipeline_mng->getOpDomainLhsPipeline());
295  CHKERR add_domain_ops_lhs(pipeline_mng->getOpDomainLhsPipeline());
296  CHKERR add_domain_base_ops(pipeline_mng->getOpDomainRhsPipeline());
297  CHKERR add_domain_ops_rhs(pipeline_mng->getOpDomainRhsPipeline());
298 
299  auto integration_rule = [](int, int, int approx_order) {
300  return 2 * (approx_order - 1);
301  };
302 
303  CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
304  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
305  CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule);
306 
308 }
309 //! [Push operators to pipeline]
310 
311 /**
312  * @brief Monitor solution
313  *
314  * This functions is called by TS solver at the end of each step. It is used
315  * to output results to the hard drive.
316  */
317 struct Monitor : public FEMethod {
318  Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc)
319  : dM(dm), postProc(post_proc){};
322  constexpr int save_every_nth_step = 1;
323  if (ts_step % save_every_nth_step == 0) {
324  CHKERR DMoFEMLoopFiniteElements(dM, "dFE", postProc);
325  CHKERR postProc->writeFile(
326  "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
327  }
329  }
330 
331 private:
333  boost::shared_ptr<PostProcEle> postProc;
334 };
335 
336 //! [Solve]
339  auto *simple = mField.getInterface<Simple>();
340  auto *pipeline_mng = mField.getInterface<PipelineManager>();
341 
342  auto dm = simple->getDM();
343  auto ts = pipeline_mng->createTSIM();
344 
345  // Setup postprocessing
346  auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
347  post_proc_fe->generateReferenceElementMesh();
348  if (SPACE_DIM == 2) {
349  auto det_ptr = boost::make_shared<VectorDouble>();
350  auto jac_ptr = boost::make_shared<MatrixDouble>();
351  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
352  post_proc_fe->getOpPtrVector().push_back(
353  new OpCalculateHOJacForFace(jac_ptr));
354  post_proc_fe->getOpPtrVector().push_back(
355  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
356  post_proc_fe->getOpPtrVector().push_back(
357  new OpSetInvJacH1ForFace(inv_jac_ptr));
358  }
359  post_proc_fe->getOpPtrVector().push_back(
361  "U", commonHenckyDataPtr->matGradPtr));
362  post_proc_fe->getOpPtrVector().push_back(
363  new OpCalculateEigenVals<SPACE_DIM>("U", commonHenckyDataPtr));
364  post_proc_fe->getOpPtrVector().push_back(
365  new OpCalculateLogC<SPACE_DIM>("U", commonHenckyDataPtr));
366  post_proc_fe->getOpPtrVector().push_back(
367  new OpCalculateLogC_dC<SPACE_DIM>("U", commonHenckyDataPtr));
368  post_proc_fe->getOpPtrVector().push_back(
369  new OpCalculateHenckyStress<SPACE_DIM>("U", commonHenckyDataPtr));
370  post_proc_fe->getOpPtrVector().push_back(
371  new OpCalculatePiolaStress<SPACE_DIM>("U", commonHenckyDataPtr));
372 
373  post_proc_fe->getOpPtrVector().push_back(new OpPostProcHencky<SPACE_DIM>(
374  "U", post_proc_fe->postProcMesh, post_proc_fe->mapGaussPts,
375  commonHenckyDataPtr));
376 
377  post_proc_fe->addFieldValuesPostProc("U");
378 
379  // Add monitor to time solver
380  boost::shared_ptr<FEMethod> null_fe;
381  auto monitor_ptr = boost::make_shared<Monitor>(dm, post_proc_fe);
382  CHKERR DMMoFEMTSSetMonitor(dm, ts, simple->getDomainFEName(), null_fe,
383  null_fe, monitor_ptr);
384 
385  double ftime = 1;
386  CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
387  CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
388 
389  auto D = smartCreateDMVector(simple->getDM());
390 
391  CHKERR TSSetSolution(ts, D);
392  CHKERR TSSetFromOptions(ts);
393 
394  CHKERR TSSolve(ts, NULL);
395  CHKERR TSGetTime(ts, &ftime);
396 
397  PetscInt steps, snesfails, rejects, nonlinits, linits;
398  CHKERR TSGetTimeStepNumber(ts, &steps);
399  CHKERR TSGetSNESFailures(ts, &snesfails);
400  CHKERR TSGetStepRejections(ts, &rejects);
401  CHKERR TSGetSNESIterations(ts, &nonlinits);
402  CHKERR TSGetKSPIterations(ts, &linits);
403  MOFEM_LOG_C("EXAMPLE", Sev::inform,
404  "steps %D (%D rejected, %D SNES fails), ftime %g, nonlinits "
405  "%D, linits %D\n",
406  steps, rejects, snesfails, ftime, nonlinits, linits);
407 
409 }
410 //! [Solve]
411 
412 //! [Postprocess results]
415  PetscBool test_flg = PETSC_FALSE;
416  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-test", &test_flg, PETSC_NULL);
417  if (test_flg) {
418  auto *simple = mField.getInterface<Simple>();
419  auto T = smartCreateDMVector(simple->getDM());
420  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
421  SCATTER_FORWARD);
422  double nrm2;
423  CHKERR VecNorm(T, NORM_2, &nrm2);
424  MOFEM_LOG("EXAMPLE", Sev::inform) << "Regression norm " << nrm2;
425  constexpr double regression_value = 1.09572;
426  if (fabs(nrm2 - regression_value) > 1e-2)
427  SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
428  "Regression test faileed; wrong norm value.");
429  }
431 }
432 //! [Postprocess results]
433 
434 //! [Check]
438 }
439 //! [Check]
440 
441 static char help[] = "...\n\n";
442 
443 int main(int argc, char *argv[]) {
444 
445  // Initialisation of MoFEM/PETSc and MOAB data structures
446  const char param_file[] = "param_file.petsc";
447  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
448 
449  // Add logging channel for example
450  auto core_log = logging::core::get();
451  core_log->add_sink(
453  LogManager::setLog("EXAMPLE");
454  MOFEM_LOG_TAG("EXAMPLE", "example");
455 
456  try {
457 
458  //! [Register MoFEM discrete manager in PETSc]
459  DMType dm_name = "DMMOFEM";
460  CHKERR DMRegister_MoFEM(dm_name);
461  //! [Register MoFEM discrete manager in PETSc
462 
463  //! [Create MoAB]
464  moab::Core mb_instance; ///< mesh database
465  moab::Interface &moab = mb_instance; ///< mesh database interface
466  //! [Create MoAB]
467 
468  //! [Create MoFEM]
469  MoFEM::Core core(moab); ///< finite element database
470  MoFEM::Interface &m_field = core; ///< finite element database insterface
471  //! [Create MoFEM]
472 
473  //! [Example]
474  Example ex(m_field);
475  CHKERR ex.runProblem();
476  //! [Example]
477  }
478  CATCH_ERRORS;
479 
481 }
std::string param_file
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
int main(int argc, char *argv[])
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Kronecker Delta class symmetric.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, SPACE_DIM > OpBodyForce
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ H1
continuous field
Definition: definitions.h:98
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ BLOCKSET
Definition: definitions.h:161
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:53
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
auto integration_rule
constexpr auto t_kd
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:956
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:481
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:544
FaceElementForcesAndSourcesCoreSwitch< 0 > FaceElementForcesAndSourcesCore
Face finite element default.
VolumeElementForcesAndSourcesCoreSwitch< 0 > VolumeElementForcesAndSourcesCore
Volume finite element default.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
#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.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
const double T
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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:994
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
CoreTmp< 0 > Core
Definition: Core.hpp:1096
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
const double D
diffusivity
double A
int save_every_nth_step
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 0 > OpBoundaryVec
Definition: plastic.cpp:95
static constexpr int approx_order
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
PipelineManager::FaceEle DomainEle
[Example]
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
MoFEMErrorCode checkResults()
MoFEMErrorCode solveSystem()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
MoFEMErrorCode runProblem()
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
Definition: HenckyOps.hpp:512
Simple interface for fast problem set-up.
Definition: BcManager.hpp:27
Core (interface) class.
Definition: Core.hpp:92
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:85
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:125
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
structure for User Loop Methods on finite elements
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
MoFEM::FaceElementForcesAndSourcesCore FaceEle
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
Simple interface for fast problem set-up.
Definition: Simple.hpp:33
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
VolumeElementForcesAndSourcesCoreBase::UserDataOperator UserDataOperator
[Push operators to pipeline]
Monitor(SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
MoFEMErrorCode postProcess()
function is run at the end of loop
Postprocess on face.
Post processing.
EntitiesFieldData::EntData EntData
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
static char help[]
[Check]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradTensorGrad< 1, SPACE_DIM, SPACE_DIM, 1 > OpK
constexpr int SPACE_DIM
constexpr double poisson_ratio
constexpr double shear_modulus_G
constexpr double bulk_modulus_K
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 0 > OpBoundaryVec
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, SPACE_DIM > OpBodyForce
constexpr double young_modulus