v0.13.1
child_and_parent.cpp
Go to the documentation of this file.
1/**
2 * \example child_and_parent.cpp
3 *
4 * Testing projection child and parent.
5 *
6 */
7
8/* This file is part of MoFEM.
9 * MoFEM is free software: you can redistribute it and/or modify it under
10 * the terms of the GNU Lesser General Public License as published by the
11 * Free Software Foundation, either version 3 of the License, or (at your
12 * option) any later version.
13 *
14 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
15 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
16 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17 * License for more details.
18 *
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
21
22#include <MoFEM.hpp>
23
24using namespace MoFEM;
25
26static char help[] = "...\n\n";
27
28constexpr char FIELD_NAME[] = "U";
29constexpr int FIELD_DIM = 1;
30constexpr int SPACE_DIM = 2;
31
32template <int DIM> struct ElementsAndOps {};
33
34template <> struct ElementsAndOps<2> {
41};
42
43template <> struct ElementsAndOps<3> {
46};
47
52
56
57template <int FIELD_DIM> struct ApproxFieldFunction;
58
59template <> struct ApproxFieldFunction<1> {
60 double operator()(const double x, const double y, const double z) {
61 return x * x + y * y + x * y + pow(x, 3) + pow(y, 3) + pow(x, 4) +
62 pow(y, 4);
63 }
64};
65
70
71struct AtomTest {
72
73 AtomTest(MoFEM::Interface &m_field) : mField(m_field) {}
74
76
77private:
80
82
88 checkResults(boost::function<bool(FEMethod *fe_method_ptr)> test_bit);
90 struct CommonData {
91 boost::shared_ptr<VectorDouble> approxVals;
94 };
95
96 template <int FIELD_DIM> struct OpError;
97};
98
101
102template <> struct AtomTest::OpError<1> : public DomainEleOp {
103 boost::shared_ptr<CommonData> commonDataPtr;
104 OpError(boost::shared_ptr<CommonData> &common_data_ptr)
105 : DomainEleOp(FIELD_NAME, OPROW), commonDataPtr(common_data_ptr) {}
108
109 if (const size_t nb_dofs = data.getIndices().size()) {
110
111 const int nb_integration_pts = getGaussPts().size2();
112 auto t_w = getFTensor0IntegrationWeight();
113 auto t_val = getFTensor0FromVec(*(commonDataPtr->approxVals));
114 auto t_coords = getFTensor1CoordsAtGaussPts();
115
116 VectorDouble nf(nb_dofs, false);
117 nf.clear();
118
120 const double volume = getMeasure();
121
122 auto t_row_base = data.getFTensor0N();
123 double error = 0;
124 for (int gg = 0; gg != nb_integration_pts; ++gg) {
125
126 const double alpha = t_w * volume;
127 double diff = t_val - AtomTest::approxFunction(t_coords(0), t_coords(1),
128 t_coords(2));
129 error += alpha * pow(diff, 2);
130
131 for (size_t r = 0; r != nb_dofs; ++r) {
132 nf[r] += alpha * t_row_base * diff;
133 ++t_row_base;
134 }
135
136 ++t_w;
137 ++t_val;
138 ++t_coords;
139 }
140
141 const int index = 0;
142 CHKERR VecSetValue(commonDataPtr->L2Vec, index, error, ADD_VALUES);
143 CHKERR VecSetValues(commonDataPtr->resVec, data, &nf[0], ADD_VALUES);
144 }
145
147 }
148};
149
150template <typename ELE_OP, typename PARENT_ELE>
151struct OpCheckGaussCoords : public ELE_OP {
153
157
158 MatrixDouble parent_coords;
159
160 PARENT_ELE parent_fe(this->getPtrFE()->mField);
161 auto op = new ELE_OP(NOSPACE, ELE_OP::OPSPACE);
162 op->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
166
167 MOFEM_LOG("SELF", Sev::noisy)
168 << "parent_coords in op "
169 << static_cast<ELE_OP *>(op_ptr)->getCoordsAtGaussPts();
170
171 parent_coords = static_cast<ELE_OP *>(op_ptr)->getCoordsAtGaussPts();
173 };
174 parent_fe.getOpPtrVector().push_back(op);
175
176 MOFEM_LOG("SELF", Sev::noisy) << "fe name " << this->getFEName();
177 CHKERR this->loopParent(this->getFEName(), &parent_fe);
178 MOFEM_LOG("SELF", Sev::noisy) << "parent_coords " << parent_coords;
179
180 MatrixDouble child_coords = this->getCoordsAtGaussPts();
181 MOFEM_LOG("SELF", Sev::noisy) << "child_coords " << child_coords;
182
183 child_coords -= parent_coords;
184
185 MOFEM_LOG("SELF", Sev::noisy) << "Corrds diffs" << child_coords;
186
187 double n = 0;
188 for (auto d : child_coords.data())
189 n += d * d;
190
191 if (sqrt(n) > 1e-12)
192 SETERRQ1(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
193 "Parent and child global coords at integration points are "
194 "diffrent norm = %3.2e",
195 sqrt(n));
196
198 }
199};
200
201//! [Run programme]
208
209 auto test_bit_child = [](FEMethod *fe_ptr) {
210 const auto &bit = fe_ptr->numeredEntFiniteElementPtr->getBitRefLevel();
211 MOFEM_LOG("SELF", Sev::noisy) << bit << " " << bit.test(0);
212 return bit.test(1);
213 };
214
218}
219//! [Run programme]
220
221//! [Read mesh]
224
228
229 MOFEM_LOG("WORLD", Sev::verbose) << "Dim " << simpleInterface->getDim();
230
231 auto bit_level0 = simpleInterface->getBitRefLevel();
232
233 auto &moab = mField.get_moab();
234
235 auto refine_mesh = [&](auto bit_level1) {
237
238 auto refine = mField.getInterface<MeshRefinement>();
239
240 auto meshset_level0_ptr = get_temp_meshset_ptr(moab);
241 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByRefLevel(
242 bit_level0, BitRefLevel().set(), *meshset_level0_ptr);
243
244 // random mesh refinement
245 auto meshset_ref_edges_ptr = get_temp_meshset_ptr(moab);
246 Range edges_to_refine;
247 CHKERR moab.get_entities_by_type(*meshset_level0_ptr, MBEDGE,
248 edges_to_refine);
249 int ii = 0;
250 for (Range::iterator eit = edges_to_refine.begin();
251 eit != edges_to_refine.end(); eit++, ii++) {
252 int numb = ii % 2;
253 if (numb == 0) {
254 CHKERR moab.add_entities(*meshset_ref_edges_ptr, &*eit, 1);
255 }
256 }
257 CHKERR refine->addVerticesInTheMiddleOfEdges(*meshset_ref_edges_ptr,
258 bit_level1, false, VERBOSE);
259 if (simpleInterface->getDim() == 3) {
260 CHKERR refine->refineTets(*meshset_level0_ptr, bit_level1, VERBOSE);
261 } else if (simpleInterface->getDim() == 2) {
262 CHKERR refine->refineTris(*meshset_level0_ptr, bit_level1, VERBOSE);
263 } else {
264 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
265 "Dimension not handled by test");
266 }
267
269 };
270
271 BitRefLevel bit_level1;
272 bit_level1.set(1);
273 CHKERR refine_mesh(bit_level1);
276
278}
279//! [Read mesh]
280
281//! [Set up problem]
284 // Add field
289 constexpr int order = 4;
292
293 CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
295 BitRefLevel().set(0));
296
298}
299//! [Set up problem]
300
301boost::shared_ptr<DomainEle> domainChildLhs, domainChildRhs;
302
303//! [Push operators to pipeline]
307
308 auto rule = [](int, int, int p) -> int { return 2 * p; };
309
310 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule);
311 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
312
313 auto test_bit_parent = [](FEMethod *fe_ptr) {
314 const auto &bit = fe_ptr->numeredEntFiniteElementPtr->getBitRefLevel();
315 MOFEM_LOG("SELF", Sev::noisy) << bit << " " << bit.test(0);
316 return bit.test(0);
317 };
318
319 pipeline_mng->getDomainLhsFE()->exeTestHook = test_bit_parent;
320 pipeline_mng->getDomainRhsFE()->exeTestHook = test_bit_parent;
321
322 auto beta = [](const double, const double, const double) { return 1; };
323
324 // Make aliased shared pointer, and create child element
325 domainChildLhs = boost::make_shared<DomainEle>(mField);
326 domainChildLhs->getRuleHook = rule;
327 domainChildLhs->getOpPtrVector().push_back(
329
330 domainChildRhs = boost::make_shared<DomainEle>(mField);
331 domainChildLhs->getRuleHook = rule;
332 domainChildRhs->getOpPtrVector().push_back(
334
335 auto parent_op_lhs = new DomainEleOp(NOSPACE, DomainEleOp::OPSPACE);
336 parent_op_lhs->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
339 auto domain_op = static_cast<DomainEleOp *>(op_ptr);
341
342 MOFEM_LOG("SELF", Sev::noisy) << "LHS Pipeline FE";
343
344 if (!domainChildLhs)
345 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "FE not allocated");
346
347 auto &bit =
348 domain_op->getFEMethod()->numeredEntFiniteElementPtr->getBitRefLevel();
349 if (bit == BitRefLevel().set(0)) {
350 CHKERR domain_op->loopChildren(domain_op->getFEName(),
351 domainChildLhs.get(), VERBOSE, Sev::noisy);
352 } else {
353 CHKERR domain_op->loopThis(domain_op->getFEName(), domainChildLhs.get(),
354 VERBOSE, Sev::noisy);
355 }
357 };
358
359 auto parent_op_rhs = new DomainEleOp(NOSPACE, DomainEleOp::OPSPACE);
360 parent_op_rhs->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
363 auto domain_op = static_cast<DomainEleOp *>(op_ptr);
365
366 MOFEM_LOG("SELF", Sev::noisy) << "RHS Pipeline FE";
367
368 if (!domainChildRhs)
369 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "FE not allocated");
370
371 auto &bit =
372 domain_op->getFEMethod()->numeredEntFiniteElementPtr->getBitRefLevel();
373 if (bit == BitRefLevel().set(0)) {
374 CHKERR domain_op->loopChildren(domain_op->getFEName(),
375 domainChildRhs.get(), VERBOSE, Sev::noisy);
376 } else if ((bit & BitRefLevel().set(0)).any()) {
377 CHKERR domain_op->loopThis(domain_op->getFEName(), domainChildRhs.get(),
378 VERBOSE, Sev::noisy);
379 }
381 };
382
383 pipeline_mng->getOpDomainLhsPipeline().push_back(parent_op_lhs);
384 pipeline_mng->getOpDomainRhsPipeline().push_back(parent_op_rhs);
385
387}
388//! [Push operators to pipeline]
389
390//! [Solve]
394
395 MOFEM_LOG("WORLD", Sev::inform) << "Solve problem";
396
397 auto solver = pipeline_mng->createKSP();
398 CHKERR KSPSetFromOptions(solver);
399 CHKERR KSPSetUp(solver);
400
401 auto dm = simpleInterface->getDM();
402 auto D = smartCreateDMVector(dm);
403 auto F = smartVectorDuplicate(D);
404
405 CHKERR KSPSolve(solver, F, D);
406 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
407 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
408 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
410}
411
412//! [Solve]
415
416 auto &moab = mField.get_moab();
417
418 auto bit_level0 = BitRefLevel().set(0);
419 auto bit_level1 = BitRefLevel().set(1);
420 auto bit_level2 = BitRefLevel().set(2);
421
422 auto refine_mesh = [&]() {
424
425 auto refine = mField.getInterface<MeshRefinement>();
426
427 auto meshset_level1_ptr = get_temp_meshset_ptr(moab);
428 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByDimAndRefLevel(
429 bit_level1, BitRefLevel().set(), simpleInterface->getDim(),
430 *meshset_level1_ptr);
431
432 // random mesh refinement
433 auto meshset_ref_edges_ptr = get_temp_meshset_ptr(moab);
434 Range edges_to_refine;
435 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
436 bit_level1, BitRefLevel().set(), MBEDGE, edges_to_refine);
437
438 CHKERR refine->addVerticesInTheMiddleOfEdges(edges_to_refine, bit_level2,
439 VERBOSE);
440 if (simpleInterface->getDim() == 3) {
441 CHKERR refine->refineTets(*meshset_level1_ptr, bit_level2, VERBOSE);
442 } else if (simpleInterface->getDim() == 2) {
443 CHKERR refine->refineTris(*meshset_level1_ptr, bit_level2, VERBOSE);
444 } else {
445 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
446 "Dimension not handled by test");
447 }
448
449 Range meshsets;
450 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets, true);
451 for (auto m : meshsets) {
453 ->updateMeshsetByEntitiesChildren(m, bit_level2, m, MBMAXTYPE, false);
454 }
455
457 };
458
459 CHKERR refine_mesh();
460
461 simpleInterface->getBitRefLevel() = bit_level1 | bit_level2;
463
465
466 CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
468 bit_level0 | bit_level1);
469
470 auto project_data = [&]() {
472
474
475 pipeline_mng->getDomainLhsFE().reset();
476 pipeline_mng->getDomainRhsFE().reset();
477 pipeline_mng->getBoundaryRhsFE().reset();
478
479 auto rule = [](int, int, int p) -> int { return 2 * p; };
480
481 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule);
482 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
483 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(rule);
484
485 auto test_bit_ref = [](FEMethod *fe_ptr) {
486 const auto &bit = fe_ptr->numeredEntFiniteElementPtr->getBitRefLevel();
487 MOFEM_LOG("SELF", Sev::noisy) << "ref : " << bit << " " << bit.test(2);
488 return bit.test(2);
489 };
490
491 pipeline_mng->getDomainLhsFE()->exeTestHook = test_bit_ref;
492 pipeline_mng->getDomainRhsFE()->exeTestHook = test_bit_ref;
493 pipeline_mng->getBoundaryRhsFE()->exeTestHook = test_bit_ref;
494
495 auto beta = [](const double, const double, const double) { return 1; };
496 auto field_vals_ptr = boost::make_shared<VectorDouble>();
497
498 auto domainParentRhs = boost::make_shared<DomainParentEle>(mField);
499 domainParentRhs->getOpPtrVector().push_back(
500 new OpCalculateScalarFieldValues(FIELD_NAME, field_vals_ptr));
501
502 pipeline_mng->getOpDomainLhsPipeline().push_back(
504 pipeline_mng->getOpDomainRhsPipeline().push_back(
505 new OpRunParent(domainParentRhs, bit_level2, bit_level2,
506 domainParentRhs, bit_level2, BitRefLevel().set()));
507
510 pipeline_mng->getOpDomainRhsPipeline().push_back(
511 new OpDomainTimesScalarField(FIELD_NAME, field_vals_ptr, beta));
512
513 pipeline_mng->getOpDomainRhsPipeline().push_back(
515
516 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
518
520
521 simpleInterface->getBitRefLevel() = bit_level2;
524
525 CHKERR checkResults([](FEMethod *fe_ptr) { return true; });
526
528 };
529
530 CHKERR project_data();
531
533}
534//! [Postprocess results]
535
536//! [Check results]
538 boost::function<bool(FEMethod *fe_method_ptr)> test_bit) {
541 pipeline_mng->getDomainLhsFE().reset();
542 pipeline_mng->getDomainRhsFE().reset();
543 pipeline_mng->getBoundaryRhsFE().reset();
544
545 auto rule = [](int, int, int p) -> int { return 2 * p + 1; };
546 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
547 pipeline_mng->getDomainRhsFE()->exeTestHook = test_bit;
548
549 auto common_data_ptr = boost::make_shared<CommonData>();
550 common_data_ptr->resVec = smartCreateDMVector(simpleInterface->getDM());
551 common_data_ptr->L2Vec = createSmartVectorMPI(
552 mField.get_comm(), (!mField.get_comm_rank()) ? 1 : 0, 1);
553 common_data_ptr->approxVals = boost::make_shared<VectorDouble>();
554
555 pipeline_mng->getOpDomainRhsPipeline().push_back(
557 common_data_ptr->approxVals));
558 pipeline_mng->getOpDomainRhsPipeline().push_back(
559 new OpError<FIELD_DIM>(common_data_ptr));
560
561 CHKERR VecZeroEntries(common_data_ptr->L2Vec);
562 CHKERR VecZeroEntries(common_data_ptr->resVec);
563
564 CHKERR pipeline_mng->loopFiniteElements();
565
566 CHKERR VecAssemblyBegin(common_data_ptr->L2Vec);
567 CHKERR VecAssemblyEnd(common_data_ptr->L2Vec);
568 CHKERR VecAssemblyBegin(common_data_ptr->resVec);
569 CHKERR VecAssemblyEnd(common_data_ptr->resVec);
570 double nrm2;
571 CHKERR VecNorm(common_data_ptr->resVec, NORM_2, &nrm2);
572 const double *array;
573 CHKERR VecGetArrayRead(common_data_ptr->L2Vec, &array);
574 MOFEM_LOG_C("WORLD", Sev::inform, "Error %6.4e Vec norm %6.4e\n",
575 std::sqrt(array[0]), nrm2);
576 CHKERR VecRestoreArrayRead(common_data_ptr->L2Vec, &array);
577
578 constexpr double eps = 1e-8;
579 if (nrm2 > eps)
580 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
581 "Not converged solution err = %6.4e", nrm2);
583}
584//! [Check results]
585
586int main(int argc, char *argv[]) {
587
588 // Initialisation of MoFEM/PETSc and MOAB data structures
589 MoFEM::Core::Initialize(&argc, &argv, NULL, help);
590
591 try {
592
593 //! [Register MoFEM discrete manager in PETSc]
594 DMType dm_name = "DMMOFEM";
595 CHKERR DMRegister_MoFEM(dm_name);
596 //! [Register MoFEM discrete manager in PETSc
597
598 //! [Create MoAB]
599 moab::Core mb_instance; ///< mesh database
600 moab::Interface &moab = mb_instance; ///< mesh database interface
601 //! [Create MoAB]
602
603 //! [Create MoFEM]
604 MoFEM::Core core(moab); ///< finite element database
605 MoFEM::Interface &m_field = core; ///< finite element database insterface
606 //! [Create MoFEM]
607
608 //! [AtomTest]
609 AtomTest ex(m_field);
610 CHKERR ex.runProblem();
611 //! [AtomTest]
612 }
614
616}
static Index< 'n', 3 > n
static Index< 'p', 3 > p
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
MoFEM::FaceElementForcesAndSourcesCore FaceEle
static const double eps
ElementsAndOps< SPACE_DIM >::DomianParentEle DomainParentEle
int main(int argc, char *argv[])
[Check results]
EntitiesFieldData::EntData EntData
ElementsAndOps< SPACE_DIM >::BoundaryParentEle BoundaryParentEle
static char help[]
constexpr char FIELD_NAME[]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, FIELD_DIM > OpDomainSource
constexpr int SPACE_DIM
boost::shared_ptr< DomainEle > domainChildRhs
constexpr int FIELD_DIM
boost::shared_ptr< DomainEle > domainChildLhs
[Set up problem]
DomainEle::UserDataOperator DomainEleOp
Finire element operator type.
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, FIELD_DIM > OpDomainMass
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
@ VERBOSE
Definition: definitions.h:222
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ H1
continuous field
Definition: definitions.h:98
@ NOSPACE
Definition: definitions.h:96
#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
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#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
constexpr double beta
constexpr double alpha
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:482
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:977
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
boost::ptr_vector< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
SmartPetscObj< KSP > createKSP(SmartPetscObj< DM > dm=nullptr)
Create KSP (linear) solver.
boost::ptr_vector< UserDataOperator > & getOpBoundaryRhsPipeline()
Get the Op Boundary Rhs Pipeline object.
boost::ptr_vector< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
auto test_bit_child
lambda function used to select elements on which finite element pipelines are executed.
auto bit
set bit
FTensor::Index< 'i', SPACE_DIM > i
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalarField< 1 > OpDomainTimesScalarField
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
auto createSmartVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
CoreTmp< 0 > Core
Definition: Core.hpp:1096
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:149
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temprary meshset.
Definition: Templates.hpp:1479
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
const double D
diffusivity
const double r
rate factor
FTensor::Index< 'm', 3 > m
double operator()(const double x, const double y, const double z)
Collected data use d by operator to evaluate errors for the test.
boost::shared_ptr< VectorDouble > approxVals
SmartPetscObj< Vec > L2Vec
SmartPetscObj< Vec > resVec
OpError(boost::shared_ptr< CommonData > &common_data_ptr)
boost::shared_ptr< CommonData > commonDataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Operator to evaluate errors.
boost::shared_ptr< CommonData > commonDataPtr
Simple * simpleInterface
MoFEMErrorCode checkResults()
[Check results]
MoFEMErrorCode solveSystem()
[Push operators to pipeline]
static ApproxFieldFunction< FIELD_DIM > approxFunction
AtomTest(MoFEM::Interface &m_field)
MoFEMErrorCode setupProblem()
[Read mesh]
MoFEMErrorCode refineResults()
[Solve]
MoFEMErrorCode readMesh()
[Run programme]
MoFEMErrorCode assembleSystem()
[Push operators to pipeline]
MoFEM::Interface & mField
MoFEMErrorCode runProblem()
[Run programme]
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
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
base operator to do operations at Gauss Pt. level
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
const VectorInt & getIndices() const
Get global indices of dofs on entity.
structure for User Loop Methods on finite elements
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
auto getFTensor0IntegrationWeight()
Get integration weights.
double getMeasure() const
get measure of element
@ OPROW
operator doWork function is executed on FE rows
@ OPSPACE
operator do Work is execute on space data
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
Mesh refinement interface.
Get value at integration points for scalar field.
Operator to execute finite element instance on parent element. This operator is typically used to pro...
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getDomainLhsFE()
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
boost::shared_ptr< FEMethod > & getBoundaryRhsFE()
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Problem manager is used to build and partition problems.
Simple interface for fast problem set-up.
Definition: Simple.hpp:35
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:272
MoFEMErrorCode addDomainField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:378
MoFEMErrorCode reSetUp(bool only_dm=false)
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:776
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:294
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:810
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:308
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:593
MoFEMErrorCode addBoundaryField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on boundary.
Definition: Simple.cpp:396
BitRefLevel & getBitRefLevelMask()
Get the BitRefLevel.
Definition: Simple.hpp:314
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:753
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:342
BitRefLevel & getBitRefLevel()
Get the BitRefLevel.
Definition: Simple.hpp:307
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)