v0.13.1
nonlinear_dynamic_elastic.cpp
Go to the documentation of this file.
1/**
2 * \file dynamic_elastic.cpp
3 * \example dynamic_elastic.cpp
4 *
5 * Plane stress elastic dynamic problem
6 *
7 */
8
9#include <MoFEM.hpp>
10#include <MatrixFunction.hpp>
11using namespace MoFEM;
12
13template <int DIM> struct ElementsAndOps {};
14
15template <> struct ElementsAndOps<2> {
17};
18
19template <> struct ElementsAndOps<3> {
21};
22
23constexpr int SPACE_DIM =
24 EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
25
30
32 GAUSS>::OpGradTensorGrad<1, SPACE_DIM, SPACE_DIM, 1>;
38 GAUSS>::OpSource<1, SPACE_DIM>;
41
42constexpr bool is_quasi_static = false;
43constexpr double rho = 1;
44constexpr double omega = 2.4;
45constexpr double young_modulus = 1;
46constexpr double poisson_ratio = 0.25;
47constexpr double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
48constexpr double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
49
50#include <HenckyOps.hpp>
51using namespace HenckyOps;
52
53static double *ts_time_ptr;
54static double *ts_aa_ptr;
55
56struct Example {
57
58 Example(MoFEM::Interface &m_field) : mField(m_field) {}
59
61
62private:
64
73
74 boost::shared_ptr<MatrixDouble> matGradPtr;
75 boost::shared_ptr<MatrixDouble> matStrainPtr;
76 boost::shared_ptr<MatrixDouble> matStressPtr;
77 boost::shared_ptr<MatrixDouble> matAccelerationPtr;
78 boost::shared_ptr<MatrixDouble> matInertiaPtr;
79 boost::shared_ptr<MatrixDouble> matDPtr;
80
81 boost::shared_ptr<MatrixDouble> matTangentPtr;
82};
83
84//! [Create common data]
87
88 auto set_matrial_stiffness = [&]() {
95 constexpr double A =
96 (SPACE_DIM == 2) ? 2 * shear_modulus_G /
97 (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
98 : 1;
99 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*matDPtr);
100 t_D(i, j, k, l) = 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
101 A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) *
102 t_kd(i, j) * t_kd(k, l);
104 };
105
106 matGradPtr = boost::make_shared<MatrixDouble>();
107 matStrainPtr = boost::make_shared<MatrixDouble>();
108 matStressPtr = boost::make_shared<MatrixDouble>();
109 matAccelerationPtr = boost::make_shared<MatrixDouble>();
110 matInertiaPtr = boost::make_shared<MatrixDouble>();
111 matDPtr = boost::make_shared<MatrixDouble>();
112
113 matTangentPtr = boost::make_shared<MatrixDouble>();
114
115 constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
116 matDPtr->resize(size_symm * size_symm, 1);
117
118 CHKERR set_matrial_stiffness();
119
121}
122//! [Create common data]
123
124//! [Run problem]
136}
137//! [Run problem]
138
139//! [Read mesh]
143 CHKERR simple->getOptions();
144 CHKERR simple->loadFile();
146}
147//! [Read mesh]
148
149//! [Set up problem]
153 // Add field
154 CHKERR simple->addDomainField("U", H1, AINSWORTH_BERNSTEIN_BEZIER_BASE,
155 SPACE_DIM);
156 int order = 2;
157 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
158 CHKERR simple->setFieldOrder("U", order);
159 CHKERR simple->setUp();
161}
162//! [Set up problem]
163
164//! [Boundary condition]
167
169 auto bc_mng = mField.getInterface<BcManager>();
170
171 CHKERR bc_mng->removeBlockDOFsOnEntities<DisplacementCubitBcData>(
172 simple->getProblemName(), "U");
173
175}
176//! [Boundary condition]
177
178//! [Push operators to pipeline]
181 auto *simple = mField.getInterface<Simple>();
182 auto *pipeline_mng = mField.getInterface<PipelineManager>();
183
184 auto det_ptr = boost::make_shared<VectorDouble>();
185 auto jac_ptr = boost::make_shared<MatrixDouble>();
186 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
187
188 auto add_base_transform = [&](auto &pipeline) {
189 pipeline.push_back(new OpCalculateHOJac<SPACE_DIM>(jac_ptr));
190 pipeline.push_back(
191 new OpInvertMatrix<SPACE_DIM>(jac_ptr, det_ptr, inv_jac_ptr));
192 pipeline.push_back(
193 new OpSetHOInvJacToScalarBases<SPACE_DIM>(H1, inv_jac_ptr));
194 pipeline.push_back(new OpSetHOWeights(det_ptr));
195 };
196 add_base_transform(pipeline_mng->getOpDomainRhsPipeline());
197 add_base_transform(pipeline_mng->getOpDomainLhsPipeline());
198
199 // Get pointer to U_tt shift in domain element
200 auto get_rho = [this](const double, const double, const double) {
201 auto *pipeline_mng = mField.getInterface<PipelineManager>();
202 auto &fe_domain_lhs = pipeline_mng->getDomainLhsFE();
203 return rho * fe_domain_lhs->ts_aa;
204 };
205
206 auto get_body_force = [this](const double, const double, const double) {
207 auto *pipeline_mng = mField.getInterface<PipelineManager>();
208 auto fe_domain_rhs = pipeline_mng->getDomainRhsFE();
211 t_source(i) = 0.;
212 t_source(0) = 0.1;
213 t_source(1) = 1.;
214 const auto time = fe_domain_rhs->ts_t;
215 t_source(i) *= sin(time * omega * M_PI);
216 return t_source;
217 };
218
219 auto henky_common_data_ptr = boost::make_shared<HenckyOps::CommonData>();
220 henky_common_data_ptr->matGradPtr = matGradPtr;
221 henky_common_data_ptr->matDPtr = matDPtr;
222
223 pipeline_mng->getOpDomainLhsPipeline().push_back(
225 matGradPtr));
226
227 pipeline_mng->getOpDomainLhsPipeline().push_back(
228 new OpCalculateEigenVals<SPACE_DIM>("U", henky_common_data_ptr));
229 pipeline_mng->getOpDomainLhsPipeline().push_back(
230 new OpCalculateLogC<SPACE_DIM>("U", henky_common_data_ptr));
231 pipeline_mng->getOpDomainLhsPipeline().push_back(
232 new OpCalculateLogC_dC<SPACE_DIM>("U", henky_common_data_ptr));
233 pipeline_mng->getOpDomainLhsPipeline().push_back(
234 new OpCalculateHenckyStress<SPACE_DIM>("U", henky_common_data_ptr));
235 pipeline_mng->getOpDomainLhsPipeline().push_back(
236 new OpCalculatePiolaStress<SPACE_DIM>("U", henky_common_data_ptr));
237 pipeline_mng->getOpDomainLhsPipeline().push_back(
238 new OpHenckyTangent<SPACE_DIM>("U", henky_common_data_ptr));
239 pipeline_mng->getOpDomainLhsPipeline().push_back(
240 new OpK("U", "U", henky_common_data_ptr->getMatTangent()));
241
242 if (!is_quasi_static)
243 pipeline_mng->getOpDomainLhsPipeline().push_back(
244 new OpMass("U", "U", get_rho));
245
246 pipeline_mng->getOpDomainRhsPipeline().push_back(
247 new OpBodyForce("U", get_body_force));
248 pipeline_mng->getOpDomainRhsPipeline().push_back(
250 matGradPtr));
251
252 pipeline_mng->getOpDomainRhsPipeline().push_back(
253 new OpCalculateEigenVals<SPACE_DIM>("U", henky_common_data_ptr));
254 pipeline_mng->getOpDomainRhsPipeline().push_back(
255 new OpCalculateLogC<SPACE_DIM>("U", henky_common_data_ptr));
256 pipeline_mng->getOpDomainRhsPipeline().push_back(
257 new OpCalculateLogC_dC<SPACE_DIM>("U", henky_common_data_ptr));
258 pipeline_mng->getOpDomainRhsPipeline().push_back(
259 new OpCalculateHenckyStress<SPACE_DIM>("U", henky_common_data_ptr));
260 pipeline_mng->getOpDomainRhsPipeline().push_back(
261 new OpCalculatePiolaStress<SPACE_DIM>("U", henky_common_data_ptr));
262
263 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInternalForce(
264 "U", henky_common_data_ptr->getMatFirstPiolaStress()));
265 if (!is_quasi_static) {
266 pipeline_mng->getOpDomainRhsPipeline().push_back(
269 pipeline_mng->getOpDomainRhsPipeline().push_back(
271 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInertiaForce(
272 "U", matInertiaPtr, [](double, double, double) { return 1.; }));
273 }
274
275 auto integration_rule = [](int, int, int approx_order) {
276 return 2 * approx_order;
277 };
278 CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
279 CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
280
281 auto get_bc_hook = [&]() {
283 mField, pipeline_mng->getDomainRhsFE(),
284 {boost::make_shared<TimeScale>()});
285 return hook;
286 };
287
288 pipeline_mng->getDomainRhsFE()->preProcessHook = get_bc_hook();
289
291}
292//! [Push operators to pipeline]
293
294/**
295 * @brief Monitor solution
296 *
297 * This functions is called by TS solver at the end of each step. It is used
298 * to output results to the hard drive.
299 */
300struct Monitor : public FEMethod {
301 Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc)
302 : dM(dm), postProc(post_proc){};
305 constexpr int save_every_nth_step = 1;
306 if (ts_step % save_every_nth_step == 0) {
308 CHKERR postProc->writeFile(
309 "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
310 }
312 }
313
314private:
316 boost::shared_ptr<PostProcEle> postProc;
317};
318
319//! [Solve]
322 auto *simple = mField.getInterface<Simple>();
323 auto *pipeline_mng = mField.getInterface<PipelineManager>();
324
325 auto dm = simple->getDM();
327 if (is_quasi_static)
328 ts = pipeline_mng->createTSIM();
329 else
330 ts = pipeline_mng->createTSIM2();
331
332 // Setup postprocessing
333 auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
334
335 auto det_ptr = boost::make_shared<VectorDouble>();
336 auto jac_ptr = boost::make_shared<MatrixDouble>();
337 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
338 post_proc_fe->getOpPtrVector().push_back(
339 new OpCalculateHOJac<SPACE_DIM>(jac_ptr));
340 post_proc_fe->getOpPtrVector().push_back(
341 new OpInvertMatrix<SPACE_DIM>(jac_ptr, det_ptr, inv_jac_ptr));
342 post_proc_fe->getOpPtrVector().push_back(
343 new OpSetHOInvJacToScalarBases<SPACE_DIM>(H1, inv_jac_ptr));
344
345 post_proc_fe->getOpPtrVector().push_back(
347 matGradPtr));
348 post_proc_fe->getOpPtrVector().push_back(
350 post_proc_fe->getOpPtrVector().push_back(
353
354 auto u_ptr = boost::make_shared<MatrixDouble>();
355 post_proc_fe->getOpPtrVector().push_back(
357
359
360 post_proc_fe->getOpPtrVector().push_back(
361
362 new OpPPMap(
363
364 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
365
366 {},
367
368 {{"U", u_ptr}},
369
370 {},
371
372 {{"STRAIN", matStrainPtr}, {"STRESS", matStressPtr}}
373
374 )
375
376 );
377
378 // Add monitor to time solver
379 boost::shared_ptr<FEMethod> null_fe;
380 auto monitor_ptr = boost::make_shared<Monitor>(dm, post_proc_fe);
381 CHKERR DMMoFEMTSSetMonitor(dm, ts, simple->getDomainFEName(), null_fe,
382 null_fe, monitor_ptr);
383
384 double ftime = 1;
385 // CHKERR TSSetMaxTime(ts, ftime);
386 CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
387
388 auto T = smartCreateDMVector(simple->getDM());
389 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
390 SCATTER_FORWARD);
391 if (is_quasi_static) {
392 CHKERR TSSetSolution(ts, T);
393 CHKERR TSSetFromOptions(ts);
394 } else {
395 auto TT = smartVectorDuplicate(T);
396 CHKERR TS2SetSolution(ts, T, TT);
397 CHKERR TSSetFromOptions(ts);
398 }
399
400 CHKERR TSSolve(ts, NULL);
401 CHKERR TSGetTime(ts, &ftime);
402
403 PetscInt steps, snesfails, rejects, nonlinits, linits;
404#if PETSC_VERSION_GE(3, 8, 0)
405 CHKERR TSGetStepNumber(ts, &steps);
406#else
407 CHKERR TSGetTimeStepNumber(ts, &steps);
408#endif
409 CHKERR TSGetSNESFailures(ts, &snesfails);
410 CHKERR TSGetStepRejections(ts, &rejects);
411 CHKERR TSGetSNESIterations(ts, &nonlinits);
412 CHKERR TSGetKSPIterations(ts, &linits);
413 MOFEM_LOG_C("EXAMPLE", Sev::inform,
414 "steps %D (%D rejected, %D SNES fails), ftime %g, nonlinits "
415 "%D, linits %D\n",
416 steps, rejects, snesfails, ftime, nonlinits, linits);
417
419}
420//! [Solve]
421
422//! [Postprocess results]
425 PetscBool test_flg = PETSC_FALSE;
426 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-test", &test_flg, PETSC_NULL);
427 if (test_flg) {
428 auto *simple = mField.getInterface<Simple>();
429 auto T = smartCreateDMVector(simple->getDM());
430 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
431 SCATTER_FORWARD);
432 double nrm2;
433 CHKERR VecNorm(T, NORM_2, &nrm2);
434 MOFEM_LOG("EXAMPLE", Sev::inform) << "Regression norm " << nrm2;
435 constexpr double regression_value = 1.09572;
436 if (fabs(nrm2 - regression_value) > 1e-2)
437 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
438 "Regression test faileed; wrong norm value.");
439 }
441}
442//! [Postprocess results]
443
444//! [Check]
448}
449//! [Check]
450
451static char help[] = "...\n\n";
452
453int main(int argc, char *argv[]) {
454
455 // Initialisation of MoFEM/PETSc and MOAB data structures
456 const char param_file[] = "param_file.petsc";
458
459 // Add logging channel for example
460 auto core_log = logging::core::get();
461 core_log->add_sink(
462 LogManager::createSink(LogManager::getStrmWorld(), "EXAMPLE"));
463 LogManager::setLog("EXAMPLE");
464 MOFEM_LOG_TAG("EXAMPLE", "example");
465
466 try {
467
468 //! [Register MoFEM discrete manager in PETSc]
469 DMType dm_name = "DMMOFEM";
470 CHKERR DMRegister_MoFEM(dm_name);
471 //! [Register MoFEM discrete manager in PETSc
472
473 //! [Create MoAB]
474 moab::Core mb_instance; ///< mesh database
475 moab::Interface &moab = mb_instance; ///< mesh database interface
476 //! [Create MoAB]
477
478 //! [Create MoFEM]
479 MoFEM::Core core(moab); ///< finite element database
480 MoFEM::Interface &m_field = core; ///< finite element database insterface
481 //! [Create MoFEM]
482
483 //! [Example]
484 Example ex(m_field);
485 CHKERR ex.runProblem();
486 //! [Example]
487 }
489
491}
std::string param_file
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:304
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Kronecker Delta class symmetric.
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
@ H1
continuous field
Definition: definitions.h:85
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
#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
FormsIntegrators< DomainEleOp >::Assembly< A >::LinearForm< I >::OpGradTimesSymTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
Definition: elastic.cpp:44
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
auto integration_rule
constexpr auto t_kd
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:470
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:47
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:533
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:965
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
const double T
FormsIntegrators< DomainEleOp >::Assembly< USER_ASSEMBLE >::LinearForm< GAUSS >::OpBaseTimesVector< 1, 3, 1 > OpInertiaForce
FormsIntegrators< DomainEleOp >::Assembly< USER_ASSEMBLE >::LinearForm< GAUSS >::OpSource< 1, 3 > OpBodyForce
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: MoFEM.hpp:24
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:1003
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)
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
CoreTmp< 0 > Core
Definition: Core.hpp:1086
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
int main(int argc, char *argv[])
EntitiesFieldData::EntData EntData
constexpr double rho
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
static char help[]
[Check]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradTensorGrad< 1, SPACE_DIM, SPACE_DIM, 1 > OpK
static double * ts_aa_ptr
constexpr int SPACE_DIM
constexpr double poisson_ratio
constexpr double shear_modulus_G
constexpr double bulk_modulus_K
constexpr double omega
static double * ts_time_ptr
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, SPACE_DIM > OpBodyForce
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 1 > OpInertiaForce
constexpr bool is_quasi_static
constexpr double young_modulus
double A
int save_every_nth_step
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:10
static constexpr int approx_order
[Example]
Definition: plastic.cpp:120
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
boost::shared_ptr< MatrixDouble > matAccelerationPtr
boost::shared_ptr< MatrixDouble > matInertiaPtr
boost::shared_ptr< MatrixDouble > matGradPtr
boost::shared_ptr< MatrixDouble > matStressPtr
boost::shared_ptr< MatrixDouble > matStrainPtr
MoFEMErrorCode checkResults()
MoFEMErrorCode solveSystem()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
boost::shared_ptr< MatrixDouble > matTangentPtr
MoFEMErrorCode runProblem()
MoFEM::Interface & mField
Definition: plastic.cpp:127
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
boost::shared_ptr< MatrixDouble > matDPtr
Simple interface for fast problem set-up.
Definition: BcManager.hpp:23
Core (interface) class.
Definition: Core.hpp:82
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
Deprecated interface functions.
Definition of the displacement bc data structure.
Definition: BCData.hpp:72
Data on single entity (This is passed as argument to DataOperator::doWork)
structure for User Loop Methods on finite elements
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Approximate field valuse for given petsc vector.
Get values at integration pts for tensor filed rank 1, i.e. vector field.
Post post-proc data at points from hash maps.
Set inverse jacobian to base functions.
Set inverse jacobian to base functions.
PipelineManager interface.
Simple interface for fast problem set-up.
Definition: Simple.hpp:26
PetscInt ts_step
time step
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Monitor solution.
Monitor(SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
SmartPetscObj< DM > dM
MoFEMErrorCode postProcess()
function is run at the end of loop
boost::shared_ptr< PostProcEle > postProc