v0.9.0
Public Member Functions | Public Attributes | List of all members
MyTransport Struct Reference

Application of mix transport data structure. More...

Inheritance diagram for MyTransport:
[legend]
Collaboration diagram for MyTransport:
[legend]

Public Member Functions

 MyTransport (MoFEM::Interface &m_field)
 
MoFEMErrorCode getSource (EntityHandle ent, const double x, const double y, const double z, double &flux)
 set source term More...
 
MoFEMErrorCode getBcOnValues (const EntityHandle ent, const double x, const double y, const double z, double &value)
 
MoFEMErrorCode getBcOnFluxes (const EntityHandle ent, const double x, const double y, const double z, double &flux)
 essential (Neumann) boundary condition (set fluxes) More...
 
 MyTransport (MoFEM::Interface &m_field, BcFluxMap &bc_flux_map)
 
MoFEMErrorCode getSource (EntityHandle ent, const double x, const double y, const double z, double &flux)
 set source term More...
 
MoFEMErrorCode getBcOnValues (const EntityHandle ent, const double x, const double y, const double z, double &value)
 natural (Dirihlet) boundary conditions (set values) More...
 
MoFEMErrorCode getBcOnFluxes (const EntityHandle ent, const double x, const double y, const double z, double &flux)
 essential (Neumann) boundary condition (set fluxes) More...
 
MoFEMErrorCode addBoundaryElements (BitRefLevel &ref_level)
 set-up boundary conditions More...
 
MoFEMErrorCode refineMesh (MixTransportElement &ufe, const int nb_levels, const int order)
 Refine mesh. More...
 
MoFEMErrorCode squashBits ()
 Squash bits of entities. More...
 
MoFEMErrorCode updateMeshsetsFieldsAndElements (const int nb_levels)
 update meshsets with new entities after mesh refinement More...
 
- Public Member Functions inherited from MixTransport::MixTransportElement
 MixTransportElement (MoFEM::Interface &m_field)
 construction of this data structure More...
 
virtual ~MixTransportElement ()
 destructor More...
 
MoFEMErrorCode getDirichletBCIndices (IS *is)
 get dof indices where essential boundary conditions are applied More...
 
virtual MoFEMErrorCode getResistivity (const EntityHandle ent, const double x, const double y, const double z, MatrixDouble3by3 &inv_k)
 natural (Dirichlet) boundary conditions (set values) More...
 
virtual MoFEMErrorCode getBcOnValues (const EntityHandle ent, const int gg, const double x, const double y, const double z, double &value)
 evaluate natural (Dirichlet) boundary conditions More...
 
MoFEMErrorCode addFields (const std::string &values, const std::string &fluxes, const int order)
 Add fields to database. More...
 
MoFEMErrorCode addFiniteElements (const std::string &fluxes_name, const std::string &values_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add finite elements More...
 
MoFEMErrorCode buildProblem (BitRefLevel &ref_level)
 Build problem. More...
 
MoFEMErrorCode postProc (const string out_file)
 Post process results. More...
 
MoFEMErrorCode createMatrices ()
 create matrices More...
 
MoFEMErrorCode solveLinearProblem ()
 solve problem More...
 
MoFEMErrorCode calculateResidual ()
 calculate residual More...
 
MoFEMErrorCode evaluateError ()
 Calculate error on elements. More...
 
MoFEMErrorCode destroyMatrices ()
 destroy matrices More...
 

Public Attributes

BcFluxMapbcFluxMap
 
EntityHandle lastEnt
 
double lastFlux
 
- Public Attributes inherited from MixTransport::MixTransportElement
MoFEM::InterfacemField
 
MyVolumeFE feVol
 Instance of volume element. More...
 
MyTriFE feTri
 Instance of surface element. More...
 
VectorDouble valuesAtGaussPts
 values at integration points on element More...
 
MatrixDouble valuesGradientAtGaussPts
 gradients at integration points on element More...
 
VectorDouble divergenceAtGaussPts
 divergence at integration points on element More...
 
MatrixDouble fluxesAtGaussPts
 fluxes at integration points on element More...
 
set< intbcIndices
 
std::map< int, BlockDatasetOfBlocks
 maps block set id with appropriate BlockData More...
 
Vec D
 
Vec D0
 
Vec F
 
Mat Aij
 
map< double, EntityHandleerrorMap
 
double sumErrorFlux
 
double sumErrorDiv
 
double sumErrorJump
 

Detailed Description

Application of mix transport data structure.

define sources and other stuff

MixTransportElement is a class collecting functions, operators and data for mix implementation of transport element. See there to learn how elements are created or how operators look like.

Some methods in MixTransportElement are abstract, f.e. user need to implement own source therm.

MixTransportElement is a class collecting functions, operators and data for mix implementation of transport element. See there to learn how elements are created or how operators look like.

Some methods in MixTransportElement are abstract, f.e. user need to implement own source therm.

Definition at line 36 of file mix_transport.cpp.

Constructor & Destructor Documentation

◆ MyTransport() [1/2]

MyTransport::MyTransport ( MoFEM::Interface m_field)

Definition at line 38 of file mix_transport.cpp.

38 : MixTransportElement(m_field){};
MixTransportElement(MoFEM::Interface &m_field)
construction of this data structure

◆ MyTransport() [2/2]

MyTransport::MyTransport ( MoFEM::Interface m_field,
BcFluxMap bc_flux_map 
)

Definition at line 65 of file h_adaptive_transport.cpp.

66  : MixTransportElement(m_field), bcFluxMap(bc_flux_map), lastEnt(0),
67  lastFlux(0) {}
BcFluxMap & bcFluxMap
EntityHandle lastEnt
MixTransportElement(MoFEM::Interface &m_field)
construction of this data structure

Member Function Documentation

◆ addBoundaryElements()

MoFEMErrorCode MyTransport::addBoundaryElements ( BitRefLevel ref_level)

set-up boundary conditions

Parameters
ref_levelmesh refinement level
Note
It is assumed that user would like to something non-standard with boundary conditions, have a own type of data structures to pass to functions calculating values and fluxes on boundary. For example BcFluxMap. That way this function is implemented here not in generic class MixTransportElement.
Returns
error code

Definition at line 141 of file h_adaptive_transport.cpp.

141  {
143  Range tets;
144  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
145  ref_level, BitRefLevel().set(), MBTET, tets);
146  Skinner skin(&mField.get_moab());
147  Range skin_faces; // skin faces from 3d ents
148  CHKERR skin.find_skin(0, tets, false, skin_faces);
149  // note: what is essential (dirichlet) is natural (neumann) for mix-FE
150  // compared to classical FE
151  Range natural_bc;
153  mField, NODESET | TEMPERATURESET, it)) {
154  Range tris;
155  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, tris,
156  true);
157  natural_bc.insert(tris.begin(), tris.end());
158  }
160  mField, SIDESET | HEATFLUXSET, it)) {
161  HeatFluxCubitBcData mydata;
162  CHKERR it->getBcDataStructure(mydata);
163  if (mydata.data.flag1 == 1) {
164  Range tris;
165  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, tris,
166  true);
167  bcFluxMap[it->getMeshsetId()].eNts = tris;
168  bcFluxMap[it->getMeshsetId()].fLux = mydata.data.value1;
169  // cerr << bcFluxMap[it->getMeshsetId()].eNts << endl;
170  // cerr << bcFluxMap[it->getMeshsetId()].fLux << endl;
171  }
172  }
173  Range essential_bc = subtract(skin_faces, natural_bc);
174  Range bit_tris;
175  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
176  ref_level, BitRefLevel().set(), MBTRI, bit_tris);
177  essential_bc = intersect(bit_tris, essential_bc);
178  natural_bc = intersect(bit_tris, natural_bc);
180  "MIX_BCFLUX");
182  "MIX_BCVALUE");
183  // CHKERR
184  // mField.add_ents_to_finite_element_by_type(skin_faces,MBTRI,"MIX_BCVALUE");
186  }
BcFluxMap & bcFluxMap
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
Managing BitRefLevels.
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
Definition of the heat flux bc data structure.
Definition: BCData.hpp:415

◆ getBcOnFluxes() [1/2]

MoFEMErrorCode MyTransport::getBcOnFluxes ( const EntityHandle  ent,
const double  x,
const double  y,
const double  z,
double flux 
)
virtual

essential (Neumann) boundary condition (set fluxes)

Parameters
enthandle to finite element entity
xcoord
ycoord
zcoord
fluxreference to flux which is set by function
Returns
[description]

Reimplemented from MixTransport::MixTransportElement.

Definition at line 55 of file mix_transport.cpp.

56  {
58  flux = 0;
60  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ getBcOnFluxes() [2/2]

MoFEMErrorCode MyTransport::getBcOnFluxes ( const EntityHandle  ent,
const double  x,
const double  y,
const double  z,
double flux 
)
virtual

essential (Neumann) boundary condition (set fluxes)

Parameters
enthandle to finite element entity
xcoord
ycoord
zcoord
fluxreference to flux which is set by function
Returns
[description]

Reimplemented from MixTransport::MixTransportElement.

Definition at line 110 of file h_adaptive_transport.cpp.

111  {
113  if (lastEnt == ent) {
114  flux = lastFlux;
115  } else {
116  flux = 0;
117  for (BcFluxMap::iterator mit = bcFluxMap.begin(); mit != bcFluxMap.end();
118  mit++) {
119  Range &tris = mit->second.eNts;
120  if (tris.find(ent) != tris.end()) {
121  flux = mit->second.fLux;
122  }
123  }
124  lastEnt = ent;
125  lastFlux = flux;
126  }
128  }
BcFluxMap & bcFluxMap
EntityHandle lastEnt
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ getBcOnValues() [1/2]

MoFEMErrorCode MyTransport::getBcOnValues ( const EntityHandle  ent,
const double  x,
const double  y,
const double  z,
double value 
)

Definition at line 48 of file mix_transport.cpp.

49  {
51  value = 1;
53  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ getBcOnValues() [2/2]

MoFEMErrorCode MyTransport::getBcOnValues ( const EntityHandle  ent,
const double  x,
const double  y,
const double  z,
double value 
)

natural (Dirihlet) boundary conditions (set values)

Parameters
enthandle to finite element entity
xcoord
ycoord
zcoord
valuereference to value set by function
Returns
error code

Definition at line 94 of file h_adaptive_transport.cpp.

95  {
97  value = 0;
99  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ getSource() [1/2]

MoFEMErrorCode MyTransport::getSource ( EntityHandle  ent,
const double  x,
const double  y,
const double  z,
double flux 
)
virtual

set source term

Parameters
enthandle to entity on which function is evaluated
xcoord
ycoord
zcoord
fluxreference to source term set by function
Returns
error code

Reimplemented from MixTransport::MixTransportElement.

Definition at line 40 of file mix_transport.cpp.

41  {
43  // double d = sqrt(x*x+y*y+z*z);
44  flux = 1; //-pow(d,5./4.);
46  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ getSource() [2/2]

MoFEMErrorCode MyTransport::getSource ( EntityHandle  ent,
const double  x,
const double  y,
const double  z,
double flux 
)
virtual

set source term

Parameters
enthandle to entity on which function is evaluated
xcoord
ycoord
zcoord
fluxreference to source term set by function
Returns
error code

Reimplemented from MixTransport::MixTransportElement.

Definition at line 78 of file h_adaptive_transport.cpp.

79  {
81  flux = 0;
83  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ refineMesh()

MoFEMErrorCode MyTransport::refineMesh ( MixTransportElement ufe,
const int  nb_levels,
const int  order 
)

Refine mesh.

Parameters
ufegeneral data structure
nb_levelsnumber of refinement levels
orderset order of approximation
Returns
errpr code

Refinement of could result in distorted mesh, for example, imagine when you have two levels of non-uniform refinement. Some tetrahedra on the mesh at first refinement instance are only refined by splitting subset of edges on it. Then refined child tetrahedra usually will have worse quality than quality of parent element. Refining such element in subsequent mesh refinement, potentially will deteriorate elements quality even worse. To prevent that adding new refinement level, recreate whole hierarchy of meshes.

Note on subsequent improvement could include refinement of tetrahedra from different levels, including initial mesh. So refinement two could split elements created during refinement one and also split elements from an initial mesh.

That adding the new refinement level creates refinement hierarchy of meshes from a scratch, not adding to existing one.

Entities from previous hierarchy are used in that process, but bit levels on those entities are squashed.

Definition at line 215 of file h_adaptive_transport.cpp.

216  {
217  MeshRefinement *refine_ptr;
219  // get refined edges having child vertex
220  const RefEntity_multiIndex *ref_ents_ptr;
221  CHKERR mField.get_ref_ents(&ref_ents_ptr);
222  typedef RefEntity_multiIndex::index<
223  Composite_EntType_and_ParentEntType_mi_tag>::type RefEntsByComposite;
224  const RefEntsByComposite &ref_ents =
225  ref_ents_ptr->get<Composite_EntType_and_ParentEntType_mi_tag>();
226  RefEntsByComposite::iterator rit, hi_rit;
227  rit = ref_ents.lower_bound(boost::make_tuple(MBVERTEX, MBEDGE));
228  hi_rit = ref_ents.upper_bound(boost::make_tuple(MBVERTEX, MBEDGE));
229  Range refined_edges;
230  // thist loop is over vertices which parent is edge
231  for (; rit != hi_rit; rit++) {
232  refined_edges.insert((*rit)->getParentEnt()); // get parent edge
233  }
234  // get tets which has large error
235  Range tets_to_refine;
236  const double max_error = ufe.errorMap.rbegin()->first;
237  // int size = ((double)5/6)*ufe.errorMap.size();
238  for (map<double, EntityHandle>::iterator mit = ufe.errorMap.begin();
239  mit != ufe.errorMap.end(); mit++) {
240  // cerr << mit->first << " " << mit->second << endl;
241  // if((size--)>0) continue;
242  if (mit->first < 0.25 * max_error)
243  continue;
244  tets_to_refine.insert(mit->second);
245  }
246  Range tets_to_refine_edges;
247  CHKERR mField.get_moab().get_adjacencies(
248  tets_to_refine, 1, false, tets_to_refine_edges, moab::Interface::UNION);
249  refined_edges.merge(tets_to_refine_edges);
250  CHKERR mField.getInterface(refine_ptr);
251  for (int ll = 0; ll != nb_levels; ll++) {
252  Range edges;
253  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
254  BitRefLevel().set(ll), BitRefLevel().set(), MBEDGE, edges);
255  edges = intersect(edges, refined_edges);
256  // add edges to refine at current level edges (some of the where refined
257  // before)
259  edges, BitRefLevel().set(ll + 1));
260  // get tets at current level
261  Range tets;
262  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
263  BitRefLevel().set(ll), BitRefLevel().set(), MBTET, tets);
264  CHKERR refine_ptr->refine_TET(tets, BitRefLevel().set(ll + 1));
266  }
267 
268  // update fields and elements
269  EntityHandle ref_meshset;
270  CHKERR mField.get_moab().create_meshset(MESHSET_SET, ref_meshset);
271  {
272  // cerr << BitRefLevel().set(nb_levels) << endl;
273  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
274  BitRefLevel().set(nb_levels), BitRefLevel().set(), MBTET,
275  ref_meshset);
276 
277  Range ref_tets;
278  CHKERR mField.get_moab().get_entities_by_type(ref_meshset, MBTET,
279  ref_tets);
280 
281  // add entities to field
282  CHKERR mField.add_ents_to_field_by_type(ref_meshset, MBTET, "FLUXES");
283  CHKERR mField.add_ents_to_field_by_type(ref_meshset, MBTET, "VALUES");
284  CHKERR mField.set_field_order(0, MBTET, "FLUXES", order + 1);
285  CHKERR mField.set_field_order(0, MBTRI, "FLUXES", order + 1);
286  CHKERR mField.set_field_order(0, MBTET, "VALUES", order);
287 
288  // add entities to skeleton
289  Range ref_tris;
290  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
291  BitRefLevel().set(nb_levels), BitRefLevel().set(), MBTRI, ref_tris);
293  "MIX_SKELETON");
294 
295  // add entities to finite elements
297  mField, BLOCKSET | MAT_THERMALSET, it)) {
298  Mat_Thermal temp_data;
299  CHKERR it->getAttributeDataStructure(temp_data);
300  setOfBlocks[it->getMeshsetId()].cOnductivity =
301  temp_data.data.Conductivity;
302  setOfBlocks[it->getMeshsetId()].cApacity = temp_data.data.HeatCapacity;
303  CHKERR mField.get_moab().get_entities_by_type(
304  it->meshset, MBTET, setOfBlocks[it->getMeshsetId()].tEts, true);
305  setOfBlocks[it->getMeshsetId()].tEts =
306  intersect(ref_tets, setOfBlocks[it->getMeshsetId()].tEts);
308  setOfBlocks[it->getMeshsetId()].tEts, MBTET, "MIX");
309  }
310  }
311  CHKERR mField.get_moab().delete_entities(&ref_meshset, 1);
313  }
virtual MoFEMErrorCode get_ref_ents(const RefEntity_multiIndex **refined_ents_ptr) const =0
Get ref entities multi-index from database.
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode refine_TET(const EntityHandle meshset, const BitRefLevel &bit, const bool respect_interface=false, int verb=QUIET, Range *ref_edges=NULL, const bool debug=false)
refine TET in the meshset
virtual MoFEMErrorCode add_vertices_in_the_middle_of_edges(const EntityHandle meshset, const BitRefLevel &bit, const bool recursive=false, int verb=QUIET, EntityHandle start_v=0)
make vertices in the middle of edges in meshset and add them to refinement levels defined by bit
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
Thermal material data structure.
MoFEMErrorCode updateMeshsetsFieldsAndElements(const int nb_levels)
update meshsets with new entities after mesh refinement
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
Managing BitRefLevels.
block name is "MAT_THERMAL"
Definition: definitions.h:224
Mesh refinement interface.
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, member< RefEntity::BasicEntity, EntityHandle, &RefEntity::ent > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType > >, ordered_non_unique< tag< ParentEntType_mi_tag >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
map< double, EntityHandle > errorMap
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ squashBits()

MoFEMErrorCode MyTransport::squashBits ( )

Squash bits of entities.

Information about hierarchy of meshsets is lost, but entities are not deleted from the mesh. After squash entities bits, new hierarchy can be created.

Returns
error code

Definition at line 323 of file h_adaptive_transport.cpp.

323  {
325  BitRefLevel all_but_0;
326  all_but_0.set(0);
327  all_but_0.flip();
328  BitRefLevel garbage_bit;
329  garbage_bit.set(BITREFLEVEL_SIZE - 1); // Garbage level
330  const RefEntity_multiIndex *refined_ents_ptr;
331  CHKERR mField.get_ref_ents(&refined_ents_ptr);
332  RefEntity_multiIndex::iterator mit = refined_ents_ptr->begin();
333  for (; mit != refined_ents_ptr->end(); mit++) {
334  if (mit->get()->getEntType() == MBENTITYSET)
335  continue;
336  BitRefLevel bit = mit->get()->getBitRefLevel();
337  if ((all_but_0 & bit) == bit) {
338  *(const_cast<RefEntity *>(mit->get())->getBitRefLevelPtr()) =
339  garbage_bit;
340  } else {
341  *(const_cast<RefEntity *>(mit->get())->getBitRefLevelPtr()) =
342  BitRefLevel().set(0);
343  }
344  }
346  }
virtual MoFEMErrorCode get_ref_ents(const RefEntity_multiIndex **refined_ents_ptr) const =0
Get ref entities multi-index from database.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define BITREFLEVEL_SIZE
max number of refinements
Definition: definitions.h:279
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, member< RefEntity::BasicEntity, EntityHandle, &RefEntity::ent > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType > >, ordered_non_unique< tag< ParentEntType_mi_tag >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ updateMeshsetsFieldsAndElements()

MoFEMErrorCode MyTransport::updateMeshsetsFieldsAndElements ( const int  nb_levels)

update meshsets with new entities after mesh refinement

Parameters
nb_levelsnb_levels
orderappropriate order
Returns
error code

Definition at line 354 of file h_adaptive_transport.cpp.

354  {
355  BitRefLevel ref_level;
357  ref_level.set(nb_levels);
359  EntityHandle meshset = it->meshset;
361  ->updateMeshsetByEntitiesChildren(meshset, ref_level, meshset,
362  MBVERTEX, true);
364  ->updateMeshsetByEntitiesChildren(meshset, ref_level, meshset, MBEDGE,
365  true);
367  ->updateMeshsetByEntitiesChildren(meshset, ref_level, meshset, MBTRI,
368  true);
370  ->updateMeshsetByEntitiesChildren(meshset, ref_level, meshset, MBTET,
371  true);
372  }
374  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
Managing BitRefLevels.
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

Member Data Documentation

◆ bcFluxMap

BcFluxMap& MyTransport::bcFluxMap

Definition at line 61 of file h_adaptive_transport.cpp.

◆ lastEnt

EntityHandle MyTransport::lastEnt

Definition at line 62 of file h_adaptive_transport.cpp.

◆ lastFlux

double MyTransport::lastFlux

Definition at line 63 of file h_adaptive_transport.cpp.


The documentation for this struct was generated from the following files: