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<VectorAdaptor> &&field_data_adaptor_ptr,
152  boost::shared_ptr<const int> &&t_max_order_ptr)
154  ref_ent_ptr),
155  fieldDataAdaptorPtr(field_data_adaptor_ptr),
156  tagMaxOrderPtr(t_max_order_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<VectorAdaptor> 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<double *>(MoFEM::get_tag_ptr(field_ptr->moab,
175  field_ptr->th_FieldDataVerts,
176  ref_ent_ptr->ent, &size));
177  break;
178  default:
179  ptr = static_cast<double *>(MoFEM::get_tag_ptr(
180  field_ptr->moab, field_ptr->th_FieldData, ref_ent_ptr->ent, &size));
181  }
182  return boost::make_shared<VectorAdaptor>(
183  size, ublas::shallow_array_adaptor<double>(size, ptr));
184 }
185 
187 std::ostream &operator<<(std::ostream &os, const FieldEntity &e) {
188  os << "ent_global_uid "
189  << (UId)e.getGlobalUniqueId()
190  // << " ent_local_uid " << (UId)e.get_local_unique_id()
191  << " entity " << e.getEnt() << " type " << e.getEntType() << " pstatus "
192  << std::bitset<8>(e.getPStatus()) << " owner handle " << e.getOwnerEnt()
193  << " owner proc " << e.getOwnerProc() << " order " << e.getMaxOrder()
194  << " " << *e.sFieldPtr;
195  return os;
196 }
197 
199 
200  moab::Interface &moab = e->sPtr->basicDataPtr->moab;
201  const EntityHandle ent = e->getEnt();
202  *const_cast<ApproximationOrder *>(e->getMaxOrderPtr()) = order;
203  std::size_t nb_dofs = e->getOrderNbDofs(order) * e->getNbOfCoeffs();
204 
205  double *tag_field_data;
206  int tag_field_data_size;
207 
208  auto set_verts = [&]() {
209  if (e->sFieldPtr->th_FieldDataVertsType == MB_TAG_SPARSE) {
210  // Get pointer and size of field values tag
211  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
212  (const void **)&tag_field_data,
213  &tag_field_data_size);
214  if (nb_dofs) {
215  if (nb_dofs != tag_field_data_size) {
216  rval = moab.tag_set_data(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
217  &*data.begin());
218  MOAB_THROW(rval);
219  }
220  } else if (rval == MB_SUCCESS) {
221  rval = moab.tag_delete_data(e->sFieldPtr->th_FieldDataVerts, &ent, 1);
222  MOAB_THROW(rval);
223  }
224  } else {
225  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
226  (const void **)&tag_field_data);
227  MOAB_THROW(rval);
228  rval = moab.tag_set_data(e->sFieldPtr->th_FieldDataVerts, &ent, 1,
229  tag_field_data);
230  MOAB_THROW(rval);
231  }
232  };
233 
234  auto set_default = [&]() {
235  // Get pointer and size of field values tag
236  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldData, &ent, 1,
237  (const void **)&tag_field_data,
238  &tag_field_data_size);
239  // Tag exist and are some data on it
240  if (rval == MB_SUCCESS) {
241  // Check if size of filed values tag is correct
242  if (nb_dofs <= (unsigned int)tag_field_data_size)
243  return;
244  else if (nb_dofs == 0) {
245  // Delete data on this entity
246  rval = moab.tag_delete_data(e->sFieldPtr->th_FieldData, &ent, 1);
247  MOAB_THROW(rval);
248  return;
249  }
250  // Size of tag is different than new seize, so copy data to new
251  // container
252  data.resize(tag_field_data_size);
253  FieldData *ptr_begin = static_cast<FieldData *>(tag_field_data);
254  FieldData *ptr_end =
255  static_cast<FieldData *>(tag_field_data) + tag_field_data_size;
256  std::copy(ptr_begin, ptr_end, data.begin());
257  }
258  // Set new data
259  if (nb_dofs > 0) {
260  // Set field dof data
261  data.resize(nb_dofs, 0);
262  int tag_size[1];
263  tag_size[0] = data.size();
264  void const *tag_data[] = {&data[0]};
265  rval = moab.tag_set_by_ptr(e->sFieldPtr->th_FieldData, &ent, 1, tag_data,
266  tag_size);
267  MOAB_THROW(rval);
268  rval = moab.tag_get_by_ptr(e->sFieldPtr->th_FieldData, &ent, 1,
269  (const void **)&tag_field_data,
270  &tag_field_data_size);
271  MOAB_THROW(rval);
272  }
273  };
274 
275  switch (e->getEntType()) {
276  case MBVERTEX:
277  set_verts();
278  break;
279  default:
280  set_default();
281  }
282 
283  if (nb_dofs)
285  tag_field_data_size, ublas::shallow_array_adaptor<double>(
286  tag_field_data_size, tag_field_data));
287  else
288  *(e->getEntFieldDataPtr()) =
289  VectorAdaptor(0, ublas::shallow_array_adaptor<double>(0, nullptr));
290 }
291 
292 } // 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
boost::shared_ptr< const int > tagMaxOrderPtr
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< VectorAdaptor > 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:476
UId getGlobalUniqueIdCalculate() const
Calculate global UId.
FieldCoefficientsNumber getNbOfCoeffs() const
BasicEntity(const boost::shared_ptr< BasicEntityData > &basic_data_ptr, const EntityHandle ent)
boost::shared_ptr< VectorAdaptor > fieldDataAdaptorPtr
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:619
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:601
Pointer interface for MoFEM::Field.
const UId & getGlobalUniqueId() const
Get global unique id.
boost::shared_ptr< VectorAdaptor > & getEntFieldDataPtr() const
Get shared ptr to vector adaptor pointing to the field tag data on entity.
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.
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)
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:104
EntityType getParentEntType() const
Get patent entity.
interface to RefEntity
static moab::Error error
FieldEntity(const boost::shared_ptr< Field > &field_ptr, const boost::shared_ptr< RefEntity > &ref_ent_ptr, boost::shared_ptr< VectorAdaptor > &&field_data_adaptor_ptr, boost::shared_ptr< const int > &&t_max_order_ptr)
EntityHandle moab_owner_handle
#define CHKERR
Inline error check.
Definition: definitions.h:595
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
void * get_tag_ptr(moab::Interface &moab, Tag th, EntityHandle ent, int *tag_size)
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:406
boost::shared_ptr< BasicEntityData > basicDataPtr
static BitRefEdges DummyBitRefEdges
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.