v0.13.1
Loading...
Searching...
No Matches
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#include <MoFEM.hpp>
9
10using namespace MoFEM;
11
12static char help[] = "...\n\n";
13
14constexpr char FIELD_NAME[] = "U";
15constexpr int FIELD_DIM = 1;
16constexpr int SPACE_DIM = 2;
17
18template <int DIM> struct ElementsAndOps {};
19
20template <> struct ElementsAndOps<2> {
27};
28
29template <> struct ElementsAndOps<3> {
32};
33
38
42
43template <int FIELD_DIM> struct ApproxFieldFunction;
44
45template <> struct ApproxFieldFunction<1> {
46 double operator()(const double x, const double y, const double z) {
47 return x * x + y * y + x * y + pow(x, 3) + pow(y, 3) + pow(x, 4) +
48 pow(y, 4);
49 }
50};
51
56
57struct AtomTest {
58
59 AtomTest(MoFEM::Interface &m_field) : mField(m_field) {}
60
62
63private:
66
68
74 checkResults(boost::function<bool(FEMethod *fe_method_ptr)> test_bit);
76 struct CommonData {
77 boost::shared_ptr<VectorDouble> approxVals;
80 };
81
82 template <int FIELD_DIM> struct OpError;
83};
84
87
88template <> struct AtomTest::OpError<1> : public DomainEleOp {
89 boost::shared_ptr<CommonData> commonDataPtr;
90 OpError(boost::shared_ptr<CommonData> &common_data_ptr)
91 : DomainEleOp(FIELD_NAME, OPROW), commonDataPtr(common_data_ptr) {}
92 MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
94
95 if (const size_t nb_dofs = data.getIndices().size()) {
96
97 const int nb_integration_pts = getGaussPts().size2();
99 auto t_val = getFTensor0FromVec(*(commonDataPtr->approxVals));
100 auto t_coords = getFTensor1CoordsAtGaussPts();
101
102 VectorDouble nf(nb_dofs, false);
103 nf.clear();
104
106 const double volume = getMeasure();
107
108 auto t_row_base = data.getFTensor0N();
109 double error = 0;
110 for (int gg = 0; gg != nb_integration_pts; ++gg) {
111
112 const double alpha = t_w * volume;
113 double diff = t_val - AtomTest::approxFunction(t_coords(0), t_coords(1),
114 t_coords(2));
115 error += alpha * pow(diff, 2);
116
117 for (size_t r = 0; r != nb_dofs; ++r) {
118 nf[r] += alpha * t_row_base * diff;
119 ++t_row_base;
120 }
121
122 ++t_w;
123 ++t_val;
124 ++t_coords;
125 }
126
127 const int index = 0;
128 CHKERR VecSetValue(commonDataPtr->L2Vec, index, error, ADD_VALUES);
129 CHKERR VecSetValues(commonDataPtr->resVec, data, &nf[0], ADD_VALUES);
130 }
131
133 }
134};
135
136template <typename ELE_OP, typename PARENT_ELE>
137struct OpCheckGaussCoords : public ELE_OP {
139
143
144 MatrixDouble parent_coords;
145
146 PARENT_ELE parent_fe(this->getPtrFE()->mField);
147 auto op = new ELE_OP(NOSPACE, ELE_OP::OPSPACE);
148 op->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
149 EntityType type,
152
153 MOFEM_LOG("SELF", Sev::noisy)
154 << "parent_coords in op "
155 << static_cast<ELE_OP *>(op_ptr)->getCoordsAtGaussPts();
156
157 parent_coords = static_cast<ELE_OP *>(op_ptr)->getCoordsAtGaussPts();
159 };
160 parent_fe.getOpPtrVector().push_back(op);
161
162 MOFEM_LOG("SELF", Sev::noisy) << "fe name " << this->getFEName();
163 CHKERR this->loopParent(this->getFEName(), &parent_fe);
164 MOFEM_LOG("SELF", Sev::noisy) << "parent_coords " << parent_coords;
165
166 MatrixDouble child_coords = this->getCoordsAtGaussPts();
167 MOFEM_LOG("SELF", Sev::noisy) << "child_coords " << child_coords;
168
169 child_coords -= parent_coords;
170
171 MOFEM_LOG("SELF", Sev::noisy) << "Corrds diffs" << child_coords;
172
173 double n = 0;
174 for (auto d : child_coords.data())
175 n += d * d;
176
177 if (sqrt(n) > 1e-12)
178 SETERRQ1(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
179 "Parent and child global coords at integration points are "
180 "diffrent norm = %3.2e",
181 sqrt(n));
182
184 }
185};
186
187//! [Run programme]
194
195 auto test_bit_child = [](FEMethod *fe_ptr) {
196 const auto &bit = fe_ptr->numeredEntFiniteElementPtr->getBitRefLevel();
197 MOFEM_LOG("SELF", Sev::noisy) << bit << " " << bit.test(0);
198 return bit.test(1);
199 };
200
204}
205//! [Run programme]
206
207//! [Read mesh]
210
214
215 MOFEM_LOG("WORLD", Sev::verbose) << "Dim " << simpleInterface->getDim();
216
217 auto bit_level0 = simpleInterface->getBitRefLevel();
218
219 auto &moab = mField.get_moab();
220
221 auto refine_mesh = [&](auto bit_level1) {
223
224 auto refine = mField.getInterface<MeshRefinement>();
225
226 auto meshset_level0_ptr = get_temp_meshset_ptr(moab);
227 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByRefLevel(
228 bit_level0, BitRefLevel().set(), *meshset_level0_ptr);
229
230 // random mesh refinement
231 auto meshset_ref_edges_ptr = get_temp_meshset_ptr(moab);
232 Range edges_to_refine;
233 CHKERR moab.get_entities_by_type(*meshset_level0_ptr, MBEDGE,
234 edges_to_refine);
235 int ii = 0;
236 for (Range::iterator eit = edges_to_refine.begin();
237 eit != edges_to_refine.end(); eit++, ii++) {
238 int numb = ii % 2;
239 if (numb == 0) {
240 CHKERR moab.add_entities(*meshset_ref_edges_ptr, &*eit, 1);
241 }
242 }
243 CHKERR refine->addVerticesInTheMiddleOfEdges(*meshset_ref_edges_ptr,
244 bit_level1, false, VERBOSE);
245 if (simpleInterface->getDim() == 3) {
246 CHKERR refine->refineTets(*meshset_level0_ptr, bit_level1, VERBOSE);
247 } else if (simpleInterface->getDim() == 2) {
248 CHKERR refine->refineTris(*meshset_level0_ptr, bit_level1, VERBOSE);
249 } else {
250 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
251 "Dimension not handled by test");
252 }
253
255 };
256
257 BitRefLevel bit_level1;
258 bit_level1.set(1);
259 CHKERR refine_mesh(bit_level1);
262
264}
265//! [Read mesh]
266
267//! [Set up problem]
270 // Add field
275 constexpr int order = 4;
278
279 CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
281 BitRefLevel().set(0));
282
284}
285//! [Set up problem]
286
287boost::shared_ptr<DomainEle> domainChildLhs, domainChildRhs;
288
289//! [Push operators to pipeline]
293
294 auto rule = [](int, int, int p) -> int { return 2 * p; };
295
296 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule);
297 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
298
299 auto test_bit_parent = [](FEMethod *fe_ptr) {
300 const auto &bit = fe_ptr->numeredEntFiniteElementPtr->getBitRefLevel();
301 MOFEM_LOG("SELF", Sev::noisy) << bit << " " << bit.test(0);
302 return bit.test(0);
303 };
304
305 pipeline_mng->getDomainLhsFE()->exeTestHook = test_bit_parent;
306 pipeline_mng->getDomainRhsFE()->exeTestHook = test_bit_parent;
307
308 auto beta = [](const double, const double, const double) { return 1; };
309
310 // Make aliased shared pointer, and create child element
311 domainChildLhs = boost::make_shared<DomainEle>(mField);
312 domainChildLhs->getRuleHook = rule;
313 domainChildLhs->getOpPtrVector().push_back(
315
316 domainChildRhs = boost::make_shared<DomainEle>(mField);
317 domainChildLhs->getRuleHook = rule;
318 domainChildRhs->getOpPtrVector().push_back(
320
321 auto parent_op_lhs = new DomainEleOp(NOSPACE, DomainEleOp::OPSPACE);
322 parent_op_lhs->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
323 EntityType type,
325 auto domain_op = static_cast<DomainEleOp *>(op_ptr);
327
328 MOFEM_LOG("SELF", Sev::noisy) << "LHS Pipeline FE";
329
330 if (!domainChildLhs)
331 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "FE not allocated");
332
333 auto &bit =
334 domain_op->getFEMethod()->numeredEntFiniteElementPtr->getBitRefLevel();
335 if (bit == BitRefLevel().set(0)) {
336 CHKERR domain_op->loopChildren(domain_op->getFEName(),
337 domainChildLhs.get(), VERBOSE, Sev::noisy);
338 } else {
339 CHKERR domain_op->loopThis(domain_op->getFEName(), domainChildLhs.get(),
340 VERBOSE, Sev::noisy);
341 }
343 };
344
345 auto parent_op_rhs = new DomainEleOp(NOSPACE, DomainEleOp::OPSPACE);
346 parent_op_rhs->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
347 EntityType type,
349 auto domain_op = static_cast<DomainEleOp *>(op_ptr);
351
352 MOFEM_LOG("SELF", Sev::noisy) << "RHS Pipeline FE";
353
354 if (!domainChildRhs)
355 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "FE not allocated");
356
357 auto &bit =
358 domain_op->getFEMethod()->numeredEntFiniteElementPtr->getBitRefLevel();
359 if (bit == BitRefLevel().set(0)) {
360 CHKERR domain_op->loopChildren(domain_op->getFEName(),
361 domainChildRhs.get(), VERBOSE, Sev::noisy);
362 } else if ((bit & BitRefLevel().set(0)).any()) {
363 CHKERR domain_op->loopThis(domain_op->getFEName(), domainChildRhs.get(),
364 VERBOSE, Sev::noisy);
365 }
367 };
368
369 pipeline_mng->getOpDomainLhsPipeline().push_back(parent_op_lhs);
370 pipeline_mng->getOpDomainRhsPipeline().push_back(parent_op_rhs);
371
373}
374//! [Push operators to pipeline]
375
376//! [Solve]
380
381 MOFEM_LOG("WORLD", Sev::inform) << "Solve problem";
382
383 auto solver = pipeline_mng->createKSP();
384 CHKERR KSPSetFromOptions(solver);
385 CHKERR KSPSetUp(solver);
386
387 auto dm = simpleInterface->getDM();
388 auto D = smartCreateDMVector(dm);
389 auto F = smartVectorDuplicate(D);
390
391 CHKERR KSPSolve(solver, F, D);
392 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
393 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
394 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
396}
397
398//! [Solve]
401
402 auto &moab = mField.get_moab();
403
404 auto bit_level0 = BitRefLevel().set(0);
405 auto bit_level1 = BitRefLevel().set(1);
406 auto bit_level2 = BitRefLevel().set(2);
407
408 auto refine_mesh = [&]() {
410
411 auto refine = mField.getInterface<MeshRefinement>();
412
413 auto meshset_level1_ptr = get_temp_meshset_ptr(moab);
414 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByDimAndRefLevel(
415 bit_level1, BitRefLevel().set(), simpleInterface->getDim(),
416 *meshset_level1_ptr);
417
418 // random mesh refinement
419 auto meshset_ref_edges_ptr = get_temp_meshset_ptr(moab);
420 Range edges_to_refine;
421 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
422 bit_level1, BitRefLevel().set(), MBEDGE, edges_to_refine);
423
424 CHKERR refine->addVerticesInTheMiddleOfEdges(edges_to_refine, bit_level2,
425 VERBOSE);
426 if (simpleInterface->getDim() == 3) {
427 CHKERR refine->refineTets(*meshset_level1_ptr, bit_level2, VERBOSE);
428 } else if (simpleInterface->getDim() == 2) {
429 CHKERR refine->refineTris(*meshset_level1_ptr, bit_level2, VERBOSE);
430 } else {
431 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
432 "Dimension not handled by test");
433 }
434
435 Range meshsets;
436 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets, true);
437 for (auto m : meshsets) {
439 ->updateMeshsetByEntitiesChildren(m, bit_level2, m, MBMAXTYPE, false);
440 }
441
443 };
444
445 CHKERR refine_mesh();
446
447 simpleInterface->getBitRefLevel() = bit_level1 | bit_level2;
449
451
452 CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
454 bit_level0 | bit_level1);
455
456 auto project_data = [&]() {
458
460
461 pipeline_mng->getDomainLhsFE().reset();
462 pipeline_mng->getDomainRhsFE().reset();
463 pipeline_mng->getBoundaryRhsFE().reset();
464
465 auto rule = [](int, int, int p) -> int { return 2 * p; };
466
467 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule);
468 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
469 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(rule);
470
471 auto test_bit_ref = [](FEMethod *fe_ptr) {
472 const auto &bit = fe_ptr->numeredEntFiniteElementPtr->getBitRefLevel();
473 MOFEM_LOG("SELF", Sev::noisy) << "ref : " << bit << " " << bit.test(2);
474 return bit.test(2);
475 };
476
477 pipeline_mng->getDomainLhsFE()->exeTestHook = test_bit_ref;
478 pipeline_mng->getDomainRhsFE()->exeTestHook = test_bit_ref;
479 pipeline_mng->getBoundaryRhsFE()->exeTestHook = test_bit_ref;
480
481 auto beta = [](const double, const double, const double) { return 1; };
482 auto field_vals_ptr = boost::make_shared<VectorDouble>();
483
484 auto domainParentRhs = boost::make_shared<DomainParentEle>(mField);
485 domainParentRhs->getOpPtrVector().push_back(
486 new OpCalculateScalarFieldValues(FIELD_NAME, field_vals_ptr));
487
488 pipeline_mng->getOpDomainLhsPipeline().push_back(
490 pipeline_mng->getOpDomainRhsPipeline().push_back(
491 new OpRunParent(domainParentRhs, bit_level2, bit_level2,
492 domainParentRhs, bit_level2, BitRefLevel().set()));
493
496 pipeline_mng->getOpDomainRhsPipeline().push_back(
497 new OpDomainTimesScalarField(FIELD_NAME, field_vals_ptr, beta));
498
499 pipeline_mng->getOpDomainRhsPipeline().push_back(
501
502 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
504
506
507 simpleInterface->getBitRefLevel() = bit_level2;
510
511 CHKERR checkResults([](FEMethod *fe_ptr) { return true; });
512
514 };
515
516 CHKERR project_data();
517
519}
520//! [Postprocess results]
521
522//! [Check results]
524 boost::function<bool(FEMethod *fe_method_ptr)> test_bit) {
527 pipeline_mng->getDomainLhsFE().reset();
528 pipeline_mng->getDomainRhsFE().reset();
529 pipeline_mng->getBoundaryRhsFE().reset();
530
531 auto rule = [](int, int, int p) -> int { return 2 * p + 1; };
532 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
533 pipeline_mng->getDomainRhsFE()->exeTestHook = test_bit;
534
535 auto common_data_ptr = boost::make_shared<CommonData>();
536 common_data_ptr->resVec = smartCreateDMVector(simpleInterface->getDM());
537 common_data_ptr->L2Vec = createSmartVectorMPI(
538 mField.get_comm(), (!mField.get_comm_rank()) ? 1 : 0, 1);
539 common_data_ptr->approxVals = boost::make_shared<VectorDouble>();
540
541 pipeline_mng->getOpDomainRhsPipeline().push_back(
543 common_data_ptr->approxVals));
544 pipeline_mng->getOpDomainRhsPipeline().push_back(
545 new OpError<FIELD_DIM>(common_data_ptr));
546
547 CHKERR VecZeroEntries(common_data_ptr->L2Vec);
548 CHKERR VecZeroEntries(common_data_ptr->resVec);
549
550 CHKERR pipeline_mng->loopFiniteElements();
551
552 CHKERR VecAssemblyBegin(common_data_ptr->L2Vec);
553 CHKERR VecAssemblyEnd(common_data_ptr->L2Vec);
554 CHKERR VecAssemblyBegin(common_data_ptr->resVec);
555 CHKERR VecAssemblyEnd(common_data_ptr->resVec);
556 double nrm2;
557 CHKERR VecNorm(common_data_ptr->resVec, NORM_2, &nrm2);
558 const double *array;
559 CHKERR VecGetArrayRead(common_data_ptr->L2Vec, &array);
560 MOFEM_LOG_C("WORLD", Sev::inform, "Error %6.4e Vec norm %6.4e\n",
561 std::sqrt(array[0]), nrm2);
562 CHKERR VecRestoreArrayRead(common_data_ptr->L2Vec, &array);
563
564 constexpr double eps = 1e-8;
565 if (nrm2 > eps)
566 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
567 "Not converged solution err = %6.4e", nrm2);
569}
570//! [Check results]
571
572int main(int argc, char *argv[]) {
573
574 // Initialisation of MoFEM/PETSc and MOAB data structures
575 MoFEM::Core::Initialize(&argc, &argv, NULL, help);
576
577 try {
578
579 //! [Register MoFEM discrete manager in PETSc]
580 DMType dm_name = "DMMOFEM";
581 CHKERR DMRegister_MoFEM(dm_name);
582 //! [Register MoFEM discrete manager in PETSc
583
584 //! [Create MoAB]
585 moab::Core mb_instance; ///< mesh database
586 moab::Interface &moab = mb_instance; ///< mesh database interface
587 //! [Create MoAB]
588
589 //! [Create MoFEM]
590 MoFEM::Core core(moab); ///< finite element database
591 MoFEM::Interface &m_field = core; ///< finite element database insterface
592 //! [Create MoFEM]
593
594 //! [AtomTest]
595 AtomTest ex(m_field);
596 CHKERR ex.runProblem();
597 //! [AtomTest]
598 }
600
602}
static Index< 'p', 3 > p
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:304
int main()
Definition: adol-c_atom.cpp:46
static const double eps
ElementsAndOps< SPACE_DIM >::DomianParentEle DomainParentEle
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]
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, FIELD_DIM > OpDomainMass
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
@ VERBOSE
Definition: definitions.h:209
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ H1
continuous field
Definition: definitions.h:85
@ NOSPACE
Definition: definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#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
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
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
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:965
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:301
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
constexpr char FIELD_NAME[]
constexpr int FIELD_DIM
DomainEle::UserDataOperator DomainEleOp
Finire element operator type.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< BASE_DIM, FIELD_DIM > OpDomainMass
OPerator to integrate mass matrix for least square approximation.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< BASE_DIM, FIELD_DIM > OpDomainSource
Operator to integrate the right hand side matrix for the problem.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1 > OpDomainTimesScalarField
double D
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
auto createSmartVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temprary meshset.
Definition: Templates.hpp:1584
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
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: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
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()
MoFEM::FaceElementForcesAndSourcesCore FaceEle
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
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:26
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:264
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:373
MoFEMErrorCode reSetUp(bool only_dm=false)
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:746
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:289
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:780
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:303
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:588
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:391
BitRefLevel & getBitRefLevelMask()
Get the BitRefLevel.
Definition: Simple.hpp:306
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:723
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:334
BitRefLevel & getBitRefLevel()
Get the BitRefLevel.
Definition: Simple.hpp:299
intrusive_ptr for managing petsc objects
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)