v0.10.0
Classes | Typedefs | Functions | Variables
boundary_marker.cpp File Reference
#include <MoFEM.hpp>

Go to the source code of this file.

Classes

struct  OpFace
 

Typedefs

using FaceEle = MoFEM::FaceElementForcesAndSourcesCore
 
using FaceEleOp = FaceEle::UserDataOperator
 

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 

Typedef Documentation

◆ FaceEle

◆ FaceEleOp

Examples
boundary_marker.cpp.

Definition at line 31 of file boundary_marker.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)
Examples
boundary_marker.cpp.

Definition at line 85 of file boundary_marker.cpp.

85  {
86 
87  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
88 
89  try {
90 
91  moab::Core mb_instance;
92  moab::Interface &moab = mb_instance;
93  ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
94  if (pcomm == NULL)
95  pcomm = new ParallelComm(&moab, PETSC_COMM_WORLD);
96 
97  // Create MoFEM instance
98  MoFEM::Core core(moab);
99  MoFEM::Interface &m_field = core;
100 
101  int order = 4;
102  CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
103 
104  DMType dm_name = "DMMOFEM";
105  CHKERR DMRegister_MoFEM(dm_name);
106 
107  auto *simple_interface = m_field.getInterface<Simple>();
108  CHKERR simple_interface->getOptions();
109  CHKERR simple_interface->loadFile("");
110  CHKERR simple_interface->addDomainField("FIELD1", H1,
112  CHKERR simple_interface->setFieldOrder("FIELD1", order);
113  CHKERR simple_interface->setUp();
114 
115  auto get_ents_on_mesh_skin = [&]() {
116  Range faces;
117  CHKERR m_field.get_moab().get_entities_by_type(0, MBTRI, faces);
118  Skinner skin(&m_field.get_moab());
119  Range skin_edges;
120  CHKERR skin.find_skin(0, faces, false, skin_edges);
121  Range skin_verts;
122  CHKERR moab.get_connectivity(skin_edges, skin_verts, true);
123  skin_edges.merge(skin_verts);
124  return skin_edges;
125  };
126 
127  auto mark_boundary_dofs = [&](Range &skin_edges) {
128  auto problem_manager = m_field.getInterface<ProblemsManager>();
129  std::vector<bool> marker;
130  problem_manager->markDofs(simple_interface->getProblemName(), ROW,
131  skin_edges, marker);
132  return marker;
133  };
134 
135  auto skin_ents = get_ents_on_mesh_skin();
136 
137  // Get global local vector of marked DOFs. Is global, since is set for all
138  // DOFs on processor. Is local since only DOFs on processor are in the
139  // vector. To access DOFs use local indices.
140  auto marker = mark_boundary_dofs(skin_ents);
141 
142  boost::shared_ptr<FaceEle> fe(new FaceEle(m_field));
143  fe->getOpPtrVector().push_back(new OpFace(skin_ents, marker));
144 
145  auto dm = simple_interface->getDM();
146  CHKERR DMoFEMLoopFiniteElements(dm, simple_interface->getDomainFEName(),
147  fe);
148  }
149  CATCH_ERRORS;
150 
152 }
Deprecated interface functions.
Problem manager is used to build and partition problems.
virtual moab::Interface & get_moab()=0
CoreTmp< 0 > Core
Definition: Core.hpp:1129
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:60
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:507
Simple interface for fast problem set-up.
Definition: Simple.hpp:36
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
constexpr int order
MoFEM::FaceElementForcesAndSourcesCore FaceEle
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:48
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
#define CHKERR
Inline error check.
Definition: definitions.h:604
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:292
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1943
MoFEMErrorCode markDofs(const std::string problem_name, RowColData rc, const Range ents, std::vector< bool > &marker)
Create vector with marked indices.
Core (interface) class.
Definition: Core.hpp:77
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:441
static char help[]
continuous field
Definition: definitions.h:177
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:100

Variable Documentation

◆ help

char help[] = "...\n\n"
static
Examples
boundary_marker.cpp.

Definition at line 28 of file boundary_marker.cpp.