v0.8.23
EntsMultiIndices.cpp
Go to the documentation of this file.
1 /** \file CoreDataStructures.cpp
2  * \brief Multi-index containers, data structures and other low-level functions
3  */
4 
5 /* MoFEM is free software: you can redistribute it and/or modify it under
6  * the terms of the GNU Lesser General Public License as published by the
7  * Free Software Foundation, either version 3 of the License, or (at your
8  * option) any later version.
9  *
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 #define IS_BUILDING_MB
20 #include <moab/Error.hpp>
21 
22 namespace MoFEM {
23 
24 static moab::Error error;
25 
26 inline void *get_tag_ptr(moab::Interface &moab, Tag th, EntityHandle ent,
27  int *tag_size) {
28  ApproximationOrder *ret_val;
29  rval = moab.tag_get_by_ptr(th, &ent, 1, (const void **)&ret_val, tag_size);
30  if (rval != MB_SUCCESS) {
31  *tag_size = 0;
32  return NULL;
33  } else {
34  return ret_val;
35  }
36 }
37 
38 BasicEntityData::BasicEntityData(const moab::Interface &moab,
39  const int pcomm_id)
40  : moab(const_cast<moab::Interface &>(moab)), pcommID(pcomm_id),
41  distributedMesh(true) {
42  rval = moab.tag_get_handle("_RefParentHandle", th_RefParentHandle);
44  rval = moab.tag_get_handle("_RefBitLevel", th_RefBitLevel);
46 }
48 
49 // basic moab ent
51  const boost::shared_ptr<BasicEntityData> &basic_data_ptr,
52  const EntityHandle ent)
53  : basicDataPtr(basic_data_ptr), ent(ent) {
54  switch (getEntType()) {
55  case MBVERTEX:
56  case MBEDGE:
57  case MBTRI:
58  case MBQUAD:
59  case MBTET:
60  case MBPRISM:
61  case MBENTITYSET:
62  break;
63  default:
64  THROW_MESSAGE("this entity type is currently not implemented");
65  }
66  ParallelComm *pcomm =
67  ParallelComm::get_pcomm(&basicDataPtr->moab, basicDataPtr->pcommID);
68  if (pcomm == NULL)
69  THROW_MESSAGE("pcomm is null");
70  rval = pcomm->get_owner_handle(ent, owner_proc, moab_owner_handle);
73 }
74 
75 unsigned char BasicEntity::getPStatus() const {
76  ParallelComm *pcomm =
77  ParallelComm::get_pcomm(&basicDataPtr->moab, basicDataPtr->pcommID);
78  return *((unsigned char *)MoFEM::get_tag_ptr(
79  basicDataPtr->moab, pcomm->pstatus_tag(), ent, NULL));
80 }
81 
82 // ref moab ent
84 RefEntity::RefEntity(const boost::shared_ptr<BasicEntityData> &basic_data_ptr,
85  const EntityHandle ent)
86  : BasicEntity(basic_data_ptr, ent) {}
87 
89  return static_cast<EntityHandle *>(get_tag_ptr(
90  basicDataPtr->moab, basicDataPtr->th_RefParentHandle, ent, NULL));
91 }
92 
94  return static_cast<BitRefLevel *>(
95  get_tag_ptr(basicDataPtr->moab, basicDataPtr->th_RefBitLevel, ent, NULL));
96 }
97 
98 MoFEMErrorCode getParentEnt(moab::Interface &moab, Range ents,
99  std::vector<EntityHandle> vec_patent_ent) {
100 
102  Tag th_ref_parent_handle;
103  CHKERR moab.tag_get_handle("_RefParentHandle", th_ref_parent_handle);
104  vec_patent_ent.resize(ents.size());
105  CHKERR moab.tag_get_data(th_ref_parent_handle, ents,
106  &*vec_patent_ent.begin());
108 }
109 
111 RefEntity::getBitRefLevel(moab::Interface &moab, Range ents,
112  std::vector<BitRefLevel> &vec_bit_ref_level) {
113 
115  Tag th_ref_bit_level;
116  CHKERR moab.tag_get_handle("_RefBitLevel", th_ref_bit_level);
117  vec_bit_ref_level.resize(ents.size());
118  CHKERR moab.tag_get_data(th_ref_bit_level, ents, &*vec_bit_ref_level.begin());
120 }
121 
123  moab::Interface &moab, Range ents,
124  std::vector<const BitRefLevel *> &vec_ptr_bit_ref_level) {
126  Tag th_ref_bit_level;
127  CHKERR moab.tag_get_handle("_RefBitLevel", th_ref_bit_level);
128  vec_ptr_bit_ref_level.resize(ents.size());
129  CHKERR moab.tag_get_by_ptr(
130  th_ref_bit_level, ents,
131  reinterpret_cast<const void **>(&*vec_ptr_bit_ref_level.begin()));
133 }
134 
135 std::ostream &operator<<(std::ostream &os, const RefEntity &e) {
136  os << "ent " << e.ent;
137  os << " pstatus " << std::bitset<8>(e.getPStatus());
138  os << " owner ent " << e.getOwnerEnt();
139  os << " owner proc " << e.getOwnerProc();
140  os << " parent ent " << e.getParentEnt();
141  // os << " BitRefLevel " << e.getBitRefLevel();
142  os << " ent type " << e.getEntType();
143  os << " ent parent type " << e.getParentEntType();
144  return os;
145 }
146 
147 // moab ent
149  const boost::shared_ptr<Field> &field_ptr,
150  const boost::shared_ptr<RefEntity> &ref_ent_ptr,
151  boost::shared_ptr<double *const> &&field_data_adaptor_ptr,
152  boost::shared_ptr<const int> &&t_max_order_ptr)
154  ref_ent_ptr),
155  tagMaxOrderPtr(t_max_order_ptr),
156  fieldDataAdaptorPtr(field_data_adaptor_ptr) {
158 
159  if (PetscUnlikely(!fieldDataAdaptorPtr))
160  THROW_MESSAGE("Pointer to field data adaptor not set");
161 
162  if (PetscUnlikely(!tagMaxOrderPtr))
163  THROW_MESSAGE("Pointer to max order not set");
164 }
165 
166 boost::shared_ptr<FieldData *const> FieldEntity::makeSharedFieldDataAdaptorPtr(
167  const boost::shared_ptr<Field> &field_ptr,
168  const boost::shared_ptr<RefEntity> &ref_ent_ptr) {
169  int size;
170  double *ptr;
171  switch (ref_ent_ptr->getEntType()) {
172  case MBVERTEX:
173  size = field_ptr->getNbOfCoeffs();
174  ptr = static_cast<FieldData *>(
175  MoFEM::get_tag_ptr(field_ptr->moab, field_ptr->th_FieldDataVerts,
176  ref_ent_ptr->ent, &size));
177  break;
178  default:
179  ptr = static_cast<FieldData *>(MoFEM::get_tag_ptr(
180  field_ptr->moab, field_ptr->th_FieldData, ref_ent_ptr->ent, &size));
181  }
182  return boost::make_shared<FieldData *const>(ptr);
183 }
184 
186 std::ostream &operator<<(std::ostream &os, const FieldEntity &e) {
187  os << "ent_global_uid "
188  << (UId)e.getGlobalUniqueId()
189  // << " ent_local_uid " << (UId)e.get_local_unique_id()
190  << " entity " << e.getEnt() << " type " << e.getEntType() << " pstatus "
191  << std::bitset<8>(e.getPStatus()) << " owner handle " << e.getOwnerEnt()
192  << " owner proc " << e.getOwnerProc() << " order " << e.getMaxOrder()
193  << " " << *e.sFieldPtr;
194  return os;
195 }
196 
198 
199  moab::Interface &moab = e->sPtr->basicDataPtr->moab;
200  const EntityHandle ent = e->getEnt();
201  *const_cast<ApproximationOrder *>(e->getMaxOrderPtr()) = order;
202  std::size_t nb_dofs = e->getOrderNbDofs(order) * e->getNbOfCoeffs();
203 
204  double *tag_field_data = nullptr;
205  int tag_field_data_size = 0;
206 
207  auto set_verts = [&]() {
208  if (e->sFieldPtr->th_FieldDataVertsType == MB_TAG_SPARSE) {
209  // Get pointer and size of field values tag
210  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
211  (const void **)&tag_field_data,
212  &tag_field_data_size);
213  if (nb_dofs) {
214  if (nb_dofs != tag_field_data_size) {
215  rval = moab.tag_set_data(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
216  &*data.begin());
217  MOAB_THROW(rval);
218  }
219  } else if (rval == MB_SUCCESS) {
220  rval = moab.tag_delete_data(e->sFieldPtr->th_FieldDataVerts, &ent, 1);
221  MOAB_THROW(rval);
222  }
223  } else {
224  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
225  (const void **)&tag_field_data);
226  MOAB_THROW(rval);
227  rval = moab.tag_set_data(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
228  tag_field_data);
229  MOAB_THROW(rval);
230  }
231  };
232 
233  auto set_default = [&]() {
234  if (reduceTagSize || nb_dofs) {
235 
236  // Get pointer and size of field values tag
237  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldData, &ent, 1,
238  (const void **)&tag_field_data,
239  &tag_field_data_size);
240 
241  if ((reduceTagSize && nb_dofs != tag_field_data_size) ||
242  nb_dofs > tag_field_data_size) {
243 
244  // Tag exist and are some data on it
245  if (rval == MB_SUCCESS) {
246  // Size of tag is different than new size, so copy data to new
247  // container
248  data.resize(tag_field_data_size);
249  FieldData *ptr_begin = static_cast<FieldData *>(tag_field_data);
250  FieldData *ptr_end =
251  static_cast<FieldData *>(tag_field_data) + tag_field_data_size;
252  std::copy(ptr_begin, ptr_end, data.begin());
253  }
254 
255  if (rval != MB_SUCCESS || nb_dofs) {
256 
257  // Set field dof data
258  data.resize(nb_dofs, 0);
259  int tag_size[1];
260  tag_size[0] = data.size();
261  void const *tag_data[] = {&data[0]};
262  rval = moab.tag_set_by_ptr(e->sFieldPtr->th_FieldData, &ent, 1,
263  tag_data, tag_size);
264  MOAB_THROW(rval);
265  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldData, &ent, 1,
266  (const void **)&tag_field_data,
267  &tag_field_data_size);
268  MOAB_THROW(rval);
269 
270  } else {
271 
272  rval = moab.tag_delete_data(e->sFieldPtr->th_FieldData, &ent, 1);
273  MOAB_THROW(rval);
274  }
275  }
276  }
277  };
278 
279  switch (e->getEntType()) {
280  case MBVERTEX:
281  set_verts();
282  break;
283  default:
284  set_default();
285  }
286 
287  if (nb_dofs)
288  const_cast<double *&>(*(e->getEntFieldDataPtr())) = tag_field_data;
289  else
290  const_cast<double *&>(*(e->getEntFieldDataPtr())) = nullptr;
291 }
292 
293 } // namespace MoFEM
const BitRefLevel & getBitRefLevel() const
Get entity ref bit refinement signature.
unsigned char getPStatus() const
get pstatus This tag stores various aspects of parallel status in bits; see also define following,...
boost::shared_ptr< T > sPtr
RefEntity(const boost::shared_ptr< BasicEntityData > &basic_data_ptr, const EntityHandle ent)
int getOrderNbDofs(int order) const
Get number of DOFs on entity for given order of approximation.
this struct keeps basic methods for moab entity
Provide data structure for (tensor) field approximation.The Field is intended to provide support for ...
EntityHandle * getParentEntPtr() const
Get pointer to parent entity tag.
static boost::shared_ptr< FieldData *const > makeSharedFieldDataAdaptorPtr(const boost::shared_ptr< Field > &field_ptr, const boost::shared_ptr< RefEntity > &ref_ent_ptr)
Return shared pointer to entity field data vector adaptor.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
UId getGlobalUniqueIdCalculate() const
Calculate global UId.
FieldCoefficientsNumber getNbOfCoeffs() const
BasicEntity(const boost::shared_ptr< BasicEntityData > &basic_data_ptr, const EntityHandle ent)
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:620
Struct keeps handle to entity in the field.
UId globalUId
Global unique id for this entity.
#define MOAB_THROW(a)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:602
Pointer interface for MoFEM::Field.
const UId & getGlobalUniqueId() const
Get global unique id.
Struct keeps handle to refined handle.
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
moab::Interface & moab
MoFEMErrorCode getParentEnt(moab::Interface &moab, Range ents, std::vector< EntityHandle > vec_patent_ent)
EntityHandle getOwnerEnt() const
std::bitset< BITREFEDGES_SIZE > BitRefEdges
Definition: Types.hpp:45
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:84
std::ostream & operator<<(std::ostream &os, const DataForcesAndSourcesCore::EntData &e)
const ApproximationOrder * getMaxOrderPtr() const
Get pinter to Tag keeping approximation order.
ApproximationOrder getMaxOrder() const
Get order set to the entity (Allocated tag size for such number)
int part_proc
this can be changed on distributed
EntityType getEntType() const
Get entity type.
int owner_proc
this never can not be changed if distributed mesh
EntityHandle getOwnerEnt() const
Owner handle on this or other processors.
boost::shared_ptr< FieldData *const > & getEntFieldDataPtr() const
Get shared ptr to vector adaptor pointing to the field tag data on entity.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
BasicEntityData(const moab::Interface &mfield, const int pcomm_id=MYPCOMM_INDEX)
static MoFEMErrorCode getParentEnt(Interface &moab, Range ents, std::vector< EntityHandle > vec_patent_ent)
EntityType getParentEntType() const
Get patent entity.
interface to RefEntity
static moab::Error error
EntityHandle moab_owner_handle
#define CHKERR
Inline error check.
Definition: definitions.h:596
unsigned char getPStatus() const
void operator()(boost::shared_ptr< FieldEntity > &e)
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
const ApproximationOrder order
FieldEntity(const boost::shared_ptr< Field > &field_ptr, const boost::shared_ptr< RefEntity > &ref_ent_ptr, boost::shared_ptr< double *const > &&field_data_adaptor_ptr, boost::shared_ptr< const int > &&t_max_order_ptr)
void * get_tag_ptr(moab::Interface &moab, Tag th, EntityHandle ent, int *tag_size)
boost::shared_ptr< FieldData *const > fieldDataAdaptorPtr
EntityHandle getEnt() const
Get entity handle.
int getOwnerProc() const
Get processor owning entity.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
boost::shared_ptr< BasicEntityData > basicDataPtr
static BitRefEdges DummyBitRefEdges
boost::shared_ptr< const ApproximationOrder > tagMaxOrderPtr
uint128_t UId
Unique Id.
Definition: Types.hpp:41
boost::shared_ptr< T > sFieldPtr
std::vector< FieldData > data
BitRefLevel * getBitRefLevelPtr() const
Get pointer to bit ref level tag.