v0.13.0
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 /* 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 using namespace MoFEM;
26 
27 template <int DIM> struct ElementsAndOps {};
28 
29 template <> struct ElementsAndOps<2> {
33 };
34 
35 template <> struct ElementsAndOps<3> {
39 };
40 
41 constexpr int SPACE_DIM =
42  EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
43 
48 
50  GAUSS>::OpGradTensorGrad<1, SPACE_DIM, SPACE_DIM, 1>;
56  GAUSS>::OpSource<1, SPACE_DIM>;
59 
60 constexpr bool is_quasi_static = false;
61 constexpr double rho = 1;
62 constexpr double omega = 2.4;
63 constexpr double young_modulus = 1;
64 constexpr double poisson_ratio = 0.25;
65 constexpr double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
66 constexpr double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
67 
68 #include <HenckyOps.hpp>
69 #include <OpPostProcElastic.hpp>
70 using namespace Tutorial;
71 using namespace HenckyOps;
72 
73 static double *ts_time_ptr;
74 static double *ts_aa_ptr;
75 
76 struct Example {
77 
78  Example(MoFEM::Interface &m_field) : mField(m_field) {}
79 
81 
82 private:
83  MoFEM::Interface &mField;
84 
93 
94  boost::shared_ptr<MatrixDouble> matGradPtr;
95  boost::shared_ptr<MatrixDouble> matStrainPtr;
96  boost::shared_ptr<MatrixDouble> matStressPtr;
97  boost::shared_ptr<MatrixDouble> matAccelerationPtr;
98  boost::shared_ptr<MatrixDouble> matInertiaPtr;
99  boost::shared_ptr<MatrixDouble> matDPtr;
100 
101  boost::shared_ptr<MatrixDouble> matTangentPtr;
102 };
103 
104 //! [Create common data]
107 
108  auto set_matrial_stiffness = [&]() {
115  constexpr double A =
116  (SPACE_DIM == 2) ? 2 * shear_modulus_G /
117  (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
118  : 1;
119  auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*matDPtr);
120  t_D(i, j, k, l) = 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
121  A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) *
122  t_kd(i, j) * t_kd(k, l);
124  };
125 
126  matGradPtr = boost::make_shared<MatrixDouble>();
127  matStrainPtr = boost::make_shared<MatrixDouble>();
128  matStressPtr = boost::make_shared<MatrixDouble>();
129  matAccelerationPtr = boost::make_shared<MatrixDouble>();
130  matInertiaPtr = boost::make_shared<MatrixDouble>();
131  matDPtr = boost::make_shared<MatrixDouble>();
132 
133  matTangentPtr = boost::make_shared<MatrixDouble>();
134 
135  constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
136  matDPtr->resize(size_symm * size_symm, 1);
137 
138  CHKERR set_matrial_stiffness();
139 
141 }
142 //! [Create common data]
143 
144 //! [Run problem]
147  CHKERR readMesh();
148  CHKERR setupProblem();
149  CHKERR createCommonData();
150  CHKERR boundaryCondition();
151  CHKERR assembleSystem();
152  CHKERR solveSystem();
153  CHKERR outputResults();
154  CHKERR checkResults();
156 }
157 //! [Run problem]
158 
159 //! [Read mesh]
162  auto simple = mField.getInterface<Simple>();
163  CHKERR simple->getOptions();
164  CHKERR simple->loadFile();
166 }
167 //! [Read mesh]
168 
169 //! [Set up problem]
172  Simple *simple = mField.getInterface<Simple>();
173  // Add field
174  CHKERR simple->addDomainField("U", H1, AINSWORTH_BERNSTEIN_BEZIER_BASE,
175  SPACE_DIM);
176  int order = 2;
177  CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
178  CHKERR simple->setFieldOrder("U", order);
179  CHKERR simple->setUp();
181 }
182 //! [Set up problem]
183 
184 //! [Boundary condition]
187 
188  auto simple = mField.getInterface<Simple>();
189  auto bc_mng = mField.getInterface<BcManager>();
190 
191  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_X",
192  "U", 0, 0);
193  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_Y",
194  "U", 1, 1);
195  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_Z",
196  "U", 2, 2);
197  CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "FIX_ALL",
198  "U", 0, 3);
199 
201 }
202 //! [Boundary condition]
203 
204 //! [Push operators to pipeline]
207  auto *simple = mField.getInterface<Simple>();
208  auto *pipeline_mng = mField.getInterface<PipelineManager>();
209 
210  if (SPACE_DIM == 2) {
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_mng->getOpDomainLhsPipeline().push_back(
215  new OpCalculateHOJacForFace(jac_ptr));
216  pipeline_mng->getOpDomainLhsPipeline().push_back(
217  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
218  pipeline_mng->getOpDomainLhsPipeline().push_back(
219  new OpSetInvJacH1ForFace(inv_jac_ptr));
220  pipeline_mng->getOpDomainRhsPipeline().push_back(
221  new OpCalculateHOJacForFace(jac_ptr));
222  pipeline_mng->getOpDomainRhsPipeline().push_back(
223  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
224  pipeline_mng->getOpDomainRhsPipeline().push_back(
225  new OpSetInvJacH1ForFace(inv_jac_ptr));
226  }
227 
228  // Get pointer to U_tt shift in domain element
229  auto get_rho = [this](const double, const double, const double) {
230  auto *pipeline_mng = mField.getInterface<PipelineManager>();
231  auto &fe_domain_lhs = pipeline_mng->getDomainLhsFE();
232  return rho * fe_domain_lhs->ts_aa;
233  };
234 
235  auto get_body_force = [this](const double, const double, const double) {
236  auto *pipeline_mng = mField.getInterface<PipelineManager>();
237  auto fe_domain_rhs = pipeline_mng->getDomainRhsFE();
240  t_source(i) = 0.;
241  t_source(0) = 0.1;
242  t_source(1) = 1.;
243  const auto time = fe_domain_rhs->ts_t;
244  t_source(i) *= sin(time * omega * M_PI);
245  return t_source;
246  };
247 
248  auto henky_common_data_ptr = boost::make_shared<HenckyOps::CommonData>();
249  henky_common_data_ptr->matGradPtr = matGradPtr;
250  henky_common_data_ptr->matDPtr = matDPtr;
251 
252  pipeline_mng->getOpDomainLhsPipeline().push_back(
254  matGradPtr));
255 
256  pipeline_mng->getOpDomainLhsPipeline().push_back(
257  new OpCalculateEigenVals<SPACE_DIM>("U", henky_common_data_ptr));
258  pipeline_mng->getOpDomainLhsPipeline().push_back(
259  new OpCalculateLogC<SPACE_DIM>("U", henky_common_data_ptr));
260  pipeline_mng->getOpDomainLhsPipeline().push_back(
261  new OpCalculateLogC_dC<SPACE_DIM>("U", henky_common_data_ptr));
262  pipeline_mng->getOpDomainLhsPipeline().push_back(
263  new OpCalculateHenckyStress<SPACE_DIM>("U", henky_common_data_ptr));
264  pipeline_mng->getOpDomainLhsPipeline().push_back(
265  new OpCalculatePiolaStress<SPACE_DIM>("U", henky_common_data_ptr));
266  pipeline_mng->getOpDomainLhsPipeline().push_back(
267  new OpHenckyTangent<SPACE_DIM>("U", henky_common_data_ptr));
268  pipeline_mng->getOpDomainLhsPipeline().push_back(
269  new OpK("U", "U", henky_common_data_ptr->getMatTangent()));
270 
271  if (!is_quasi_static)
272  pipeline_mng->getOpDomainLhsPipeline().push_back(
273  new OpMass("U", "U", get_rho));
274 
275  pipeline_mng->getOpDomainRhsPipeline().push_back(
276  new OpBodyForce("U", get_body_force));
277  pipeline_mng->getOpDomainRhsPipeline().push_back(
279  matGradPtr));
280 
281  pipeline_mng->getOpDomainRhsPipeline().push_back(
282  new OpCalculateEigenVals<SPACE_DIM>("U", henky_common_data_ptr));
283  pipeline_mng->getOpDomainRhsPipeline().push_back(
284  new OpCalculateLogC<SPACE_DIM>("U", henky_common_data_ptr));
285  pipeline_mng->getOpDomainRhsPipeline().push_back(
286  new OpCalculateLogC_dC<SPACE_DIM>("U", henky_common_data_ptr));
287  pipeline_mng->getOpDomainRhsPipeline().push_back(
288  new OpCalculateHenckyStress<SPACE_DIM>("U", henky_common_data_ptr));
289  pipeline_mng->getOpDomainRhsPipeline().push_back(
290  new OpCalculatePiolaStress<SPACE_DIM>("U", henky_common_data_ptr));
291 
292  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInternalForce(
293  "U", henky_common_data_ptr->getMatFirstPiolaStress()));
294  if (!is_quasi_static) {
295  pipeline_mng->getOpDomainRhsPipeline().push_back(
297  matAccelerationPtr));
298  pipeline_mng->getOpDomainRhsPipeline().push_back(
299  new OpScaleMatrix("U", rho, matAccelerationPtr, matInertiaPtr));
300  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInertiaForce(
301  "U", matInertiaPtr, [](double, double, double) { return 1.; }));
302  }
303 
304  auto integration_rule = [](int, int, int approx_order) {
305  return 2 * approx_order;
306  };
307  CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
308  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
309 
311 }
312 //! [Push operators to pipeline]
313 
314 /**
315  * @brief Monitor solution
316  *
317  * This functions is called by TS solver at the end of each step. It is used
318  * to output results to the hard drive.
319  */
320 struct Monitor : public FEMethod {
321  Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc)
322  : dM(dm), postProc(post_proc){};
325  constexpr int save_every_nth_step = 1;
326  if (ts_step % save_every_nth_step == 0) {
327  CHKERR DMoFEMLoopFiniteElements(dM, "dFE", postProc);
328  CHKERR postProc->writeFile(
329  "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
330  }
332  }
333 
334 private:
336  boost::shared_ptr<PostProcEle> postProc;
337 };
338 
339 //! [Solve]
342  auto *simple = mField.getInterface<Simple>();
343  auto *pipeline_mng = mField.getInterface<PipelineManager>();
344 
345  auto dm = simple->getDM();
347  if (is_quasi_static)
348  ts = pipeline_mng->createTSIM();
349  else
350  ts = pipeline_mng->createTSIM2();
351 
352  // Setup postprocessing
353  auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
354  post_proc_fe->generateReferenceElementMesh();
355  if (SPACE_DIM == 2) {
356  auto det_ptr = boost::make_shared<VectorDouble>();
357  auto jac_ptr = boost::make_shared<MatrixDouble>();
358  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
359  post_proc_fe->getOpPtrVector().push_back(
360  new OpCalculateHOJacForFace(jac_ptr));
361  post_proc_fe->getOpPtrVector().push_back(
362  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
363  post_proc_fe->getOpPtrVector().push_back(
364  new OpSetInvJacH1ForFace(inv_jac_ptr));
365  }
366  post_proc_fe->getOpPtrVector().push_back(
368  matGradPtr));
369  post_proc_fe->getOpPtrVector().push_back(
370  new OpSymmetrizeTensor<SPACE_DIM>("U", matGradPtr, matStrainPtr));
371  post_proc_fe->getOpPtrVector().push_back(
373  "U", matStrainPtr, matStressPtr, matDPtr));
374  post_proc_fe->getOpPtrVector().push_back(new OpPostProcElastic<SPACE_DIM>(
375  "U", post_proc_fe->postProcMesh, post_proc_fe->mapGaussPts, matStrainPtr,
376  matStressPtr));
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 TSSetMaxTime(ts, ftime);
387  CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
388 
389  auto T = smartCreateDMVector(simple->getDM());
390  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
391  SCATTER_FORWARD);
392  if (is_quasi_static) {
393  CHKERR TSSetSolution(ts, T);
394  CHKERR TSSetFromOptions(ts);
395  } else {
396  auto TT = smartVectorDuplicate(T);
397  CHKERR TS2SetSolution(ts, T, TT);
398  CHKERR TSSetFromOptions(ts);
399  }
400 
401  CHKERR TSSolve(ts, NULL);
402  CHKERR TSGetTime(ts, &ftime);
403 
404  PetscInt steps, snesfails, rejects, nonlinits, linits;
405 #if PETSC_VERSION_GE(3, 8, 0)
406  CHKERR TSGetStepNumber(ts, &steps);
407 #else
408  CHKERR TSGetTimeStepNumber(ts, &steps);
409 #endif
410  CHKERR TSGetSNESFailures(ts, &snesfails);
411  CHKERR TSGetStepRejections(ts, &rejects);
412  CHKERR TSGetSNESIterations(ts, &nonlinits);
413  CHKERR TSGetKSPIterations(ts, &linits);
414  MOFEM_LOG_C("EXAMPLE", Sev::inform,
415  "steps %D (%D rejected, %D SNES fails), ftime %g, nonlinits "
416  "%D, linits %D\n",
417  steps, rejects, snesfails, ftime, nonlinits, linits);
418 
420 }
421 //! [Solve]
422 
423 //! [Postprocess results]
426  PetscBool test_flg = PETSC_FALSE;
427  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-test", &test_flg, PETSC_NULL);
428  if (test_flg) {
429  auto *simple = mField.getInterface<Simple>();
430  auto T = smartCreateDMVector(simple->getDM());
431  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
432  SCATTER_FORWARD);
433  double nrm2;
434  CHKERR VecNorm(T, NORM_2, &nrm2);
435  MOFEM_LOG("EXAMPLE", Sev::inform) << "Regression norm " << nrm2;
436  constexpr double regression_value = 1.09572;
437  if (fabs(nrm2 - regression_value) > 1e-2)
438  SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
439  "Regression test faileed; wrong norm value.");
440  }
442 }
443 //! [Postprocess results]
444 
445 //! [Check]
449 }
450 //! [Check]
451 
452 static char help[] = "...\n\n";
453 
454 int main(int argc, char *argv[]) {
455 
456  // Initialisation of MoFEM/PETSc and MOAB data structures
457  const char param_file[] = "param_file.petsc";
458  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
459 
460  // Add logging channel for example
461  auto core_log = logging::core::get();
462  core_log->add_sink(
464  LogManager::setLog("EXAMPLE");
465  MOFEM_LOG_TAG("EXAMPLE", "example");
466 
467  try {
468 
469  //! [Register MoFEM discrete manager in PETSc]
470  DMType dm_name = "DMMOFEM";
471  CHKERR DMRegister_MoFEM(dm_name);
472  //! [Register MoFEM discrete manager in PETSc
473 
474  //! [Create MoAB]
475  moab::Core mb_instance; ///< mesh database
476  moab::Interface &moab = mb_instance; ///< mesh database interface
477  //! [Create MoAB]
478 
479  //! [Create MoFEM]
480  MoFEM::Core core(moab); ///< finite element database
481  MoFEM::Interface &m_field = core; ///< finite element database insterface
482  //! [Create MoFEM]
483 
484  //! [Example]
485  Example ex(m_field);
486  CHKERR ex.runProblem();
487  //! [Example]
488  }
489  CATCH_ERRORS;
490 
492 }
std::string param_file
#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
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 >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, SPACE_DIM > OpBodyForce
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 1 > OpInertiaForce
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:77
@ 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
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:53
#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
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
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
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
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
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:24
static constexpr int approx_order
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()
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.
Approximate field valuse for given petsc vector.
PipelineManager interface.
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.
[Class definition]