v0.14.0
Loading...
Searching...
No Matches
Problems

Adding and managing problems. More...

Collaboration diagram for Problems:

Make elemnts multishared

MoFEMErrorCode MoFEM::CommInterface::resolveSharedFiniteElements (const Problem *problem_ptr, const std::string &fe_name, int verb=DEFAULT_VERBOSITY)
 resolve shared entities for finite elements in the problem
 
MoFEMErrorCode MoFEM::CommInterface::resolveSharedFiniteElements (const std::string name, const std::string &fe_name, int verb=DEFAULT_VERBOSITY)
 resolve shared entities for finite elements in the problem
 

Problems

virtual MoFEMErrorCode MoFEM::CoreInterface::add_problem (const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
 Add problem.
 
virtual MoFEMErrorCode MoFEM::CoreInterface::delete_problem (const std::string name)=0
 Delete problem.
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_add_finite_element (const std::string name_problem, const std::string &fe_name)=0
 add finite element to problem, this add entities assigned to finite element to a particular problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_unset_finite_element (const std::string name_problem, const std::string &fe_name)=0
 unset finite element from problem, this remove entities assigned to finite element to a particular problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_ref_level_add_bit (const std::string &name_problem, const BitRefLevel &bit)=0
 add ref level to problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_mask_ref_level_add_bit (const std::string &name_problem, const BitRefLevel &bit)=0
 set dof mask ref level for problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_ref_level_set_bit (const std::string &name_problem, const BitRefLevel &bit)=0
 set ref level for problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_mask_ref_level_set_bit (const std::string &name_problem, const BitRefLevel &bit)=0
 set dof mask ref level for problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::list_problem () const =0
 list problems
 
virtual MoFEMErrorCode MoFEM::CoreInterface::clear_problem (const std::string name, int verb=DEFAULT_VERBOSITY)=0
 clear problem
 
virtual MoFEMErrorCode MoFEM::CoreInterface::clear_problems (int verb=DEFAULT_VERBOSITY)=0
 clear problems
 
virtual MoFEMErrorCode MoFEM::CoreInterface::get_problem_finite_elements_entities (const std::string name, const std::string &fe_name, const EntityHandle meshset)=0
 add finite elements to the meshset
 
virtual bool MoFEM::CoreInterface::check_problem (const std::string name)=0
 check if problem exist
 

Detailed Description

Adding and managing problems.

Function Documentation

◆ add_problem()

virtual MoFEMErrorCode MoFEM::CoreInterface::add_problem ( const std::string & name,
enum MoFEMTypes bh = MF_EXCL,
int verb = DEFAULT_VERBOSITY )
pure virtual

◆ check_problem()

virtual bool MoFEM::CoreInterface::check_problem ( const std::string name)
pure virtual

check if problem exist

Parameters
nameproblem name
Returns
true if problem is in database

Implemented in MoFEM::CoreTmp< 0 >.

◆ clear_problem()

virtual MoFEMErrorCode MoFEM::CoreInterface::clear_problem ( const std::string name,
int verb = DEFAULT_VERBOSITY )
pure virtual

clear problem

Implemented in MoFEM::CoreTmp< 0 >.

◆ clear_problems()

virtual MoFEMErrorCode MoFEM::CoreInterface::clear_problems ( int verb = DEFAULT_VERBOSITY)
pure virtual

clear problems

Implemented in MoFEM::CoreTmp< 0 >.

◆ delete_problem()

virtual MoFEMErrorCode MoFEM::CoreInterface::delete_problem ( const std::string name)
pure virtual

Delete problem.

Note
add_file is a collective, should be executed on all processors. Otherwise could lead to deadlock.

Implemented in MoFEM::CoreTmp< 0 >.

◆ get_problem_finite_elements_entities()

virtual MoFEMErrorCode MoFEM::CoreInterface::get_problem_finite_elements_entities ( const std::string name,
const std::string & fe_name,
const EntityHandle meshset )
pure virtual

add finite elements to the meshset

Parameters
nameis problem name
fe_name
meshset

Implemented in MoFEM::CoreTmp< 0 >.

◆ list_problem()

virtual MoFEMErrorCode MoFEM::CoreInterface::list_problem ( ) const
pure virtual

list problems

Implemented in MoFEM::CoreTmp< 0 >.

◆ modify_problem_add_finite_element()

virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_add_finite_element ( const std::string name_problem,
const std::string & fe_name )
pure virtual

◆ modify_problem_mask_ref_level_add_bit()

virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_mask_ref_level_add_bit ( const std::string & name_problem,
const BitRefLevel & bit )
pure virtual

set dof mask ref level for problem

Implemented in MoFEM::CoreTmp< 0 >.

◆ modify_problem_mask_ref_level_set_bit()

virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_mask_ref_level_set_bit ( const std::string & name_problem,
const BitRefLevel & bit )
pure virtual

set dof mask ref level for problem

Note
add_file is a collective, should be executed on all processors. Otherwise could lead to deadlock.

Implemented in MoFEM::CoreTmp< 0 >.

◆ modify_problem_ref_level_add_bit()

virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_ref_level_add_bit ( const std::string & name_problem,
const BitRefLevel & bit )
pure virtual

add ref level to problem

Note
add_file is a collective, should be executed on all processors. Otherwise could lead to deadlock.

if same finite element is solved using different level of refinements, than the level of refinement has to be specificied to problem in query

Parameters
nameProblem name
BitRefLevelbitLevel Example:
mField.modify_problem_add_finite_element("BEAM_BENDING_ON_MESH_REF1","ELASTIC");
mField.modify_problem_add_finite_element("BEAM_BENDING_ON_MESH_REF2","ELASTIC");
mField.modify_problem_ref_level_add_bit("BEAM_BENDING_ON_MESH_REF1",bit_level1);
mField.modify_problem_ref_level_add_bit("BEAM_BENDING_ON_MESH_REF2",bit_level2);
#define CHKERR
Inline error check.
Two Problems exist and solved independently, both are elastic, but solved using different mesh refinement

Implemented in MoFEM::CoreTmp< 0 >.

Examples
build_large_problem.cpp, build_problems.cpp, continuity_check_on_contact_prism_side_ele.cpp, continuity_check_on_skeleton_3d.cpp, forces_and_sources_testing_edge_element.cpp, forces_and_sources_testing_flat_prism_element.cpp, forces_and_sources_testing_users_base.cpp, hcurl_divergence_operator_2d.cpp, mesh_insert_interface_atom.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, nonlinear_dynamics.cpp, prism_elements_from_surface.cpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, remove_entities_from_problem.cpp, remove_entities_from_problem_not_partitioned.cpp, simple_contact.cpp, simple_contact_thermal.cpp, and test_jacobian_of_simple_contact_element.cpp.

◆ modify_problem_ref_level_set_bit()

virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_ref_level_set_bit ( const std::string & name_problem,
const BitRefLevel & bit )
pure virtual

set ref level for problem

Note
add_file is a collective, should be executed on all processors. Otherwise could lead to deadlock.

if same finite element is solved using different level of refinements, than the level of refinement has to be specificied to problem in query

Parameters
nameProblem name
BitRefLevelbitLevel Example:
mField.modify_problem_add_finite_element("BEAM_BENDING_ON_MESH_REF1","ELASTIC");
mField.modify_problem_add_finite_element("BEAM_BENDING_ON_MESH_REF2","ELASTIC");
mField.modify_problem_ref_level_set_bit("BEAM_BENDING_ON_MESH_REF1",bit_level1);
mField.modify_problem_ref_level_set_bit("BEAM_BENDING_ON_MESH_REF2",bit_level2);
Two Problems exist and solved independently, both are elastic, but solved using different mesh refinement

Implemented in MoFEM::CoreTmp< 0 >.

◆ modify_problem_unset_finite_element()

virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_unset_finite_element ( const std::string name_problem,
const std::string & fe_name )
pure virtual

unset finite element from problem, this remove entities assigned to finite element to a particular problem

Note: If problem is build, it need to be cleaned to make this effective

Note
add_file is a collective, should be executed on all processors. Otherwise could lead to deadlock.
Parameters
nameProblem name
nameFinite Element name

Implemented in MoFEM::CoreTmp< 0 >.

◆ resolveSharedFiniteElements() [1/2]

MoFEMErrorCode MoFEM::CommInterface::resolveSharedFiniteElements ( const Problem * problem_ptr,
const std::string & fe_name,
int verb = DEFAULT_VERBOSITY )

resolve shared entities for finite elements in the problem

Parameters
problem_ptrproblem pointer
fe_namefinite element name
verbverbosity level
Returns
error code

This allows for tag reduction or tag exchange, f.e.

CHKERR m_field.resolveSharedFiniteElements(problem_ptr,"SHELL_ELEMENT");
Tag th;
CHKERR mField.get_moab().tag_get_handle("ADAPT_ORDER",th);
CHKERR ParallelComm* pcomm =
ParallelComm::get_pcomm(&mField.get_moab(),MYPCOMM_INDEX);
CHKERR pcomm->exchange_tags(th,prisms);
#define MYPCOMM_INDEX
default communicator number PCOMM

Definition at line 550 of file CommInterface.cpp.

551 {
552 MoFEM::Interface &m_field = cOre;
554 ParallelComm *pcomm = ParallelComm::get_pcomm(
555 &m_field.get_moab(), m_field.get_basic_entity_data_ptr()->pcommID);
556 std::vector<int> shprocs(MAX_SHARING_PROCS, 0);
557 std::vector<EntityHandle> shhandles(MAX_SHARING_PROCS, 0);
558 Range ents;
559 Tag th_gid = m_field.get_moab().globalId_tag();
560 PetscLayout layout;
561 CHKERR problem_ptr->getNumberOfElementsByNameAndPart(m_field.get_comm(),
562 fe_name, &layout);
563 int gid, last_gid;
564 CHKERR PetscLayoutGetRange(layout, &gid, &last_gid);
565 CHKERR PetscLayoutDestroy(&layout);
566
567 const FiniteElement_multiIndex *fes_ptr;
568 CHKERR m_field.get_finite_elements(&fes_ptr);
569
570 auto fe_miit = fes_ptr->get<FiniteElement_name_mi_tag>().find(fe_name);
571 if (fe_miit != fes_ptr->get<FiniteElement_name_mi_tag>().end()) {
572 auto fit =
573 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>().lower_bound(
575 0, (*fe_miit)->getFEUId()));
576 auto hi_fe_it =
577 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>().upper_bound(
579 get_id_for_max_type<MBENTITYSET>(), (*fe_miit)->getFEUId()));
580 for (; fit != hi_fe_it; ++fit) {
581
582 auto ent = (*fit)->getEnt();
583 auto part = (*fit)->getPart();
584
585 ents.insert(ent);
586 CHKERR m_field.get_moab().tag_set_data(pcomm->part_tag(), &ent, 1, &part);
587 if (part == pcomm->rank()) {
588 CHKERR m_field.get_moab().tag_set_data(th_gid, &ent, 1, &gid);
589 gid++;
590 }
591 shprocs.clear();
592 shhandles.clear();
593
594 if (pcomm->size() > 1) {
595
596 unsigned char pstatus = 0;
597 if (pcomm->rank() != part) {
598 pstatus = PSTATUS_NOT_OWNED;
599 pstatus |= PSTATUS_GHOST;
600 }
601
602 if (pcomm->size() > 2) {
603 pstatus |= PSTATUS_SHARED;
604 pstatus |= PSTATUS_MULTISHARED;
605 } else {
606 pstatus |= PSTATUS_SHARED;
607 }
608
609 size_t rrr = 0;
610 for (size_t rr = 0; rr < pcomm->size(); ++rr) {
611 if (rr != pcomm->rank()) {
612 shhandles[rrr] = ent;
613 shprocs[rrr] = rr;
614 ++rrr;
615 }
616 }
617 for (; rrr != pcomm->size(); ++rrr)
618 shprocs[rrr] = -1;
619
620 if (pstatus & PSTATUS_SHARED) {
621 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedp_tag(), &ent, 1,
622 &shprocs[0]);
623 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedh_tag(), &ent, 1,
624 &shhandles[0]);
625 }
626
627 if (pstatus & PSTATUS_MULTISHARED) {
628 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedps_tag(), &ent, 1,
629 &shprocs[0]);
630 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedhs_tag(), &ent, 1,
631 &shhandles[0]);
632 }
633 CHKERR m_field.get_moab().tag_set_data(pcomm->pstatus_tag(), &ent, 1,
634 &pstatus);
635 }
636 }
637 }
638
639 CHKERR pcomm->exchange_tags(th_gid, ents);
641}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual boost::shared_ptr< BasicEntityData > & get_basic_entity_data_ptr()=0
Get pointer to basic entity data.
Deprecated interface functions.
UId getLocalUniqueIdCalculate() const
Generate UId for finite element entity.

◆ resolveSharedFiniteElements() [2/2]

MoFEMErrorCode MoFEM::CommInterface::resolveSharedFiniteElements ( const std::string name,
const std::string & fe_name,
int verb = DEFAULT_VERBOSITY )

resolve shared entities for finite elements in the problem

Parameters
nameproblem name
fe_namefinite element name
verbverbosity level
Returns
error code

This allows for tag reduction or tag exchange, f.e.

CHKERR m_field.resolveSharedFiniteElements(problem_ptr,"SHELL_ELEMENT");
Tag th;
CHKERR mField.get_moab().tag_get_handle("ADAPT_ORDER",th);
ParallelComm* pcomm =
ParallelComm::get_pcomm(&mField.get_moab(),MYPCOMM_INDEX);
// CHKERR pcomm->reduce_tags(th,MPI_SUM,prisms);
CHKERR pcomm->exchange_tags(th,prisms);

Definition at line 643 of file CommInterface.cpp.

644 {
645 MoFEM::Interface &m_field = cOre;
647 const Problem *problem_ptr;
648 CHKERR m_field.get_problem(name, &problem_ptr);
649 CHKERR resolveSharedFiniteElements(problem_ptr, fe_name, verb);
651}
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
MoFEMErrorCode resolveSharedFiniteElements(const Problem *problem_ptr, const std::string &fe_name, int verb=DEFAULT_VERBOSITY)
resolve shared entities for finite elements in the problem