v0.13.1
Loading...
Searching...
No Matches
heat_equation.cpp
Go to the documentation of this file.
1/**
2 * \file heat_equation.cpp
3 * \example heat_equation.cpp
4 *
5 * \brief Solve the time-dependent Heat Equation
6 \f[
7 \begin{aligned}
8\frac{\partial u(\mathbf{x}, t)}{\partial t}-\Delta u(\mathbf{x}, t)
9&=f(\mathbf{x}, t) & & \forall \mathbf{x} \in \Omega, t \in(0, T), \\
10u(\mathbf{x}, 0) &=u_{0}(\mathbf{x}) & & \forall \mathbf{x} \in \Omega, \\
11u(\mathbf{x}, t) &=g(\mathbf{x}, t) & & \forall \mathbf{x} \in \partial \Omega,
12t \in(0, T). \end{aligned}
13 \f]
14 **/
15
16
17
18#include <stdlib.h>
19#include <cmath>
21
22using namespace MoFEM;
23
24static char help[] = "...\n\n";
25
26template <int DIM> struct ElementsAndOps {};
27
28//! [Define dimension]
29constexpr int SPACE_DIM = 2; //< Space dimension of problem, mesh
30//! [Define dimension]
31
38
49
56
57// Capacity
58constexpr double c = 1;
59constexpr double k = 1;
60constexpr double init_u = 0.;
61
62/**
63 * @brief Monitor solution
64 *
65 * This functions is called by TS solver at the end of each step. It is used
66 * to output results to the hard drive.
67 */
68struct Monitor : public FEMethod {
69
70 Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc)
71 : dM(dm), postProc(post_proc){};
72
73 MoFEMErrorCode preProcess() { return 0; }
74 MoFEMErrorCode operator()() { return 0; }
75
76 static constexpr int saveEveryNthStep = 1;
77
80 if (ts_step % saveEveryNthStep == 0) {
82 CHKERR postProc->writeFile(
83 "out_level_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
84 }
86 }
87
88private:
90 boost::shared_ptr<PostProcEle> postProc;
91};
92
94public:
96
97 // Declaration of the main function to run analysis
99
100private:
101 // Declaration of other main functions called in runProgram()
110
111 // Main interfaces
113
114 // Object to mark boundary entities for the assembling of domain elements
115 boost::shared_ptr<std::vector<unsigned char>> boundaryMarker;
116};
117
118HeatEquation::HeatEquation(MoFEM::Interface &m_field) : mField(m_field) {}
119
122
123 auto *simple = mField.getInterface<Simple>();
125 CHKERR simple->getOptions();
126 CHKERR simple->loadFile();
127
129}
130
133
134 auto *simple = mField.getInterface<Simple>();
135 CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 1);
136 CHKERR simple->addBoundaryField("U", H1, AINSWORTH_LEGENDRE_BASE, 1);
137
138 int order = 3;
139 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
140 CHKERR simple->setFieldOrder("U", order);
141
142 CHKERR simple->setUp();
143
145}
146
149
150 auto integration_rule = [](int o_row, int o_col, int approx_order) {
151 return 2 * approx_order;
152 };
153
154 auto *pipeline_mng = mField.getInterface<PipelineManager>();
155 CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
156 CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
157 CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(integration_rule);
158 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule);
159
161}
162
165
166 // Get surface entities form blockset, set initial values in those
167 // blocksets. To keep it simple, it is assumed that inital values are on
168 // blockset 1
170 Range inner_surface;
171 CHKERR mField.getInterface<MeshsetsManager>()->getEntitiesByDimension(
172 1, BLOCKSET, 2, inner_surface, true);
173 if (!inner_surface.empty()) {
174 Range inner_surface_verts;
175 CHKERR mField.get_moab().get_connectivity(inner_surface,
176 inner_surface_verts, false);
177 CHKERR mField.getInterface<FieldBlas>()->setField(
178 init_u, MBVERTEX, inner_surface_verts, "U");
179 }
180 }
181
183}
184
187
188 auto bc_mng = mField.getInterface<BcManager>();
189 auto *simple = mField.getInterface<Simple>();
190 CHKERR bc_mng->pushMarkDOFsOnEntities(simple->getProblemName(), "ESSENTIAL",
191 "U", 0, 0);
192
193 auto &bc_map = bc_mng->getBcMapByBlockName();
194 boundaryMarker = boost::make_shared<std::vector<char unsigned>>();
195 for (auto b : bc_map) {
196 if (std::regex_match(b.first, std::regex("(.*)ESSENTIAL(.*)"))) {
197 boundaryMarker->resize(b.second->bcMarkers.size(), 0);
198 for (int i = 0; i != b.second->bcMarkers.size(); ++i) {
199 (*boundaryMarker)[i] |= b.second->bcMarkers[i];
200 }
201 }
202 }
203
205}
206
209
210 auto add_domain_base_ops = [&](auto &pipeline) {
211 auto det_ptr = boost::make_shared<VectorDouble>();
212 auto jac_ptr = boost::make_shared<MatrixDouble>();
213 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
214 pipeline.push_back(new OpCalculateHOJac<2>(jac_ptr));
215 pipeline.push_back(new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
216 pipeline.push_back(new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
217 pipeline.push_back(new OpSetHOWeightsOnFace());
218 };
219
220 auto add_domain_lhs_ops = [&](auto &pipeline) {
221 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
222 pipeline.push_back(new OpDomainGradGrad(
223 "U", "U", [](double, double, double) -> double { return k; }));
224 auto get_c = [this](const double, const double, const double) {
225 auto pipeline_mng = mField.getInterface<PipelineManager>();
226 auto &fe_domain_lhs = pipeline_mng->getDomainLhsFE();
227 return c * fe_domain_lhs->ts_a;
228 };
229 pipeline.push_back(new OpDomainMass("U", "U", get_c));
230 pipeline.push_back(new OpUnSetBc("U"));
231 };
232
233 auto add_domain_rhs_ops = [&](auto &pipeline) {
234 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
235 auto grad_u_at_gauss_pts = boost::make_shared<MatrixDouble>();
236 auto dot_u_at_gauss_pts = boost::make_shared<VectorDouble>();
237 pipeline.push_back(new OpCalculateScalarFieldGradient<SPACE_DIM>(
238 "U", grad_u_at_gauss_pts));
239 pipeline.push_back(
240 new OpCalculateScalarFieldValuesDot("U", dot_u_at_gauss_pts));
241 pipeline.push_back(new OpDomainGradTimesVec(
242 "U", grad_u_at_gauss_pts,
243 [](double, double, double) -> double { return k; }));
244 pipeline.push_back(new OpDomainTimesScalarField(
245 "U", dot_u_at_gauss_pts,
246 [](const double, const double, const double) { return c; }));
247 auto source_term = [&](const double x, const double y, const double z) {
248 auto pipeline_mng = mField.getInterface<PipelineManager>();
249 auto &fe_domain_lhs = pipeline_mng->getDomainLhsFE();
250 const auto t = fe_domain_lhs->ts_t;
251 return 1e1 * pow(M_E, -M_PI * M_PI * t) * sin(1. * M_PI * x) *
252 sin(2. * M_PI * y);
253 };
254 pipeline.push_back(new OpDomainSource("U", source_term));
255 pipeline.push_back(new OpUnSetBc("U"));
256 };
257
258 auto add_boundary_base_ops = [&](auto &pipeline) {};
259
260 auto add_lhs_base_ops = [&](auto &pipeline) {
261 pipeline.push_back(new OpSetBc("U", false, boundaryMarker));
262 pipeline.push_back(new OpBoundaryMass(
263 "U", "U", [](const double, const double, const double) { return c; }));
264 pipeline.push_back(new OpUnSetBc("U"));
265 };
266 auto add_rhs_base_ops = [&](auto &pipeline) {
267 pipeline.push_back(new OpSetBc("U", false, boundaryMarker));
268 auto u_at_gauss_pts = boost::make_shared<VectorDouble>();
269 auto boundary_function = [&](const double x, const double y,
270 const double z) {
271 auto pipeline_mng = mField.getInterface<PipelineManager>();
272 auto &fe_rhs = pipeline_mng->getBoundaryRhsFE();
273 const auto t = fe_rhs->ts_t;
274 return 0;
275 // abs(0.1 * pow(M_E, -M_PI * M_PI * t) * sin(2. * M_PI * x) *
276 // sin(3. * M_PI * y));
277 };
278 pipeline.push_back(new OpCalculateScalarFieldValues("U", u_at_gauss_pts));
279 pipeline.push_back(new OpBoundaryTimeScalarField(
280 "U", u_at_gauss_pts,
281 [](const double, const double, const double) { return c; }));
282 pipeline.push_back(new OpBoundarySource("U", boundary_function));
283 pipeline.push_back(new OpUnSetBc("U"));
284 };
285
286 auto pipeline_mng = mField.getInterface<PipelineManager>();
287 add_domain_base_ops(pipeline_mng->getOpDomainLhsPipeline());
288 add_domain_base_ops(pipeline_mng->getOpDomainRhsPipeline());
289 add_domain_lhs_ops(pipeline_mng->getOpDomainLhsPipeline());
290 add_domain_rhs_ops(pipeline_mng->getOpDomainRhsPipeline());
291
292 add_boundary_base_ops(pipeline_mng->getOpBoundaryLhsPipeline());
293 add_boundary_base_ops(pipeline_mng->getOpBoundaryRhsPipeline());
294 add_lhs_base_ops(pipeline_mng->getOpBoundaryLhsPipeline());
295 add_rhs_base_ops(pipeline_mng->getOpBoundaryRhsPipeline());
296
298}
299
302
303 auto *simple = mField.getInterface<Simple>();
304 auto *pipeline_mng = mField.getInterface<PipelineManager>();
305
306 auto create_post_process_element = [&]() {
307 auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
308 auto det_ptr = boost::make_shared<VectorDouble>();
309 auto jac_ptr = boost::make_shared<MatrixDouble>();
310 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
311 post_proc_fe->getOpPtrVector().push_back(new OpCalculateHOJac<2>(jac_ptr));
312 post_proc_fe->getOpPtrVector().push_back(
313 new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
314 post_proc_fe->getOpPtrVector().push_back(
315 new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
316
317 auto u_ptr = boost::make_shared<VectorDouble>();
318 post_proc_fe->getOpPtrVector().push_back(
319 new OpCalculateScalarFieldValues("U", u_ptr));
320
322
323 post_proc_fe->getOpPtrVector().push_back(
324
325 new OpPPMap(post_proc_fe->getPostProcMesh(),
326 post_proc_fe->getMapGaussPts(),
327
328 {{"U", u_ptr}},
329
330 {}, {}, {}
331
332 )
333
334 );
335
336 return post_proc_fe;
337 };
338
339 auto set_time_monitor = [&](auto dm, auto solver) {
341 boost::shared_ptr<Monitor> monitor_ptr(
342 new Monitor(dm, create_post_process_element()));
343 boost::shared_ptr<ForcesAndSourcesCore> null;
344 CHKERR DMMoFEMTSSetMonitor(dm, solver, simple->getDomainFEName(),
345 monitor_ptr, null, null);
347 };
348
349 auto set_fieldsplit_preconditioner = [&](auto solver) {
351 SNES snes;
352 CHKERR TSGetSNES(solver, &snes);
353 KSP ksp;
354 CHKERR SNESGetKSP(snes, &ksp);
355 PC pc;
356 CHKERR KSPGetPC(ksp, &pc);
357 PetscBool is_pcfs = PETSC_FALSE;
358 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
359
360 if (is_pcfs == PETSC_TRUE) {
361 auto bc_mng = mField.getInterface<BcManager>();
362 auto name_prb = simple->getProblemName();
363 auto is_all_bc = bc_mng->getBlockIS(name_prb, "ESSENTIAL", "U", 0, 0);
364 int is_all_bc_size;
365 CHKERR ISGetSize(is_all_bc, &is_all_bc_size);
366 MOFEM_LOG("EXAMPLE", Sev::inform)
367 << "Field split block size " << is_all_bc_size;
368 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL,
369 is_all_bc); // boundary block
370 }
372 };
373
374 auto dm = simple->getDM();
375 auto D = smartCreateDMVector(dm);
376 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_FORWARD);
377
378 auto solver = pipeline_mng->createTSIM();
379 CHKERR TSSetSolution(solver, D);
380 CHKERR set_time_monitor(dm, solver);
381 CHKERR TSSetSolution(solver, D);
382 CHKERR TSSetFromOptions(solver);
383 CHKERR set_fieldsplit_preconditioner(solver);
384 CHKERR TSSetUp(solver);
385 CHKERR TSSolve(solver, NULL);
386
387 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
388 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
389 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
390
392}
393
396
397 // Processes to set output results are integrated in solveSystem()
398
400}
401
404
413
415}
416
417int main(int argc, char *argv[]) {
418
419 // Initialisation of MoFEM/PETSc and MOAB data structures
420 const char param_file[] = "param_file.petsc";
422
423 // Add logging channel for example
424 auto core_log = logging::core::get();
425 core_log->add_sink(
427 LogManager::setLog("EXAMPLE");
428 MOFEM_LOG_TAG("EXAMPLE", "example")
429
430 // Error handling
431 try {
432 // Register MoFEM discrete manager in PETSc
433 DMType dm_name = "DMMOFEM";
434 CHKERR DMRegister_MoFEM(dm_name);
435
436 // Create MOAB instance
437 moab::Core mb_instance; // mesh database
438 moab::Interface &moab = mb_instance; // mesh database interface
439
440 // Create MoFEM instance
441 MoFEM::Core core(moab); // finite element database
442 MoFEM::Interface &m_field = core; // finite element interface
443
444 // Run the main analysis
445 HeatEquation heat_problem(m_field);
446 CHKERR heat_problem.runProgram();
447 }
449
450 // Finish work: cleaning memory, getting statistics, etc.
452
453 return 0;
454}
std::string param_file
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
int main()
Definition: adol-c_atom.cpp:46
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, FIELD_DIM > OpDomainSource
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, FIELD_DIM > OpDomainMass
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ 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
@ BLOCKSET
Definition: definitions.h:148
#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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
ElementsAndOps< SPACE_DIM >::PostProcEle PostProcEle
auto integration_rule
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
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:301
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
check for CUBIT Id and CUBIT type
FTensor::Index< 'i', SPACE_DIM > i
constexpr double init_u
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1 > OpBoundaryTimeScalarField
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 1 > OpDomainMass
static char help[]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1 > OpDomainTimesScalarField
constexpr double k
constexpr int SPACE_DIM
[Define dimension]
constexpr double c
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 1 > OpBoundarySource
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, SPACE_DIM > OpDomainGradTimesVec
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, SPACE_DIM > OpDomainGradGrad
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 1 > OpBoundaryMass
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 1 > OpDomainSource
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, SPACE_DIM > OpDomainGradGrad
Definition: helmholtz.cpp:27
FormsIntegrators< EdgeEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 1 > OpBoundarySource
Definition: helmholtz.cpp:33
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1 > OpDomainTimesScalarField
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, SPACE_DIM > OpDomainGradTimesVec
double D
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
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)
OpCalculateScalarFieldValuesFromPetscVecImpl< PetscData::CTX_SET_X_T > OpCalculateScalarFieldValuesDot
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1 > OpBoundaryTimeScalarField
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::BiLinearForm< G >::OpMass< 1, SPACE_DIM > OpBoundaryMass
[Only used with Hencky/nonlinear material]
Definition: plastic.cpp:82
constexpr double t
plate stiffness
Definition: plate.cpp:59
static constexpr int approx_order
MoFEM::Interface & mField
MoFEMErrorCode solveSystem()
MoFEMErrorCode outputResults()
HeatEquation(MoFEM::Interface &m_field)
MoFEMErrorCode setupProblem()
MoFEMErrorCode setIntegrationRules()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
MoFEMErrorCode boundaryCondition()
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
MoFEMErrorCode runProgram()
MoFEMErrorCode initialCondition()
Simple interface for fast problem set-up.
Definition: BcManager.hpp:23
virtual moab::Interface & get_moab()=0
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.
Data on single entity (This is passed as argument to DataOperator::doWork)
structure for User Loop Methods on finite elements
Basic algebra on fields.
Definition: FieldBlas.hpp:21
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
Interface for managing meshsets containing materials and boundary conditions.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Get value at integration points for scalar field.
Post post-proc data at points from hash maps.
Set indices on entities on finite element.
Set inverse jacobian to base functions.
Modify integration weights on face to take in account higher-order geometry.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainLhsFE()
MoFEM::FaceElementForcesAndSourcesCore FaceEle
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
boost::shared_ptr< FEMethod > & getBoundaryRhsFE()
Simple interface for fast problem set-up.
Definition: Simple.hpp:26
intrusive_ptr for managing petsc objects
PetscInt ts_step
time step number
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Monitor solution.
Monitor(SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
static constexpr int saveEveryNthStep
SmartPetscObj< DM > dM
MoFEMErrorCode postProcess()
function is run at the end of loop
boost::shared_ptr< PostProcEle > postProc
MoFEMErrorCode operator()()
function is run for every finite element
MoFEMErrorCode preProcess()
function is run at the beginning of loop