v0.13.0
eigen_elastic.cpp
Go to the documentation of this file.
1 /**
2  * \file eigen_elastic.cpp
3  * \example eigen_elastic.cpp
4  *
5  * Calculate natural frequencies in 2d and 3d problems.
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 #undef EPS
25 #include <slepceps.h>
26 
27 using namespace MoFEM;
28 
29 template <int DIM> struct ElementsAndOps {};
30 
31 template <> struct ElementsAndOps<2> {
35 };
36 
37 template <> struct ElementsAndOps<3> {
41 };
42 
43 constexpr int SPACE_DIM =
44  EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
45 
50 
52  GAUSS>::OpGradSymTensorGrad<1, SPACE_DIM, SPACE_DIM, 0>;
55 
56 double rho = 7829e-12;
57 double young_modulus = 207e3;
58 double poisson_ratio = 0.33;
59 
60 double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
62 
63 int order = 1;
64 
65 #include <OpPostProcElastic.hpp>
66 using namespace Tutorial;
67 
68 struct Example {
69 
70  Example(MoFEM::Interface &m_field) : mField(m_field) {}
71 
73 
74 private:
75  MoFEM::Interface &mField;
76 
85 
86  boost::shared_ptr<MatrixDouble> matGradPtr;
87  boost::shared_ptr<MatrixDouble> matStrainPtr;
88  boost::shared_ptr<MatrixDouble> matStressPtr;
89  boost::shared_ptr<MatrixDouble> matDPtr;
90 
94 
95  std::array<SmartPetscObj<Vec>, 6> rigidBodyMotion;
96 };
97 
98 //! [Create common data]
101 
102  CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-rho", &rho, PETSC_NULL);
103  CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-young_modulus", &young_modulus,
104  PETSC_NULL);
105  CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-poisson_ratio", &poisson_ratio,
106  PETSC_NULL);
107 
108  auto set_matrial_stiffens = [&]() {
115  auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*matDPtr);
116 
117  const double A = (SPACE_DIM == 2)
118  ? 2 * shear_modulus_G /
119  (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
120  : 1;
121  t_D(i, j, k, l) = 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
122  A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) *
123  t_kd(i, j) * t_kd(k, l);
124 
126  };
127 
128  matGradPtr = boost::make_shared<MatrixDouble>();
129  matStrainPtr = boost::make_shared<MatrixDouble>();
130  matStressPtr = boost::make_shared<MatrixDouble>();
131  matDPtr = boost::make_shared<MatrixDouble>();
132 
133  constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
134  matDPtr->resize(size_symm * size_symm, 1);
135 
136  CHKERR set_matrial_stiffens();
137 
139 }
140 //! [Create common data]
141 
142 //! [Run problem]
145  CHKERR readMesh();
146  CHKERR setupProblem();
147  CHKERR createCommonData();
148  CHKERR boundaryCondition();
149  CHKERR assembleSystem();
150  CHKERR solveSystem();
151  CHKERR outputResults();
152  CHKERR checkResults();
154 }
155 //! [Run problem]
156 
157 //! [Read mesh]
160  auto simple = mField.getInterface<Simple>();
161 
162  MOFEM_LOG("EXAMPLE", Sev::inform)
163  << "Read mesh for problem in " << EXECUTABLE_DIMENSION;
164  CHKERR simple->getOptions();
165  CHKERR simple->loadFile();
167 }
168 //! [Read mesh]
169 
170 //! [Set up problem]
172  auto *simple = mField.getInterface<Simple>();
174  // Add field
175  CHKERR simple->addDomainField("U", H1, AINSWORTH_BERNSTEIN_BEZIER_BASE,
176  SPACE_DIM);
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]
186  auto *simple = mField.getInterface<Simple>();
188 
189  rigidBodyMotion[0] = smartCreateDMVector(simple->getDM());
190  for (int n = 1; n != 6; ++n)
191  rigidBodyMotion[n] = smartVectorDuplicate(rigidBodyMotion[0]);
192 
193  // Create space of vectors or rigid motion
194  auto problem_ptr = mField.get_problem(simple->getProblemName());
195  auto dofs = problem_ptr->getNumeredRowDofsPtr();
196 
197  // Get all vertices
198  const auto bit_number = mField.get_field_bit_number("U");
199  auto lo_uid =
200  DofEntity::getUniqueIdCalculate(0, get_id_for_min_type<MBVERTEX>());
201  auto hi_uid =
202  DofEntity::getUniqueIdCalculate(2, get_id_for_max_type<MBVERTEX>());
203 
204  auto hi = dofs->upper_bound(lo_uid);
205  std::array<double, 3> coords;
206 
207  for (auto lo = dofs->lower_bound(lo_uid); lo != hi; ++lo) {
208 
209  if ((*lo)->getPart() == mField.get_comm_rank()) {
210 
211  auto ent = (*lo)->getEnt();
212  CHKERR mField.get_moab().get_coords(&ent, 1, coords.data());
213 
214  if ((*lo)->getDofCoeffIdx() == 0) {
215  CHKERR VecSetValue(rigidBodyMotion[0], (*lo)->getPetscGlobalDofIdx(), 1,
216  INSERT_VALUES);
217  CHKERR VecSetValue(rigidBodyMotion[3], (*lo)->getPetscGlobalDofIdx(),
218  -coords[1], INSERT_VALUES);
219  if (SPACE_DIM == 3)
220  CHKERR VecSetValue(rigidBodyMotion[4], (*lo)->getPetscGlobalDofIdx(),
221  -coords[2], INSERT_VALUES);
222 
223  } else if ((*lo)->getDofCoeffIdx() == 1) {
224  CHKERR VecSetValue(rigidBodyMotion[1], (*lo)->getPetscGlobalDofIdx(), 1,
225  INSERT_VALUES);
226  CHKERR VecSetValue(rigidBodyMotion[3], (*lo)->getPetscGlobalDofIdx(),
227  coords[0], INSERT_VALUES);
228  if (SPACE_DIM == 3)
229  CHKERR VecSetValue(rigidBodyMotion[5], (*lo)->getPetscGlobalDofIdx(),
230  -coords[2], INSERT_VALUES);
231 
232  } else if ((*lo)->getDofCoeffIdx() == 2) {
233  if (SPACE_DIM == 3) {
234  CHKERR VecSetValue(rigidBodyMotion[2], (*lo)->getPetscGlobalDofIdx(),
235  1, INSERT_VALUES);
236  CHKERR VecSetValue(rigidBodyMotion[4], (*lo)->getPetscGlobalDofIdx(),
237  coords[0], INSERT_VALUES);
238  CHKERR VecSetValue(rigidBodyMotion[5], (*lo)->getPetscGlobalDofIdx(),
239  coords[1], INSERT_VALUES);
240  }
241  }
242  }
243  }
244 
245  for (int n = 0; n != rigidBodyMotion.size(); ++n) {
246  CHKERR VecAssemblyBegin(rigidBodyMotion[n]);
247  CHKERR VecAssemblyEnd(rigidBodyMotion[n]);
248  CHKERR VecGhostUpdateBegin(rigidBodyMotion[n], INSERT_VALUES,
249  SCATTER_FORWARD);
250  CHKERR VecGhostUpdateEnd(rigidBodyMotion[n], INSERT_VALUES,
251  SCATTER_FORWARD);
252  }
253 
255 }
256 //! [Boundary condition]
257 
258 //! [Push operators to pipeline]
261  auto *simple = mField.getInterface<Simple>();
262  auto *pipeline_mng = mField.getInterface<PipelineManager>();
263 
264  auto dm = simple->getDM();
266  M = smartMatDuplicate(K, MAT_SHARE_NONZERO_PATTERN);
267 
268  auto calculate_stiffness_matrix = [&]() {
270  pipeline_mng->getDomainLhsFE().reset();
271 
272  if (SPACE_DIM == 2) {
273  auto det_ptr = boost::make_shared<VectorDouble>();
274  auto jac_ptr = boost::make_shared<MatrixDouble>();
275  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
276  pipeline_mng->getOpDomainLhsPipeline().push_back(
277  new OpCalculateHOJacForFace(jac_ptr));
278  pipeline_mng->getOpDomainLhsPipeline().push_back(
279  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
280  pipeline_mng->getOpDomainLhsPipeline().push_back(
281  new OpSetInvJacH1ForFace(inv_jac_ptr));
282  }
283 
284  pipeline_mng->getOpDomainLhsPipeline().push_back(
285  new OpK("U", "U", matDPtr));
286  auto integration_rule = [](int, int, int approx_order) {
287  return 2 * (approx_order - 1);
288  };
289 
290  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
291  pipeline_mng->getDomainLhsFE()->B = K;
292  CHKERR MatZeroEntries(K);
293  CHKERR pipeline_mng->loopFiniteElements();
294  CHKERR MatAssemblyBegin(K, MAT_FINAL_ASSEMBLY);
295  CHKERR MatAssemblyEnd(K, MAT_FINAL_ASSEMBLY);
297  };
298 
299  auto calculate_mass_matrix = [&]() {
301  pipeline_mng->getDomainLhsFE().reset();
302  auto get_rho = [](const double, const double, const double) { return rho; };
303  pipeline_mng->getOpDomainLhsPipeline().push_back(
304  new OpMass("U", "U", get_rho));
305  auto integration_rule = [](int, int, int approx_order) {
306  return 2 * approx_order;
307  };
308  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
309  CHKERR MatZeroEntries(M);
310  pipeline_mng->getDomainLhsFE()->B = M;
311  CHKERR pipeline_mng->loopFiniteElements();
312  CHKERR MatAssemblyBegin(M, MAT_FINAL_ASSEMBLY);
313  CHKERR MatAssemblyEnd(M, MAT_FINAL_ASSEMBLY);
315  };
316 
317  CHKERR calculate_stiffness_matrix();
318  CHKERR calculate_mass_matrix();
319 
321 }
322 //! [Push operators to pipeline]
323 
324 //! [Solve]
326  auto simple = mField.getInterface<Simple>();
327  auto is_manager = mField.getInterface<ISManager>();
329 
330  auto create_eps = [](MPI_Comm comm) {
331  EPS eps;
332  CHKERR EPSCreate(comm, &eps);
333  return SmartPetscObj<EPS>(eps);
334  };
335 
336  auto deflate_vectors = [&]() {
338  // Deflate vectors
339  std::array<Vec, 6> deflate_vectors;
340  for (int n = 0; n != 6; ++n) {
341  deflate_vectors[n] = rigidBodyMotion[n];
342  }
343  CHKERR EPSSetDeflationSpace(ePS, 6, &deflate_vectors[0]);
345  };
346 
347  auto print_info = [&]() {
349  ST st;
350  EPSType type;
351  PetscReal tol;
352  PetscInt nev, maxit, its;
353  // Optional: Get some information from the solver and display it
354  CHKERR EPSGetIterationNumber(ePS, &its);
355  MOFEM_LOG_C("EXAMPLE", Sev::inform,
356  " Number of iterations of the method: %d", its);
357  CHKERR EPSGetST(ePS, &st);
358  CHKERR EPSGetType(ePS, &type);
359  MOFEM_LOG_C("EXAMPLE", Sev::inform, " Solution method: %s", type);
360  CHKERR EPSGetDimensions(ePS, &nev, NULL, NULL);
361  MOFEM_LOG_C("EXAMPLE", Sev::inform, " Number of requested eigenvalues: %d",
362  nev);
363  CHKERR EPSGetTolerances(ePS, &tol, &maxit);
364  MOFEM_LOG_C("EXAMPLE", Sev::inform,
365  " Stopping condition: tol=%.4g, maxit=%d", (double)tol, maxit);
366 
367  PetscScalar eigr, eigi;
368  for (int nn = 0; nn < nev; nn++) {
369  CHKERR EPSGetEigenpair(ePS, nn, &eigr, &eigi, PETSC_NULL, PETSC_NULL);
370  MOFEM_LOG_C("EXAMPLE", Sev::inform,
371  " ncov = %d eigr = %.4g eigi = %.4g (inv eigr = %.4g)", nn,
372  eigr, eigi, 1. / eigr);
373  }
374 
376  };
377 
378  auto setup_eps = [&]() {
380  CHKERR EPSSetProblemType(ePS, EPS_GHEP);
381  CHKERR EPSSetWhichEigenpairs(ePS, EPS_SMALLEST_MAGNITUDE);
382  CHKERR EPSSetFromOptions(ePS);
384  };
385 
386  // Create eigensolver context
387  ePS = create_eps(mField.get_comm());
388  CHKERR EPSSetOperators(ePS, K, M);
389 
390  // Setup eps
391  CHKERR setup_eps();
392 
393  // Deflate vectors
394  CHKERR deflate_vectors();
395 
396  // Solve problem
397  CHKERR EPSSolve(ePS);
398 
399  // Print info
400  CHKERR print_info();
401 
403 }
404 //! [Solve]
405 
406 //! [Postprocess results]
409  auto *pipeline_mng = mField.getInterface<PipelineManager>();
410  auto *simple = mField.getInterface<Simple>();
411 
412  pipeline_mng->getDomainLhsFE().reset();
413  auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
414  post_proc_fe->generateReferenceElementMesh();
415 
416  if (SPACE_DIM == 2) {
417  auto det_ptr = boost::make_shared<VectorDouble>();
418  auto jac_ptr = boost::make_shared<MatrixDouble>();
419  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
420  post_proc_fe->getOpPtrVector().push_back(
421  new OpCalculateHOJacForFace(jac_ptr));
422  post_proc_fe->getOpPtrVector().push_back(
423  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
424  post_proc_fe->getOpPtrVector().push_back(
425  new OpSetInvJacH1ForFace(inv_jac_ptr));
426  }
427 
428  post_proc_fe->getOpPtrVector().push_back(
430  matGradPtr));
431  post_proc_fe->getOpPtrVector().push_back(
432  new OpSymmetrizeTensor<SPACE_DIM>("U", matGradPtr, matStrainPtr));
433  post_proc_fe->getOpPtrVector().push_back(
435  "U", matStrainPtr, matStressPtr, matDPtr));
436  post_proc_fe->getOpPtrVector().push_back(new OpPostProcElastic<SPACE_DIM>(
437  "U", post_proc_fe->postProcMesh, post_proc_fe->mapGaussPts, matStrainPtr,
438  matStressPtr));
439 
440  post_proc_fe->addFieldValuesPostProc("U");
441  pipeline_mng->getDomainRhsFE() = post_proc_fe;
442 
443  auto dm = simple->getDM();
444  auto D = smartCreateDMVector(dm);
445 
446  PetscInt nev;
447  CHKERR EPSGetDimensions(ePS, &nev, NULL, NULL);
448  PetscScalar eigr, eigi, nrm2r;
449  for (int nn = 0; nn < nev; nn++) {
450  CHKERR EPSGetEigenpair(ePS, nn, &eigr, &eigi, D, PETSC_NULL);
451  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
452  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
453  CHKERR VecNorm(D, NORM_2, &nrm2r);
454  MOFEM_LOG_C("EXAMPLE", Sev::inform,
455  " ncov = %d omega2 = %.8g omega = %.8g frequency = %.8g", nn,
456  eigr, std::sqrt(std::abs(eigr)), std::sqrt(std::abs(eigr)) / (2 * M_PI));
457  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
458  CHKERR pipeline_mng->loopFiniteElements();
459  post_proc_fe->writeFile("out_eig_" + boost::lexical_cast<std::string>(nn) +
460  ".h5m");
461  }
462 
464 }
465 //! [Postprocess results]
466 
467 //! [Check]
470  PetscBool test_flg = PETSC_FALSE;
471  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-test", &test_flg, PETSC_NULL);
472  if (test_flg) {
473  PetscScalar eigr, eigi;
474  CHKERR EPSGetEigenpair(ePS, 0, &eigr, &eigi, PETSC_NULL, PETSC_NULL);
475  constexpr double regression_value = 12579658;
476  if (fabs(eigr - regression_value) > 1)
477  SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
478  "Regression test faileed; wrong eigen value.");
479  }
481 }
482 //! [Check]
483 
484 static char help[] = "...\n\n";
485 
486 int main(int argc, char *argv[]) {
487 
488  // Initialisation of MoFEM/PETSc and MOAB data structures
489  const char param_file[] = "param_file.petsc";
490  SlepcInitialize(&argc, &argv, param_file, help);
491  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
492 
493  // Add logging channel for example
494  auto core_log = logging::core::get();
495  core_log->add_sink(
497  LogManager::setLog("EXAMPLE");
498  MOFEM_LOG_TAG("EXAMPLE", "example");
499 
500  try {
501 
502  //! [Register MoFEM discrete manager in PETSc]
503  DMType dm_name = "DMMOFEM";
504  CHKERR DMRegister_MoFEM(dm_name);
505  //! [Register MoFEM discrete manager in PETSc
506 
507  //! [Create MoAB]
508  moab::Core mb_instance; ///< mesh database
509  moab::Interface &moab = mb_instance; ///< mesh database interface
510  //! [Create MoAB]
511 
512  //! [Create MoFEM]
513  MoFEM::Core core(moab); ///< finite element database
514  MoFEM::Interface &m_field = core; ///< finite element database insterface
515  //! [Create MoFEM]
516 
517  //! [Example]
518  Example ex(m_field);
519  CHKERR ex.runProblem();
520  //! [Example]
521  }
522  CATCH_ERRORS;
523 
524  SlepcFinalize();
526 }
static Index< 'M', 3 > M
static Index< 'n', 3 > n
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
static const double eps
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Kronecker Delta class symmetric.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
#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
double young_modulus
int main(int argc, char *argv[])
EntitiesFieldData::EntData EntData
static char help[]
[Check]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
double rho
constexpr int SPACE_DIM
double poisson_ratio
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradSymTensorGrad< 1, SPACE_DIM, SPACE_DIM, 0 > OpK
int order
double bulk_modulus_K
double shear_modulus_G
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 DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1135
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
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
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
double tol
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
CoreTmp< 0 > Core
Definition: Core.hpp:1096
SmartPetscObj< Mat > smartMatDuplicate(Mat &mat, MatDuplicateOption op)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
const double D
diffusivity
double A
#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()
std::array< SmartPetscObj< Vec >, 6 > rigidBodyMotion
Example(MoFEM::Interface &m_field)
SmartPetscObj< Mat > M
MoFEMErrorCode runProblem()
SmartPetscObj< EPS > ePS
SmartPetscObj< Mat > K
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
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.
static UId getUniqueIdCalculate(const DofIdx dof, UId ent_uid)
Data on single entity (This is passed as argument to DataOperator::doWork)
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:33
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.
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
Postprocess on face.
Post processing.
[Class definition]