v0.13.0
cornea_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 #include <PostProcOnRefMesh.hpp>
28 
29 template <int DIM> struct ElementsAndOps {};
30 
31 template <> struct ElementsAndOps<3> {
35 };
36 
37 constexpr int SPACE_DIM = 3;
38 
43 
45  GAUSS>::OpGradTensorGrad<1, SPACE_DIM, SPACE_DIM, 1>;
51  GAUSS>::OpSource<1, SPACE_DIM>;
54 
55 constexpr bool is_quasi_static = false;
56 constexpr double rho = 1;
57 constexpr double omega = 2.4;
58 constexpr double young_modulus = 1;
59 constexpr double poisson_ratio = 0.25;
60 constexpr double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
61 constexpr double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
62 
63 #include <HenckyOps.hpp>
64 #include <OpPostProcElastic.hpp>
65 using namespace Tutorial;
66 using namespace HenckyOps;
67 
68 static double *ts_time_ptr;
69 static double *ts_aa_ptr;
70 
71 struct Example {
72 
73  Example(MoFEM::Interface &m_field) : mField(m_field) {}
74 
75  MoFEMErrorCode runProblem();
76 
77 private:
79 
80  MoFEMErrorCode readMesh();
81  MoFEMErrorCode setupProblem();
82  MoFEMErrorCode createCommonData();
83  MoFEMErrorCode boundaryCondition();
84  MoFEMErrorCode assembleSystem();
85  MoFEMErrorCode solveSystem();
86  MoFEMErrorCode outputResults();
87  MoFEMErrorCode checkResults();
88 
89  boost::shared_ptr<MatrixDouble> matGradPtr;
90  boost::shared_ptr<MatrixDouble> matStrainPtr;
91  boost::shared_ptr<MatrixDouble> matStressPtr;
92  boost::shared_ptr<MatrixDouble> matAccelerationPtr;
93  boost::shared_ptr<MatrixDouble> matInertiaPtr;
94  boost::shared_ptr<MatrixDouble> matDPtr;
95 
96  boost::shared_ptr<MatrixDouble> matTangentPtr;
97 };
98 
99 //! [Create common data]
102 
103  auto set_matrial_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  matAccelerationPtr = boost::make_shared<MatrixDouble>();
125  matInertiaPtr = boost::make_shared<MatrixDouble>();
126  matDPtr = boost::make_shared<MatrixDouble>();
127 
128  matTangentPtr = boost::make_shared<MatrixDouble>();
129 
130  constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
131  matDPtr->resize(size_symm * size_symm, 1);
132 
133  CHKERR set_matrial_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_BERNSTEIN_BEZIER_BASE,
170  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 
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(),
193  "FIX_ALL", "U", 0, 3);
194 
196 }
197 //! [Boundary condition]
198 
199 //! [Push operators to pipeline]
202  auto *simple = mField.getInterface<Simple>();
203  auto *pipeline_mng = mField.getInterface<PipelineManager>();
204 
205  if (SPACE_DIM == 2) {
206  auto det_ptr = boost::make_shared<VectorDouble>();
207  auto jac_ptr = boost::make_shared<MatrixDouble>();
208  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
209  pipeline_mng->getOpDomainLhsPipeline().push_back(
210  new OpCalculateHOJacForFace(jac_ptr));
211  pipeline_mng->getOpDomainLhsPipeline().push_back(
212  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
213  pipeline_mng->getOpDomainLhsPipeline().push_back(
214  new OpSetInvJacH1ForFace(inv_jac_ptr));
215  pipeline_mng->getOpDomainRhsPipeline().push_back(
216  new OpCalculateHOJacForFace(jac_ptr));
217  pipeline_mng->getOpDomainRhsPipeline().push_back(
218  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
219  pipeline_mng->getOpDomainRhsPipeline().push_back(
220  new OpSetInvJacH1ForFace(inv_jac_ptr));
221  }
222 
223  // Get pointer to U_tt shift in domain element
224  auto get_rho = [this](const double, const double, const double) {
225  auto *pipeline_mng = mField.getInterface<PipelineManager>();
226  auto &fe_domain_lhs = pipeline_mng->getDomainLhsFE();
227  return rho * fe_domain_lhs->ts_aa;
228  };
229 
230  auto get_body_force = [this](const double, const double, const double) {
231  auto *pipeline_mng = mField.getInterface<PipelineManager>();
232  auto fe_domain_rhs = pipeline_mng->getDomainRhsFE();
235  t_source(i) = 0.;
236  t_source(0) = 0.1;
237  t_source(1) = 1.;
238  const auto time = fe_domain_rhs->ts_t;
239  t_source(i) *= sin(time * omega * M_PI);
240  return t_source;
241  };
242 
243  auto henky_common_data_ptr = boost::make_shared<HenckyOps::CommonData>();
244  henky_common_data_ptr->matGradPtr = matGradPtr;
245  henky_common_data_ptr->matDPtr = matDPtr;
246 
247  pipeline_mng->getOpDomainLhsPipeline().push_back(
249  matGradPtr));
250 
251  pipeline_mng->getOpDomainLhsPipeline().push_back(
252  new OpCalculateEigenVals<SPACE_DIM>("U", henky_common_data_ptr));
253  pipeline_mng->getOpDomainLhsPipeline().push_back(
254  new OpCalculateLogC<SPACE_DIM>("U", henky_common_data_ptr));
255  pipeline_mng->getOpDomainLhsPipeline().push_back(
256  new OpCalculateLogC_dC<SPACE_DIM>("U", henky_common_data_ptr));
257  pipeline_mng->getOpDomainLhsPipeline().push_back(
258  new OpCalculateHenckyStress<SPACE_DIM>("U", henky_common_data_ptr));
259  pipeline_mng->getOpDomainLhsPipeline().push_back(
260  new OpCalculatePiolaStress<SPACE_DIM>("U", henky_common_data_ptr));
261  pipeline_mng->getOpDomainLhsPipeline().push_back(
262  new OpHenckyTangent<SPACE_DIM>("U", henky_common_data_ptr));
263  pipeline_mng->getOpDomainLhsPipeline().push_back(
264  new OpK("U", "U", henky_common_data_ptr->getMatTangent()));
265 
266  if (!is_quasi_static)
267  pipeline_mng->getOpDomainLhsPipeline().push_back(
268  new OpMass("U", "U", get_rho));
269 
270  pipeline_mng->getOpDomainRhsPipeline().push_back(
271  new OpBodyForce("U", get_body_force));
272  pipeline_mng->getOpDomainRhsPipeline().push_back(
274  matGradPtr));
275 
276  pipeline_mng->getOpDomainRhsPipeline().push_back(
277  new OpCalculateEigenVals<SPACE_DIM>("U", henky_common_data_ptr));
278  pipeline_mng->getOpDomainRhsPipeline().push_back(
279  new OpCalculateLogC<SPACE_DIM>("U", henky_common_data_ptr));
280  pipeline_mng->getOpDomainRhsPipeline().push_back(
281  new OpCalculateLogC_dC<SPACE_DIM>("U", henky_common_data_ptr));
282  pipeline_mng->getOpDomainRhsPipeline().push_back(
283  new OpCalculateHenckyStress<SPACE_DIM>("U", henky_common_data_ptr));
284  pipeline_mng->getOpDomainRhsPipeline().push_back(
285  new OpCalculatePiolaStress<SPACE_DIM>("U", henky_common_data_ptr));
286 
287  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInternalForce(
288  "U", henky_common_data_ptr->getMatFirstPiolaStress()));
289  if (!is_quasi_static) {
290  pipeline_mng->getOpDomainRhsPipeline().push_back(
292  matAccelerationPtr));
293  pipeline_mng->getOpDomainRhsPipeline().push_back(
294  new OpScaleMatrix("U", rho, matAccelerationPtr, matInertiaPtr));
295  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInertiaForce(
296  "U", matInertiaPtr, [](double, double, double) { return 1.; }));
297  }
298 
299  auto integration_rule = [](int, int, int approx_order) {
300  return 2 * approx_order;
301  };
302  CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
303  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
304 
306 }
307 //! [Push operators to pipeline]
308 
309 /**
310  * @brief Monitor solution
311  *
312  * This functions is called by TS solver at the end of each step. It is used
313  * to output results to the hard drive.
314  */
315 struct Monitor : public FEMethod {
316  Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc)
317  : dM(dm), postProc(post_proc){};
320  constexpr int save_every_nth_step = 1;
321  if (ts_step % save_every_nth_step == 0) {
322  CHKERR DMoFEMLoopFiniteElements(dM, "dFE", postProc);
323  CHKERR postProc->writeFile(
324  "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
325  }
327  }
328 
329 private:
331  boost::shared_ptr<PostProcEle> postProc;
332 };
333 
334 //! [Solve]
337  auto *simple = mField.getInterface<Simple>();
338  auto *pipeline_mng = mField.getInterface<PipelineManager>();
339 
340  auto dm = simple->getDM();
342  if (is_quasi_static)
343  ts = pipeline_mng->createTSIM();
344  else
345  ts = pipeline_mng->createTSIM2();
346 
347  // Setup postprocessing
348  auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
349  post_proc_fe->generateReferenceElementMesh();
350  if (SPACE_DIM == 2) {
351  auto det_ptr = boost::make_shared<VectorDouble>();
352  auto jac_ptr = boost::make_shared<MatrixDouble>();
353  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
354  post_proc_fe->getOpPtrVector().push_back(
355  new OpCalculateHOJacForFace(jac_ptr));
356  post_proc_fe->getOpPtrVector().push_back(
357  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
358  post_proc_fe->getOpPtrVector().push_back(
359  new OpSetInvJacH1ForFace(inv_jac_ptr));
360  }
361  post_proc_fe->getOpPtrVector().push_back(
363  matGradPtr));
364  post_proc_fe->getOpPtrVector().push_back(
365  new OpSymmetrizeTensor<SPACE_DIM>("U", matGradPtr, matStrainPtr));
366  post_proc_fe->getOpPtrVector().push_back(
368  "U", matStrainPtr, matStressPtr, matDPtr));
369  post_proc_fe->getOpPtrVector().push_back(new OpPostProcElastic<SPACE_DIM>(
370  "U", post_proc_fe->postProcMesh, post_proc_fe->mapGaussPts, matStrainPtr,
371  matStressPtr));
372  post_proc_fe->addFieldValuesPostProc("U");
373 
374  // Add monitor to time solver
375  boost::shared_ptr<FEMethod> null_fe;
376  auto monitor_ptr = boost::make_shared<Monitor>(dm, post_proc_fe);
377  CHKERR DMMoFEMTSSetMonitor(dm, ts, simple->getDomainFEName(), null_fe,
378  null_fe, monitor_ptr);
379 
380  double ftime = 1;
381  // CHKERR TSSetMaxTime(ts, ftime);
382  CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
383 
384  auto T = smartCreateDMVector(simple->getDM());
385  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
386  SCATTER_FORWARD);
387  if (is_quasi_static) {
388  CHKERR TSSetSolution(ts, T);
389  CHKERR TSSetFromOptions(ts);
390  } else {
391  auto TT = smartVectorDuplicate(T);
392  CHKERR TS2SetSolution(ts, T, TT);
393  CHKERR TSSetFromOptions(ts);
394  }
395 
396  CHKERR TSSolve(ts, NULL);
397  CHKERR TSGetTime(ts, &ftime);
398 
399  PetscInt steps, snesfails, rejects, nonlinits, linits;
400 #if PETSC_VERSION_GE(3, 8, 0)
401  CHKERR TSGetStepNumber(ts, &steps);
402 #else
403  CHKERR TSGetTimeStepNumber(ts, &steps);
404 #endif
405  CHKERR TSGetSNESFailures(ts, &snesfails);
406  CHKERR TSGetStepRejections(ts, &rejects);
407  CHKERR TSGetSNESIterations(ts, &nonlinits);
408  CHKERR TSGetKSPIterations(ts, &linits);
409  MOFEM_LOG_C("EXAMPLE", Sev::inform,
410  "steps %D (%D rejected, %D SNES fails), ftime %g, nonlinits "
411  "%D, linits %D\n",
412  steps, rejects, snesfails, ftime, nonlinits, linits);
413 
415 }
416 //! [Solve]
417 
418 //! [Postprocess results]
421  PetscBool test_flg = PETSC_FALSE;
422  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-test", &test_flg, PETSC_NULL);
423  if (test_flg) {
424  auto *simple = mField.getInterface<Simple>();
425  auto T = smartCreateDMVector(simple->getDM());
426  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
427  SCATTER_FORWARD);
428  double nrm2;
429  CHKERR VecNorm(T, NORM_2, &nrm2);
430  MOFEM_LOG("EXAMPLE", Sev::inform) << "Regression norm " << nrm2;
431  constexpr double regression_value = 1.09572;
432  if (fabs(nrm2 - regression_value) > 1e-2)
433  SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
434  "Regression test faileed; wrong norm value.");
435  }
437 }
438 //! [Postprocess results]
439 
440 //! [Check]
444 }
445 //! [Check]
446 
447 static char help[] = "...\n\n";
448 
449 int main(int argc, char *argv[]) {
450 
451  // Initialisation of MoFEM/PETSc and MOAB data structures
452  const char param_file[] = "param_file.petsc";
453  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
454 
455  // Add logging channel for example
456  auto core_log = logging::core::get();
457  core_log->add_sink(
459  LogManager::setLog("EXAMPLE");
460  MOFEM_LOG_TAG("EXAMPLE", "example");
461 
462  try {
463 
464  //! [Register MoFEM discrete manager in PETSc]
465  DMType dm_name = "DMMOFEM";
466  CHKERR DMRegister_MoFEM(dm_name);
467  //! [Register MoFEM discrete manager in PETSc
468 
469  //! [Create MoAB]
470  moab::Core mb_instance; ///< mesh database
471  moab::Interface &moab = mb_instance; ///< mesh database interface
472  //! [Create MoAB]
473 
474  //! [Create MoFEM]
475  MoFEM::Core core(moab); ///< finite element database
476  MoFEM::Interface &m_field = core; ///< finite element database insterface
477  //! [Create MoFEM]
478 
479  //! [Example]
480  Example ex(m_field);
481  CHKERR ex.runProblem();
482  //! [Example]
483  }
484  CATCH_ERRORS;
485 
487 }
std::string param_file
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
Post-process fields on refined mesh.
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.
int main(int argc, char *argv[])
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
DataForcesAndSourcesCore::EntData EntData
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
#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.
double A
int save_every_nth_step
static constexpr int approx_order
[Example]
MoFEMErrorCode boundaryCondition()
[Set up problem]
MoFEMErrorCode assembleSystem()
[Boundary condition]
MoFEMErrorCode readMesh()
[Run problem]
boost::shared_ptr< MatrixDouble > matTangentPtr
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()
[Postprocess results]
MoFEMErrorCode solveSystem()
[Solve]
MoFEMErrorCode createCommonData()
[Create common data]
Example(MoFEM::Interface &m_field)
MoFEMErrorCode runProblem()
[Create common data]
MoFEM::Interface & mField
MoFEMErrorCode setupProblem()
[Read mesh]
MoFEMErrorCode outputResults()
[Solve]
boost::shared_ptr< MatrixDouble > matDPtr
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)
SmartPetscObj< DM > dM
MoFEMErrorCode postProcess()
function is run at the end of loop
boost::shared_ptr< PostProcEle > postProc
Post processing.
[Class definition]