v0.13.2
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 More...
 
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 More...
 

Problems

virtual MoFEMErrorCode MoFEM::CoreInterface::add_problem (const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
 Add problem. More...
 
virtual MoFEMErrorCode MoFEM::CoreInterface::delete_problem (const std::string name)=0
 Delete problem. More...
 
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 More...
 
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 More...
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_ref_level_add_bit (const std::string &name_problem, const BitRefLevel &bit)=0
 add ref level to problem More...
 
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 More...
 
virtual MoFEMErrorCode MoFEM::CoreInterface::modify_problem_ref_level_set_bit (const std::string &name_problem, const BitRefLevel &bit)=0
 set ref level for problem More...
 
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 More...
 
virtual MoFEMErrorCode MoFEM::CoreInterface::list_problem () const =0
 list problems More...
 
virtual MoFEMErrorCode MoFEM::CoreInterface::clear_problem (const std::string name, int verb=DEFAULT_VERBOSITY)=0
 clear problem More...
 
virtual MoFEMErrorCode MoFEM::CoreInterface::clear_problems (int verb=DEFAULT_VERBOSITY)=0
 clear problems More...
 
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 More...
 
virtual bool MoFEM::CoreInterface::check_problem (const std::string name)=0
 check if problem exist More...
 

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.
Definition: definitions.h:535
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: definitions.h:215

Definition at line 529 of file CommInterface.cpp.

530 {
531 MoFEM::Interface &m_field = cOre;
533 ParallelComm *pcomm = ParallelComm::get_pcomm(
534 &m_field.get_moab(), m_field.get_basic_entity_data_ptr()->pcommID);
535 std::vector<int> shprocs(MAX_SHARING_PROCS, 0);
536 std::vector<EntityHandle> shhandles(MAX_SHARING_PROCS, 0);
537 Range ents;
538 Tag th_gid = m_field.get_moab().globalId_tag();
539 PetscLayout layout;
540 CHKERR problem_ptr->getNumberOfElementsByNameAndPart(m_field.get_comm(),
541 fe_name, &layout);
542 int gid, last_gid;
543 CHKERR PetscLayoutGetRange(layout, &gid, &last_gid);
544 CHKERR PetscLayoutDestroy(&layout);
545
546 const FiniteElement_multiIndex *fes_ptr;
547 CHKERR m_field.get_finite_elements(&fes_ptr);
548
549 auto fe_miit = fes_ptr->get<FiniteElement_name_mi_tag>().find(fe_name);
550 if (fe_miit != fes_ptr->get<FiniteElement_name_mi_tag>().end()) {
551 auto fit =
552 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>().lower_bound(
554 0, (*fe_miit)->getFEUId()));
555 auto hi_fe_it =
556 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>().upper_bound(
558 get_id_for_max_type<MBENTITYSET>(), (*fe_miit)->getFEUId()));
559 for (; fit != hi_fe_it; ++fit) {
560
561 auto ent = (*fit)->getEnt();
562 auto part = (*fit)->getPart();
563
564 ents.insert(ent);
565 CHKERR m_field.get_moab().tag_set_data(pcomm->part_tag(), &ent, 1, &part);
566 if (part == pcomm->rank()) {
567 CHKERR m_field.get_moab().tag_set_data(th_gid, &ent, 1, &gid);
568 gid++;
569 }
570 shprocs.clear();
571 shhandles.clear();
572
573 if (pcomm->size() > 1) {
574
575 unsigned char pstatus = 0;
576 if (pcomm->rank() != part) {
577 pstatus = PSTATUS_NOT_OWNED;
578 pstatus |= PSTATUS_GHOST;
579 }
580
581 if (pcomm->size() > 2) {
582 pstatus |= PSTATUS_SHARED;
583 pstatus |= PSTATUS_MULTISHARED;
584 } else {
585 pstatus |= PSTATUS_SHARED;
586 }
587
588 size_t rrr = 0;
589 for (size_t rr = 0; rr < pcomm->size(); ++rr) {
590 if (rr != pcomm->rank()) {
591 shhandles[rrr] = ent;
592 shprocs[rrr] = rr;
593 ++rrr;
594 }
595 }
596 for (; rrr != pcomm->size(); ++rrr)
597 shprocs[rrr] = -1;
598
599 if (pstatus & PSTATUS_SHARED) {
600 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedp_tag(), &ent, 1,
601 &shprocs[0]);
602 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedh_tag(), &ent, 1,
603 &shhandles[0]);
604 }
605
606 if (pstatus & PSTATUS_MULTISHARED) {
607 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedps_tag(), &ent, 1,
608 &shprocs[0]);
609 CHKERR m_field.get_moab().tag_set_data(pcomm->sharedhs_tag(), &ent, 1,
610 &shhandles[0]);
611 }
612 CHKERR m_field.get_moab().tag_set_data(pcomm->pstatus_tag(), &ent, 1,
613 &pstatus);
614 }
615 }
616 }
617
618 CHKERR pcomm->exchange_tags(th_gid, ents);
620}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
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 622 of file CommInterface.cpp.

623 {
624 MoFEM::Interface &m_field = cOre;
626 const Problem *problem_ptr;
627 CHKERR m_field.get_problem(name, &problem_ptr);
628 CHKERR resolveSharedFiniteElements(problem_ptr, fe_name, verb);
630}
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