v0.14.0
ForcesAndSourcesCore.cpp
Go to the documentation of this file.
1 /** \file ForcesAndSourcesCore.cpp
2 
3 \brief Implementation of Elements on Entities for Forces and Sources
4 */
5 
6 #ifdef __cplusplus
7 extern "C" {
8 #endif
9 #include <cblas.h>
10 #include <lapack_wrap.h>
11 // #include <gm_rule.h>
12 #ifdef __cplusplus
13 }
14 #endif
15 
16 namespace MoFEM {
17 
18 static auto cmp_uid_lo(const boost::weak_ptr<FieldEntity> &a, const UId &b) {
19  if (auto a_ptr = a.lock()) {
20  if (a_ptr->getLocalUniqueId() < b)
21  return true;
22  else
23  return false;
24  } else {
25  return false;
26  }
27 }
28 
29 static auto cmp_uid_hi(const UId &b, const boost::weak_ptr<FieldEntity> &a) {
30  if (auto a_ptr = a.lock()) {
31  if (b < a_ptr->getLocalUniqueId())
32  return true;
33  else
34  return false;
35  } else {
36  return true;
37  }
38 }
39 
41  :
42 
43  mField(m_field), getRuleHook(0), setRuleHook(0),
44  dataOnElement{
45 
46  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // NOSPACE,
47  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // NOFIELD
48  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // H1
49  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HCURL
50  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HDIV
51  boost::make_shared<EntitiesFieldData>(MBENTITYSET) // L2
52 
53  },
54  derivedDataOnElement{
55 
56  nullptr,
57  boost::make_shared<DerivedEntitiesFieldData>(
58  dataOnElement[NOFIELD]), // NOFIELD
59  boost::make_shared<DerivedEntitiesFieldData>(dataOnElement[H1]), // H1
60  boost::make_shared<DerivedEntitiesFieldData>(
61  dataOnElement[HCURL]), // HCURL
62  boost::make_shared<DerivedEntitiesFieldData>(
63  dataOnElement[HDIV]), // HDIV
64  boost::make_shared<DerivedEntitiesFieldData>(dataOnElement[L2]) // L2
65 
66  },
67  dataNoField(*dataOnElement[NOFIELD].get()),
68  dataH1(*dataOnElement[H1].get()), dataHcurl(*dataOnElement[HCURL].get()),
69  dataHdiv(*dataOnElement[HDIV].get()), dataL2(*dataOnElement[L2].get()),
70  lastEvaluatedElementEntityType(MBMAXTYPE), sidePtrFE(nullptr),
71  refinePtrFE(nullptr) {
72 
73  dataOnElement[NOSPACE]->dataOnEntities[MBENTITYSET].push_back(
75 
76  dataOnElement[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
78  derivedDataOnElement[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
80 }
81 
82 // ** Sense **
83 
85  const EntityType type,
86  boost::ptr_vector<EntitiesFieldData::EntData> &data) const {
88 
89  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable().get<0>();
90  auto sit = side_table.lower_bound(get_id_for_min_type(type));
91  if (sit != side_table.end()) {
92  auto hi_sit = side_table.upper_bound(get_id_for_max_type(type));
93  for (; sit != hi_sit; ++sit) {
94  if (const auto side_number = (*sit)->side_number; side_number >= 0) {
95  const int brother_side_number = (*sit)->brother_side_number;
96  const int sense = (*sit)->sense;
97 
98  data[side_number].getSense() = sense;
99  if (brother_side_number != -1)
100  data[brother_side_number].getSense() = sense;
101  }
102  }
103  }
105 }
106 
107 // ** Order **
108 
109 template <typename ENTMULTIINDEX>
110 static inline int getMaxOrder(const ENTMULTIINDEX &multi_index) {
111  int max_order = 0;
112  for (auto ent_field_weak_ptr : multi_index)
113  if (auto e = ent_field_weak_ptr.lock()) {
114  const int order = e->getMaxOrder();
115  max_order = (max_order < order) ? order : max_order;
116  }
117  return max_order;
118 }
119 
121  int max_order = 0;
122  for (auto e : getDataFieldEnts()) {
123  if (auto ptr = e.lock()) {
124  const int order = ptr->getMaxOrder();
125  max_order = (max_order < order) ? order : max_order;
126  }
127  }
128  return max_order;
129 }
130 
132  return getMaxOrder(getRowFieldEnts());
133 }
134 
136  return getMaxOrder(getColFieldEnts());
137 }
138 
140  const EntityType type, const FieldSpace space,
141  boost::ptr_vector<EntitiesFieldData::EntData> &data) const {
143 
144  auto set_order = [&]() {
146  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable();
147 
148  for (unsigned int s = 0; s != data.size(); ++s)
149  data[s].getOrder() = 0;
150 
151  const FieldEntity_vector_view &data_field_ent = getDataFieldEnts();
152 
153  for (auto r = fieldsPtr->get<BitFieldId_space_mi_tag>().equal_range(space);
154  r.first != r.second; ++r.first) {
155 
156  const auto field_bit_number = (*r.first)->getBitNumber();
157  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
158  field_bit_number, get_id_for_min_type(type));
159  auto lo = std::lower_bound(data_field_ent.begin(), data_field_ent.end(),
160  lo_uid, cmp_uid_lo);
161  if (lo != data_field_ent.end()) {
162  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
163  field_bit_number, get_id_for_max_type(type));
164  auto hi =
165  std::upper_bound(lo, data_field_ent.end(), hi_uid, cmp_uid_hi);
166  for (; lo != hi; ++lo) {
167 
168  if (auto ptr = lo->lock()) {
169 
170  auto &e = *ptr;
171  auto sit = side_table.find(e.getEnt());
172  if (sit != side_table.end()) {
173  auto &side = *sit;
174  if (const auto side_number = side->side_number;
175  side_number >= 0) {
176  ApproximationOrder ent_order = e.getMaxOrder();
177  auto &dat = data[side_number];
178  dat.getOrder() =
179  dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
180  }
181  } else
182  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
183  "Entity on side of the element not found");
184  }
185  }
186  }
187  }
188 
190  };
191 
192  auto set_order_on_brother = [&]() {
194  auto &side_table =
195  numeredEntFiniteElementPtr->getSideNumberTable().get<0>();
196  auto sit = side_table.lower_bound(get_id_for_min_type(type));
197  if (sit != side_table.end()) {
198  auto hi_sit = side_table.upper_bound(get_id_for_max_type(type));
199  for (; sit != hi_sit; ++sit) {
200  const int brother_side_number = (*sit)->brother_side_number;
201  if (brother_side_number != -1) {
202  const int side_number = (*sit)->side_number;
203  data[brother_side_number].getOrder() = data[side_number].getOrder();
204  }
205  }
206  }
208  };
209 
210  CHKERR set_order();
211  CHKERR set_order_on_brother();
212 
214 }
215 
216 // ** Indices **
217 
218 template <typename EXTRACTOR>
220  const int bit_number, FieldEntity_vector_view &ents_field,
221  VectorInt &nodes_indices, VectorInt &local_nodes_indices,
222  EXTRACTOR &&extractor) const {
224 
225  auto field_it = fieldsPtr->get<BitFieldId_mi_tag>().find(
226  BitFieldId().set(bit_number - 1));
227  if (field_it != fieldsPtr->get<BitFieldId_mi_tag>().end()) {
228 
229 #ifndef NDEBUG
230  if ((*field_it)->getBitNumber() != bit_number)
231  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong bit number");
232 #endif
233  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
234  bit_number, get_id_for_min_type<MBVERTEX>());
235  auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
236  cmp_uid_lo);
237  if (lo != ents_field.end()) {
238  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
239  bit_number, get_id_for_max_type<MBVERTEX>());
240  auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
241 
242  const int num_nodes = getNumberOfNodes();
243  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
244  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
245 
246  int nb_dofs = 0;
247  for (auto it = lo; it != hi; ++it) {
248  if (auto e = it->lock()) {
249  if (auto cache = extractor(e).lock()) {
250  if (cache->loHi[0] != cache->loHi[1]) {
251  nb_dofs += std::distance(cache->loHi[0], cache->loHi[1]);
252  break;
253  }
254  }
255  }
256  }
257 
258  if (nb_dofs) {
259  nodes_indices.resize(max_nb_dofs, false);
260  local_nodes_indices.resize(max_nb_dofs, false);
261  } else {
262  nodes_indices.resize(0, false);
263  local_nodes_indices.resize(0, false);
264  }
265 
266  if (nb_dofs != max_nb_dofs) {
267  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
268  std::fill(local_nodes_indices.begin(), local_nodes_indices.end(), -1);
269  }
270 
271  for (auto it = lo; it != hi; ++it) {
272  if (auto e = it->lock()) {
273  auto side_ptr = e->getSideNumberPtr();
274  if (const auto side_number = side_ptr->side_number;
275  side_number >= 0) {
276  const auto brother_side_number = side_ptr->brother_side_number;
277  if (auto cache = extractor(e).lock()) {
278  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
279  auto &dof = **dit;
280  const int idx = dof.getPetscGlobalDofIdx();
281  const int local_idx = dof.getPetscLocalDofIdx();
282  const int pos =
283  side_number * nb_dofs_on_vert + dof.getDofCoeffIdx();
284  nodes_indices[pos] = idx;
285  local_nodes_indices[pos] = local_idx;
286  if (brother_side_number != -1) {
287  const int elem_idx = brother_side_number * nb_dofs_on_vert +
288  (*dit)->getDofCoeffIdx();
289  nodes_indices[elem_idx] = idx;
290  local_nodes_indices[elem_idx] = local_idx;
291  }
292  }
293  }
294  }
295  }
296  }
297  } else {
298  nodes_indices.resize(0, false);
299  local_nodes_indices.resize(0, false);
300  }
301 
302  } else {
303  nodes_indices.resize(0, false);
304  local_nodes_indices.resize(0, false);
305  }
306 
308 }
309 
312  const int bit_number) const {
313 
314  struct Extractor {
315  boost::weak_ptr<EntityCacheNumeredDofs>
316  operator()(boost::shared_ptr<FieldEntity> &e) {
317  return e->entityCacheRowDofs;
318  }
319  };
320 
321  return getNodesIndices(bit_number, getRowFieldEnts(),
322  data.dataOnEntities[MBVERTEX][0].getIndices(),
323  data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
324  Extractor());
325 }
326 
329  const int bit_number) const {
330 
331  struct Extractor {
332  boost::weak_ptr<EntityCacheNumeredDofs>
333  operator()(boost::shared_ptr<FieldEntity> &e) {
334  return e->entityCacheColDofs;
335  }
336  };
337 
338  return getNodesIndices(bit_number, getColFieldEnts(),
339  data.dataOnEntities[MBVERTEX][0].getIndices(),
340  data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
341  Extractor());
342 }
343 
344 template <typename EXTRACTOR>
346  EntitiesFieldData &data, const int bit_number,
347  FieldEntity_vector_view &ents_field, const EntityType type_lo,
348  const EntityType type_hi, EXTRACTOR &&extractor) const {
350 
351  for (EntityType t = type_lo; t != type_hi; ++t) {
352  for (auto &dat : data.dataOnEntities[t]) {
353  dat.getIndices().resize(0, false);
354  dat.getLocalIndices().resize(0, false);
355  }
356  }
357 
358  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
359  bit_number, get_id_for_min_type(type_lo));
360  auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
361  cmp_uid_lo);
362  if (lo != ents_field.end()) {
363  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
364  bit_number, get_id_for_max_type(type_hi));
365  auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
366  if (lo != hi) {
367 
368  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
369 
370  for (auto it = lo; it != hi; ++it)
371  if (auto e = it->lock()) {
372 
373  const EntityType type = e->getEntType();
374  auto side_ptr = e->getSideNumberPtr();
375  if (const auto side = side_ptr->side_number; side >= 0) {
376  const auto nb_dofs_on_ent = e->getNbDofsOnEnt();
377  const auto brother_side = side_ptr->brother_side_number;
378  auto &dat = data.dataOnEntities[type][side];
379  auto &ent_field_indices = dat.getIndices();
380  auto &ent_field_local_indices = dat.getLocalIndices();
381 
382  ent_field_indices.resize(nb_dofs_on_ent, false);
383  ent_field_local_indices.resize(nb_dofs_on_ent, false);
384  std::fill(ent_field_indices.data().begin(),
385  ent_field_indices.data().end(), -1);
386  std::fill(ent_field_local_indices.data().begin(),
387  ent_field_local_indices.data().end(), -1);
388 
389  if (auto cache = extractor(e).lock()) {
390  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
391  const int idx = (*dit)->getEntDofIdx();
392  ent_field_indices[idx] = (*dit)->getPetscGlobalDofIdx();
393  ent_field_local_indices[idx] = (*dit)->getPetscLocalDofIdx();
394  }
395  }
396 
397  if (brother_side != -1) {
398  auto &dat_brother = data.dataOnEntities[type][brother_side];
399  dat_brother.getIndices().resize(nb_dofs_on_ent, false);
400  dat_brother.getLocalIndices().resize(nb_dofs_on_ent, false);
401  noalias(dat_brother.getIndices()) = dat.getIndices();
402  noalias(dat_brother.getLocalIndices()) = dat.getLocalIndices();
403  }
404  }
405  }
406  }
407  }
408 
410 }
411 
413  EntitiesFieldData &data, const int bit_number, const EntityType type_lo,
414  const EntityType type_hi) const {
415 
416  struct Extractor {
417  boost::weak_ptr<EntityCacheNumeredDofs>
418  operator()(boost::shared_ptr<FieldEntity> &e) {
419  return e->entityCacheRowDofs;
420  }
421  };
422 
423  return getEntityIndices(data, bit_number, getRowFieldEnts(), type_lo, type_hi,
424  Extractor());
425 }
426 
428  EntitiesFieldData &data, const int bit_number, const EntityType type_lo,
429  const EntityType type_hi) const {
430 
431  struct Extractor {
432  boost::weak_ptr<EntityCacheNumeredDofs>
433  operator()(boost::shared_ptr<FieldEntity> &e) {
434  return e->entityCacheColDofs;
435  }
436  };
437 
438  return getEntityIndices(data, bit_number, getColFieldEnts(), type_lo, type_hi,
439  Extractor());
440 }
441 
443  const int bit_number, boost::shared_ptr<FENumeredDofEntity_multiIndex> dofs,
444  VectorInt &indices) const {
446  auto dit = dofs->get<Unique_mi_tag>().lower_bound(
447  FieldEntity::getLoBitNumberUId(bit_number));
448  auto hi_dit = dofs->get<Unique_mi_tag>().upper_bound(
449  FieldEntity::getHiBitNumberUId(bit_number));
450  indices.resize(std::distance(dit, hi_dit));
451  for (; dit != hi_dit; dit++) {
452  int idx = (*dit)->getPetscGlobalDofIdx();
453  indices[(*dit)->getDofCoeffIdx()] = idx;
454  }
456 }
457 
460  const int bit_number) const {
462 #ifndef NDEBUG
463  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
464  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
465  "data.dataOnEntities[MBENTITYSET] is empty");
466  }
467 #endif
468  CHKERR getNoFieldIndices(bit_number, getRowDofsPtr(),
469  data.dataOnEntities[MBENTITYSET][0].getIndices());
471 }
472 
475  const int bit_number) const {
477 #ifndef NDEBUG
478  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
479  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
480  "data.dataOnEntities[MBENTITYSET] is empty");
481  }
482 #endif
483  CHKERR getNoFieldIndices(bit_number, getColDofsPtr(),
484  data.dataOnEntities[MBENTITYSET][0].getIndices());
486 }
487 
488 // ** Indices from problem **
489 
491  const std::string &field_name, const NumeredDofEntity_multiIndex &dofs,
492  VectorInt &nodes_indices) const {
494 
495  const Field *field_struture = mField.get_field_structure(field_name);
496  if (field_struture->getSpace() == H1) {
497 
498  const int num_nodes = getNumberOfNodes();
499  nodes_indices.resize(field_struture->getNbOfCoeffs() * num_nodes, false);
500  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
501 
502  auto &side_table = const_cast<SideNumber_multiIndex &>(
503  numeredEntFiniteElementPtr->getSideNumberTable());
504  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 0));
505  auto hi_siit =
506  side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 10000));
507 
508  int nn = 0;
509  for (; siit != hi_siit; siit++, nn++) {
510  if (siit->get()->side_number >= 0) {
511  auto bit_number = mField.get_field_bit_number(field_name);
512  const EntityHandle ent = siit->get()->ent;
513  auto dit = dofs.get<Unique_mi_tag>().lower_bound(
514  FieldEntity::getLoLocalEntityBitNumber(bit_number, ent));
515  auto hi_dit = dofs.get<Unique_mi_tag>().upper_bound(
516  FieldEntity::getHiLocalEntityBitNumber(bit_number, ent));
517  for (; dit != hi_dit; dit++) {
518  nodes_indices[siit->get()->side_number * (*dit)->getNbOfCoeffs() +
519  (*dit)->getDofCoeffIdx()] =
520  (*dit)->getPetscGlobalDofIdx();
521  }
522  }
523  }
524  } else {
525  nodes_indices.resize(0, false);
526  }
527 
529 }
530 
532  const std::string &field_name, const NumeredDofEntity_multiIndex &dofs,
533  EntityType type, int side_number, VectorInt &indices) const {
535 
536  indices.resize(0);
537 
538  auto &side_table = const_cast<SideNumber_multiIndex &>(
539  numeredEntFiniteElementPtr->getSideNumberTable());
540  auto siit =
541  side_table.get<1>().lower_bound(boost::make_tuple(type, side_number));
542  auto hi_siit =
543  side_table.get<1>().upper_bound(boost::make_tuple(type, side_number));
544 
545  for (; siit != hi_siit; siit++) {
546  if (siit->get()->side_number >= 0) {
547 
548  const EntityHandle ent = siit->get()->ent;
549  auto bit_number = mField.get_field_bit_number(field_name);
550  auto dit = dofs.get<Unique_mi_tag>().lower_bound(
551  FieldEntity::getLoLocalEntityBitNumber(bit_number, ent));
552  auto hi_dit = dofs.get<Unique_mi_tag>().upper_bound(
553  FieldEntity::getHiLocalEntityBitNumber(bit_number, ent));
554  indices.resize(std::distance(dit, hi_dit));
555  for (; dit != hi_dit; dit++) {
556 
557  indices[(*dit)->getEntDofIdx()] = (*dit)->getPetscGlobalDofIdx();
558  }
559  }
560  }
561 
563 }
564 
566  const std::string &field_name, VectorInt &nodes_indices) const {
568  nodes_indices);
569 }
570 
573  EntityType type, int side_number,
574  VectorInt &indices) const {
576  type, side_number, indices);
577 }
578 
580  const std::string &field_name, VectorInt &nodes_indices) const {
582  nodes_indices);
583 }
584 
587  EntityType type, int side_number,
588  VectorInt &indices) const {
590  type, side_number, indices);
591 }
592 
593 // ** Data **
594 
597 
598  for (auto &data : dataOnElement) {
599  if (data) {
600  for (auto &dat : data->dataOnEntities) {
601  for (auto &ent_dat : dat) {
602  ent_dat.getEntDataBitRefLevel().clear();
603  }
604  }
605  }
606  }
607 
608  auto &field_ents = getDataFieldEnts();
609  for (auto it : field_ents) {
610  if (auto e = it.lock()) {
611  const FieldSpace space = e->getSpace();
612  if (space > NOFIELD) {
613  const EntityType type = e->getEntType();
614  const signed char side = e->getSideNumberPtr()->side_number;
615  if (side >= 0) {
616  if (auto &data = dataOnElement[space]) {
617  if (type == MBVERTEX) {
618  auto &dat = data->dataOnEntities[type][0];
619  dat.getEntDataBitRefLevel().resize(getNumberOfNodes(), false);
620  dat.getEntDataBitRefLevel()[side] = e->getBitRefLevel();
621  } else {
622  auto &dat = data->dataOnEntities[type][side];
623  dat.getEntDataBitRefLevel().resize(1, false);
624  dat.getEntDataBitRefLevel()[0] = e->getBitRefLevel();
625  }
626  }
627  }
628  } else {
629  if (auto &data = dataOnElement[NOFIELD]) {
630  auto &dat = data->dataOnEntities[MBENTITYSET][0];
631  dat.getEntDataBitRefLevel().resize(1, false);
632  dat.getEntDataBitRefLevel()[0] = e->getBitRefLevel();
633  }
634  }
635  }
636  }
637 
639 };
640 
643  const int bit_number) const {
644 
645  auto get_nodes_field_data = [&](VectorDouble &nodes_data,
646  VectorFieldEntities &field_entities,
647  VectorDofs &nodes_dofs, FieldSpace &space,
649  VectorInt &bb_node_order) {
651 
652  nodes_data.resize(0, false);
653  nodes_dofs.resize(0, false);
654  field_entities.resize(0, false);
655 
656  auto field_it =
657  fieldsPtr->get<BitFieldId_mi_tag>().find(BitFEId().set(bit_number - 1));
658  if (field_it != fieldsPtr->get<BitFieldId_mi_tag>().end()) {
659 
660 #ifndef NDEBUG
661  if ((*field_it)->getBitNumber() != bit_number)
662  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong bit number");
663 #endif
664  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
665  space = (*field_it)->getSpace();
666  base = (*field_it)->getApproxBase();
667 
668  auto &field_ents = getDataFieldEnts();
669  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
670  bit_number, get_id_for_min_type<MBVERTEX>());
671  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
672  cmp_uid_lo);
673  if (lo != field_ents.end()) {
674  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
675  bit_number, get_id_for_max_type<MBVERTEX>());
676  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
677  if (lo != hi) {
678 
679  int nb_dofs = 0;
680  for (auto it = lo; it != hi; ++it) {
681  if (auto e = it->lock()) {
682  if (auto cache = e->entityCacheDataDofs.lock()) {
683  if (cache->loHi[0] != cache->loHi[1]) {
684  nb_dofs += std::distance(cache->loHi[0], cache->loHi[1]);
685  break;
686  }
687  }
688  }
689  }
690 
691  if (nb_dofs) {
692 
693  const int num_nodes = getNumberOfNodes();
694  bb_node_order.resize(num_nodes, false);
695  bb_node_order.clear();
696  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
697  nodes_data.resize(max_nb_dofs, false);
698  nodes_dofs.resize(max_nb_dofs, false);
699  field_entities.resize(num_nodes, false);
700  std::fill(nodes_data.begin(), nodes_data.end(), 0);
701  std::fill(nodes_dofs.begin(), nodes_dofs.end(), nullptr);
702  std::fill(field_entities.begin(), field_entities.end(), nullptr);
703 
704  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
705 
706  for (auto it = lo; it != hi; ++it) {
707  if (auto e = it->lock()) {
708  const auto &sn = e->getSideNumberPtr();
709  // Some field entities on skeleton can have negative side
710  // number
711  if (const auto side_number = sn->side_number;
712  side_number >= 0) {
713  const int brother_side_number = sn->brother_side_number;
714 
715  field_entities[side_number] = e.get();
716  if (brother_side_number != -1) {
717  brother_ents_vec.emplace_back(e);
718  field_entities[side_number] = field_entities[side_number];
719  }
720 
721  bb_node_order[side_number] = e->getMaxOrder();
722  int pos = side_number * nb_dofs_on_vert;
723  auto ent_filed_data_vec = e->getEntFieldData();
724  if (auto cache = e->entityCacheDataDofs.lock()) {
725  for (auto dit = cache->loHi[0]; dit != cache->loHi[1];
726  ++dit) {
727  const auto dof_idx = (*dit)->getEntDofIdx();
728  nodes_data[pos + dof_idx] = ent_filed_data_vec[dof_idx];
729  nodes_dofs[pos + dof_idx] =
730  reinterpret_cast<FEDofEntity *>((*dit).get());
731  }
732  }
733  }
734  }
735  }
736 
737  for (auto &it : brother_ents_vec) {
738  if (const auto e = it.lock()) {
739  const auto &sn = e->getSideNumberPtr();
740  const int side_number = sn->side_number;
741  const int brother_side_number = sn->brother_side_number;
742  bb_node_order[brother_side_number] = bb_node_order[side_number];
743  int pos = side_number * nb_dofs_on_vert;
744  int brother_pos = brother_side_number * nb_dofs_on_vert;
745  for (int ii = 0; ii != nb_dofs_on_vert; ++ii) {
746  nodes_data[brother_pos] = nodes_data[pos];
747  nodes_dofs[brother_pos] = nodes_dofs[pos];
748  ++pos;
749  ++brother_pos;
750  }
751  }
752  }
753  }
754  }
755  }
756  }
757 
759  };
760 
761  return get_nodes_field_data(
762  data.dataOnEntities[MBVERTEX][0].getFieldData(),
763  data.dataOnEntities[MBVERTEX][0].getFieldEntities(),
764  data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
765  data.dataOnEntities[MBVERTEX][0].getSpace(),
766  data.dataOnEntities[MBVERTEX][0].getBase(),
767  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder());
768 }
769 
771  EntitiesFieldData &data, const int bit_number, const EntityType type_lo,
772  const EntityType type_hi) const {
774  for (EntityType t = type_lo; t != type_hi; ++t) {
775  for (auto &dat : data.dataOnEntities[t]) {
776  dat.getOrder() = 0;
777  dat.getBase() = NOBASE;
778  dat.getSpace() = NOSPACE;
779  dat.getFieldData().resize(0, false);
780  dat.getFieldDofs().resize(0, false);
781  dat.getFieldEntities().resize(0, false);
782  }
783  }
784 
785  auto &field_ents = getDataFieldEnts();
786  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
787  bit_number, get_id_for_min_type(type_lo));
788  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
789  cmp_uid_lo);
790  if (lo != field_ents.end()) {
791  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
792  bit_number, get_id_for_max_type(type_hi));
793  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
794  if (lo != hi) {
795 
796  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
797 
798  for (auto it = lo; it != hi; ++it)
799  if (auto e = it->lock()) {
800  auto side_ptr = e->getSideNumberPtr();
801  if (const auto side = side_ptr->side_number; side >= 0) {
802  const EntityType type = e->getEntType();
803  auto &dat = data.dataOnEntities[type][side];
804  auto &ent_field = dat.getFieldEntities();
805  auto &ent_field_dofs = dat.getFieldDofs();
806  auto &ent_field_data = dat.getFieldData();
807 
808  const int brother_side = side_ptr->brother_side_number;
809  if (brother_side != -1)
810  brother_ents_vec.emplace_back(e);
811 
812  dat.getBase() = e->getApproxBase();
813  dat.getSpace() = e->getSpace();
814  const int ent_order = e->getMaxOrder();
815  dat.getOrder() =
816  dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
817 
818  auto ent_data = e->getEntFieldData();
819  ent_field_data.resize(ent_data.size(), false);
820  noalias(ent_field_data) = ent_data;
821  ent_field_dofs.resize(ent_data.size(), false);
822  std::fill(ent_field_dofs.begin(), ent_field_dofs.end(), nullptr);
823  ent_field.resize(1, false);
824  ent_field[0] = e.get();
825  if (auto cache = e->entityCacheDataDofs.lock()) {
826  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
827  ent_field_dofs[(*dit)->getEntDofIdx()] =
828  reinterpret_cast<FEDofEntity *>((*dit).get());
829  }
830  }
831  }
832  }
833 
834  for (auto &it : brother_ents_vec) {
835  if (const auto e = it.lock()) {
836  const EntityType type = e->getEntType();
837  const int side = e->getSideNumberPtr()->side_number;
838  const int brother_side = e->getSideNumberPtr()->brother_side_number;
839  auto &dat = data.dataOnEntities[type][side];
840  auto &dat_brother = data.dataOnEntities[type][brother_side];
841  dat_brother.getBase() = dat.getBase();
842  dat_brother.getSpace() = dat.getSpace();
843  dat_brother.getOrder() = dat.getOrder();
844  dat_brother.getFieldData() = dat.getFieldData();
845  dat_brother.getFieldDofs() = dat.getFieldDofs();
846  dat_brother.getFieldEntities() = dat.getFieldEntities();
847  }
848  }
849  }
850  }
851 
853 }
854 
856  const int bit_number, VectorDouble &ent_field_data,
857  VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const {
858 
860 
861  ent_field_data.resize(0, false);
862  ent_field_dofs.resize(0, false);
863  ent_field.resize(0, false);
864 
865  auto &field_ents = getDataFieldEnts();
866  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
867  bit_number, get_id_for_min_type<MBVERTEX>());
868  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
869  cmp_uid_lo);
870  if (lo != field_ents.end()) {
871 
872  ent_field.resize(field_ents.size(), false);
873  std::fill(ent_field.begin(), ent_field.end(), nullptr);
874 
875  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
876  bit_number, get_id_for_max_type<MBENTITYSET>());
877  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
878  if (lo != hi) {
879 
880  int side = 0;
881  for (auto it = lo; it != hi; ++it, ++side)
882  if (auto e = it->lock()) {
883 
884  const auto size = e->getNbDofsOnEnt();
885  ent_field_data.resize(size, false);
886  ent_field_dofs.resize(size, false);
887  ent_field[side] = e.get();
888  noalias(ent_field_data) = e->getEntFieldData();
889 
890  if (auto cache = e->entityCacheDataDofs.lock()) {
891  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
892  ent_field_dofs[(*dit)->getEntDofIdx()] =
893  reinterpret_cast<FEDofEntity *>((*dit).get());
894  }
895  }
896  }
897  }
898  }
899 
901 }
902 
905  const int bit_number) const {
907  if (data.dataOnEntities[MBENTITYSET].size() == 0)
908  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
909  "No space to insert data");
910 
912  bit_number, data.dataOnEntities[MBENTITYSET][0].getFieldData(),
913  data.dataOnEntities[MBENTITYSET][0].getFieldDofs(),
914  data.dataOnEntities[MBENTITYSET][0].getFieldEntities());
916 }
917 
918 // ** Face **
919 
923  auto &face_nodes = data.facesNodes;
924  auto &face_nodes_order = data.facesNodesOrder;
925  auto &side_table = const_cast<SideNumber_multiIndex &>(
926  numeredEntFiniteElementPtr->getSideNumberTable());
927  const auto ent = numeredEntFiniteElementPtr->getEnt();
928  const auto type = numeredEntFiniteElementPtr->getEntType();
929  const auto nb_faces = CN::NumSubEntities(type, 2);
930  const EntityHandle *conn_ele;
931  int num_nodes_ele;
932  CHKERR mField.get_moab().get_connectivity(ent, conn_ele, num_nodes_ele, true);
933  auto side_ptr_it = side_table.get<1>().lower_bound(
934  boost::make_tuple(CN::TypeDimensionMap[2].first, 0));
935  auto hi_side_ptr_it = side_table.get<1>().upper_bound(
936  boost::make_tuple(CN::TypeDimensionMap[2].second, 100));
937 
938  for (; side_ptr_it != hi_side_ptr_it; ++side_ptr_it) {
939  const auto side = (*side_ptr_it)->side_number;
940  const auto sense = (*side_ptr_it)->sense;
941  const auto offset = (*side_ptr_it)->offset;
942 
943  EntityType face_type;
944  int nb_nodes_face;
945  auto face_indices =
946  CN::SubEntityVertexIndices(type, 2, side, face_type, nb_nodes_face);
947  face_nodes.resize(nb_faces, nb_nodes_face);
948  face_nodes_order.resize(nb_faces, nb_nodes_face);
949 
950  if (sense == 1)
951  for (int n = 0; n != nb_nodes_face; ++n)
952  face_nodes_order(side, n) = (n + offset) % nb_nodes_face;
953  else
954  for (int n = 0; n != nb_nodes_face; ++n)
955  face_nodes_order(side, n) =
956  (nb_nodes_face - (n - offset) % nb_nodes_face) % nb_nodes_face;
957 
958  for (int n = 0; n != nb_nodes_face; ++n)
959  face_nodes(side, n) = face_indices[face_nodes_order(side, n)];
960 
961 #ifndef NDEBUG
962  const auto face_ent = (*side_ptr_it)->ent;
963  auto check = [&]() {
965  const EntityHandle *conn_face;
966  // int nb_nodes_face;
967  CHKERR mField.get_moab().get_connectivity(face_ent, conn_face,
968  nb_nodes_face, true);
969  face_nodes.resize(nb_faces, nb_nodes_face);
970  for (int nn = 0; nn != nb_nodes_face; ++nn) {
971  if (face_nodes(side, nn) !=
972  std::distance(
973  conn_ele,
974  std::find(conn_ele, &conn_ele[num_nodes_ele], conn_face[nn]))) {
975  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
976  "Wrong face numeration");
977  }
978  }
980  };
981  CHKERR check();
982 #endif
983  }
984 
986 }
987 
988 // ** Space and Base **
989 
991  EntitiesFieldData &data) const {
993 
994  if (nInTheLoop == 0) {
995  data.sPace.reset();
996  data.bAse.reset();
997  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
998  data.spacesOnEntities[t].reset();
999  data.basesOnEntities[t].reset();
1000  }
1001  for (int s = 0; s != LASTSPACE; ++s) {
1002  data.basesOnSpaces[s].reset();
1003  }
1004  }
1005 
1006  auto fe_type = numeredEntFiniteElementPtr->getEntType();
1007 
1008  if (getDataFieldEntsPtr())
1009  for (auto e : getDataFieldEnts()) {
1010  if (auto ptr = e.lock()) {
1011  // get data from entity
1012  const EntityType type = ptr->getEntType();
1014  const FieldSpace space = ptr->getSpace();
1015  const FieldApproximationBase approx = ptr->getApproxBase();
1016  // set data
1017  data.sPace.set(space);
1018  data.bAse.set(approx);
1019  data.spacesOnEntities[type].set(space);
1020  data.basesOnEntities[type].set(approx);
1021  data.basesOnSpaces[space].set(approx);
1022  }
1023  }
1024  }
1025  else
1026  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1027  "data fields ents not allocated on element");
1028 
1030 }
1031 
1033  const FieldApproximationBase b) {
1035  if (dataOnElement[H1]->bAse.test(b)) {
1036  switch (static_cast<FieldApproximationBase>(b)) {
1037  case NOBASE:
1038  break;
1040  // BERNSTEIN_BEZIER_BASE is not hierarchical base
1041  break;
1044  case DEMKOWICZ_JACOBI_BASE:
1045  if (!getElementPolynomialBase())
1046  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1047  "Functions genrating approximation base not defined");
1048 
1049  for (int space = H1; space != LASTSPACE; ++space) {
1050  if (dataOnElement[H1]->sPace.test(space) &&
1051  dataOnElement[H1]->bAse.test(b) &&
1052  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1053  CHKERR getElementPolynomialBase()->getValue(
1054  gaussPts,
1055  boost::make_shared<EntPolynomialBaseCtx>(
1056  *dataOnElement[space], static_cast<FieldSpace>(space),
1057  static_cast<FieldApproximationBase>(b), NOBASE));
1058  }
1059  }
1060  break;
1061  case USER_BASE:
1062  if (!getUserPolynomialBase())
1063  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1064  "Functions generating user approximation base not defined");
1065 
1066  for (int space = H1; space != LASTSPACE; ++space)
1067  if (dataOnElement[H1]->sPace.test(space) &&
1068  dataOnElement[H1]->bAse.test(b) &&
1069  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1070  CHKERR getUserPolynomialBase()->getValue(
1071  gaussPts,
1072  boost::make_shared<EntPolynomialBaseCtx>(
1073  *dataOnElement[space], static_cast<FieldSpace>(space),
1074  static_cast<FieldApproximationBase>(b), NOBASE));
1075  }
1076  break;
1077  default:
1079  "Base <%s> not yet implemented",
1081  }
1082  }
1084 }
1085 
1088  /// Use the some node base. Node base is usually used for construction other
1089  /// bases.
1090  for (int space = HCURL; space != LASTSPACE; ++space) {
1091  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
1092  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
1093  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(
1094  NOBASE) =
1095  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(
1096  NOBASE);
1097  }
1098  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
1100  static_cast<FieldApproximationBase>(b));
1101  }
1103 }
1104 
1108 
1109  const auto ele_type = numeredEntFiniteElementPtr->getEntType();
1110 
1111  auto get_nodal_base_data = [&](EntitiesFieldData &data, auto field_ptr) {
1113  auto &space = data.dataOnEntities[MBVERTEX][0].getSpace();
1114  auto &base = data.dataOnEntities[MBVERTEX][0].getBase();
1115  auto &bb_node_order = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder();
1116 
1117  auto &field_ents = getDataFieldEnts();
1118  auto bit_number = field_ptr->getBitNumber();
1119  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1120  bit_number, get_id_for_min_type<MBVERTEX>());
1121  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
1122  cmp_uid_lo);
1123  if (lo != field_ents.end()) {
1124  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1125  bit_number, get_id_for_max_type<MBVERTEX>());
1126  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
1127  if (lo != hi) {
1128 
1129  for (auto it = lo; it != hi; ++it)
1130  if (auto first_e = it->lock()) {
1131  space = first_e->getSpace();
1132  base = first_e->getApproxBase();
1133  const int num_nodes = getNumberOfNodes();
1134  bb_node_order.resize(num_nodes, false);
1135  bb_node_order.clear();
1136 
1137  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1138 
1139  for (; it != hi; ++it) {
1140  if (auto e = it->lock()) {
1141  const auto &sn = e->getSideNumberPtr();
1142  const int side_number = sn->side_number;
1143  const int brother_side_number = sn->brother_side_number;
1144  if (brother_side_number != -1)
1145  brother_ents_vec.emplace_back(e);
1146  bb_node_order[side_number] = e->getMaxOrder();
1147  }
1148  }
1149 
1150  for (auto &it : brother_ents_vec) {
1151  if (const auto e = it.lock()) {
1152  const auto &sn = e->getSideNumberPtr();
1153  const int side_number = sn->side_number;
1154  const int brother_side_number = sn->brother_side_number;
1155  bb_node_order[brother_side_number] = bb_node_order[side_number];
1156  }
1157  }
1158 
1159  break;
1160  }
1161  }
1162  }
1163 
1165  };
1166 
1167  auto get_entity_base_data = [&](EntitiesFieldData &data, auto field_ptr,
1168  const EntityType type_lo,
1169  const EntityType type_hi) {
1171  for (EntityType t = MBEDGE; t != MBPOLYHEDRON; ++t) {
1172  for (auto &dat : data.dataOnEntities[t]) {
1173  dat.getOrder() = 0;
1174  dat.getBase() = NOBASE;
1175  dat.getSpace() = NOSPACE;
1176  dat.getFieldData().resize(0, false);
1177  dat.getFieldDofs().resize(0, false);
1178  }
1179  }
1180 
1181  auto &field_ents = getDataFieldEnts();
1182  auto bit_number = field_ptr->getBitNumber();
1183  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1184  bit_number, get_id_for_min_type(type_lo));
1185  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
1186  cmp_uid_lo);
1187  if (lo != field_ents.end()) {
1188  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1189  bit_number, get_id_for_max_type(type_hi));
1190  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
1191 
1192  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1193  for (; lo != hi; ++lo) {
1194  if (auto e = lo->lock()) {
1195  if (auto cache = e->entityCacheDataDofs.lock()) {
1196  if (cache->loHi[0] != cache->loHi[1]) {
1197  if (const auto side = e->getSideNumberPtr()->side_number;
1198  side >= 0) {
1199  const EntityType type = e->getEntType();
1200  auto &dat = data.dataOnEntities[type][side];
1201  const int brother_side =
1202  e->getSideNumberPtr()->brother_side_number;
1203  if (brother_side != -1)
1204  brother_ents_vec.emplace_back(e);
1205  dat.getBase() = e->getApproxBase();
1206  dat.getSpace() = e->getSpace();
1207  const auto ent_order = e->getMaxOrder();
1208  dat.getOrder() =
1209  dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
1210  }
1211  }
1212  }
1213  }
1214  }
1215 
1216  for (auto &ent_ptr : brother_ents_vec) {
1217  if (auto e = ent_ptr.lock()) {
1218  const EntityType type = e->getEntType();
1219  const int side = e->getSideNumberPtr()->side_number;
1220  const int brother_side = e->getSideNumberPtr()->brother_side_number;
1221  auto &dat = data.dataOnEntities[type][side];
1222  auto &dat_brother = data.dataOnEntities[type][brother_side];
1223  dat_brother.getBase() = dat.getBase();
1224  dat_brother.getSpace() = dat.getSpace();
1225  dat_brother.getOrder() = dat.getOrder();
1226  }
1227  }
1228  }
1230  };
1231 
1232  for (auto &e : getDataFieldEnts()) {
1233  if (auto ent_data_ptr = e.lock()) {
1234  if (ent_data_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1235  auto space = ent_data_ptr->getSpace();
1236  for (EntityType t = MBVERTEX; t != MBPOLYHEDRON; ++t) {
1237  for (auto &dat : (*dataOnElement[space]).dataOnEntities[t]) {
1238  for (auto &ptr : dat.getBBAlphaIndicesByOrderArray())
1239  ptr.reset();
1240  for (auto &ptr : dat.getBBNByOrderArray())
1241  ptr.reset();
1242  for (auto &ptr : dat.getBBDiffNByOrderArray())
1243  ptr.reset();
1244  }
1245  }
1246  }
1247  }
1248  }
1249 
1250  auto check_space = [&](const auto space) {
1251  switch (space) {
1252  case H1:
1253  for (auto t = MBVERTEX; t <= ele_type; ++t) {
1254  if (dataOnElement[H1]->spacesOnEntities[t].test(H1))
1255  return true;
1256  }
1257  return false;
1258  case HCURL:
1259  for (auto t = MBEDGE; t <= ele_type; ++t) {
1260  if (dataOnElement[HCURL]->spacesOnEntities[t].test(HCURL))
1261  return true;
1262  }
1263  return false;
1264  case HDIV:
1265  for (auto t = MBTRI; t <= ele_type; ++t) {
1266  if (dataOnElement[HDIV]->spacesOnEntities[t].test(HDIV))
1267  return true;
1268  }
1269  return false;
1270  case L2:
1271  return dataOnElement[L2]->spacesOnEntities[ele_type].test(L2);
1272  break;
1273  default:
1274  THROW_MESSAGE("Not implemented");
1275  }
1276  };
1277 
1278  std::set<string> fields_list;
1279  for (auto &e : getDataFieldEnts()) {
1280  if (auto ent_data_ptr = e.lock()) {
1281  if (ent_data_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1282  auto field_name = ent_data_ptr->getName();
1283  if (fields_list.find(field_name) == fields_list.end()) {
1284  auto field_ptr = ent_data_ptr->getFieldRawPtr();
1285  auto space = ent_data_ptr->getSpace();
1286  CHKERR get_nodal_base_data(*dataOnElement[space], field_ptr);
1287  CHKERR get_entity_base_data(*dataOnElement[space], field_ptr, MBEDGE,
1288  ele_type);
1289  if (check_space(space)) {
1290  CHKERR getElementPolynomialBase()->getValue(
1291  gaussPts, boost::make_shared<EntPolynomialBaseCtx>(
1292  *dataOnElement[space], field_name,
1293  static_cast<FieldSpace>(space),
1295  fields_list.insert(field_name);
1296  }
1297  }
1298  }
1299  }
1300  }
1301 
1303 };
1304 
1307 
1308  // Data on elements for proper spaces
1309  for (int space = H1; space != LASTSPACE; ++space) {
1310  dataOnElement[space]->setElementType(type);
1311  derivedDataOnElement[space]->setElementType(type);
1312  }
1313 
1315 }
1316 
1317 #define FUNCTION_NAME_WITH_OP_NAME(OP) \
1318  std::ostringstream ss; \
1319  ss << "(Calling user data operator " \
1320  << boost::typeindex::type_id_runtime(OP).pretty_name() << " rowField " \
1321  << (OP).rowFieldName << " colField " << (OP).colFieldName << ") "
1322 
1323 #define CATCH_OP_ERRORS(OP) \
1324  catch (MoFEMExceptionInitial const &ex) { \
1325  FUNCTION_NAME_WITH_OP_NAME(OP) << PETSC_FUNCTION_NAME; \
1326  return PetscError(PETSC_COMM_SELF, ex.lINE, ss.str().c_str(), __FILE__, \
1327  ex.errorCode, PETSC_ERROR_INITIAL, ex.what()); \
1328  } \
1329  catch (MoFEMExceptionRepeat const &ex) { \
1330  FUNCTION_NAME_WITH_OP_NAME(OP) << PETSC_FUNCTION_NAME; \
1331  return PetscError(PETSC_COMM_SELF, ex.lINE, ss.str().c_str(), __FILE__, \
1332  ex.errorCode, PETSC_ERROR_REPEAT, " "); \
1333  } \
1334  catch (MoFEMException const &ex) { \
1335  FUNCTION_NAME_WITH_OP_NAME(OP) << ex.errorMessage; \
1336  SETERRQ(PETSC_COMM_SELF, ex.errorCode, ss.str().c_str()); \
1337  } \
1338  catch (std::exception const &ex) { \
1339  std::string message("Error: " + std::string(ex.what()) + " at " + \
1340  boost::lexical_cast<std::string>(__LINE__) + " : " + \
1341  std::string(__FILE__) + " in " + \
1342  std::string(PETSC_FUNCTION_NAME)); \
1343  FUNCTION_NAME_WITH_OP_NAME(OP) << message; \
1344  SETERRQ(PETSC_COMM_SELF, MOFEM_STD_EXCEPTION_THROW, ss.str().c_str()); \
1345  }
1346 
1349 
1350  using UDO = UserDataOperator;
1351 
1352  std::array<std::string, 2> field_name;
1353  std::array<const Field *, 3> field_struture;
1354  std::array<int, 2>
1355  field_id; // note the this is field bit number (nor field bit)
1356  std::array<FieldSpace, 2> space;
1357  std::array<FieldApproximationBase, 2> base;
1358 
1359  constexpr std::array<UDO::OpType, 2> types = {UDO::OPROW, UDO::OPCOL};
1360  std::array<int, 2> last_eval_field_id = {0, 0};
1361 
1362  std::array<boost::shared_ptr<EntitiesFieldData>, 2> op_data;
1363 
1364  auto swap_bases = [&](auto &op) {
1366  for (size_t ss = 0; ss != 2; ++ss) {
1367  if (op.getOpType() & types[ss] || op.getOpType() & UDO::OPROWCOL) {
1368  switch (base[ss]) {
1370  CHKERR op_data[ss]->baseSwap(field_name[ss],
1372  default:
1373  break;
1374  }
1375  }
1376  }
1377 
1379  };
1380 
1381  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1382 
1383  // evaluate entity data only, no field specific data provided or known
1384  auto evaluate_op_space = [&](auto &op) {
1386 
1387  // reseat all data which all field dependent
1388  dataOnElement[op.sPace]->resetFieldDependentData();
1389  std::fill(last_eval_field_id.begin(), last_eval_field_id.end(), 0);
1390 
1391  switch (op.sPace) {
1392  case NOSPACE:
1393  try {
1394  CHKERR op.doWork(
1395  0, MBENTITYSET,
1396  dataOnElement[NOSPACE]->dataOnEntities[MBENTITYSET][0]);
1397  }
1398  CATCH_OP_ERRORS(op);
1399  break;
1400  case NOFIELD:
1401  case H1:
1402  case HCURL:
1403  case HDIV:
1404  case L2:
1405  try {
1406 
1407  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
1408  for (auto &e : dataOnElement[op.sPace]->dataOnEntities[t]) {
1409  e.getSpace() = op.sPace;
1410  }
1411  }
1412 
1413  CHKERR op.opRhs(*dataOnElement[op.sPace], false);
1414  }
1415  CATCH_OP_ERRORS(op);
1416  break;
1417  default:
1418  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1419  "Not implemented for this space", op.sPace);
1420  }
1421 
1423  };
1424 
1425  // set entity data
1426  auto set_op_entities_data = [&](auto ss, auto &op) {
1428 
1429 #ifndef NDEBUG
1430  if ((op.getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
1432  .none()) {
1433  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1434  "no data field < %s > on finite element < %s >",
1435  field_name[ss].c_str(), getFEName().c_str());
1436  }
1437 #endif
1438 
1439  op_data[ss] =
1440  !ss ? dataOnElement[space[ss]] : derivedDataOnElement[space[ss]];
1441 
1442  for (auto &data : op_data[ss]->dataOnEntities[MBENTITYSET]) {
1443  CHKERR data.resetFieldDependentData();
1444  }
1445 
1446  auto get_data_for_nodes = [&]() {
1448  if (!ss)
1449  CHKERR getRowNodesIndices(*op_data[ss], field_id[ss]);
1450  else
1451  CHKERR getColNodesIndices(*op_data[ss], field_id[ss]);
1452  CHKERR getNodesFieldData(*op_data[ss], field_id[ss]);
1454  };
1455 
1456  // get data on entities but not nodes
1457  auto get_data_for_entities = [&]() {
1459  CHKERR getEntityFieldData(*op_data[ss], field_id[ss], MBEDGE);
1460  if (!ss)
1461  CHKERR getEntityRowIndices(*op_data[ss], field_id[ss], MBEDGE);
1462  else
1463  CHKERR getEntityColIndices(*op_data[ss], field_id[ss], MBEDGE);
1465  };
1466 
1467  auto get_data_for_meshset = [&]() {
1469  if (!ss) {
1470  CHKERR getNoFieldRowIndices(*op_data[ss], field_id[ss]);
1471  } else {
1472  CHKERR getNoFieldColIndices(*op_data[ss], field_id[ss]);
1473  }
1474  CHKERR getNoFieldFieldData(*op_data[ss], field_id[ss]);
1476  };
1477 
1478  switch (space[ss]) {
1479  case H1:
1480  CHKERR get_data_for_nodes();
1481  case HCURL:
1482  case HDIV:
1483  CHKERR get_data_for_entities();
1484  break;
1485  case L2:
1486  switch (type) {
1487  case MBVERTEX:
1488  CHKERR get_data_for_nodes();
1489  break;
1490  default:
1491  CHKERR get_data_for_entities();
1492  }
1493  break;
1494  case NOFIELD:
1495  // if (!getNinTheLoop()) {
1496  // NOFIELD data are the same for each element, can be
1497  // retrieved only once
1498  CHKERR get_data_for_meshset();
1499  // }
1500  break;
1501  default:
1502  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1503  "not implemented for this space < %s >",
1504  FieldSpaceNames[space[ss]]);
1505  }
1507  };
1508 
1509  // evaluate operators with field data
1510  auto evaluate_op_for_fields = [&](auto &op) {
1512 
1513  if (op.getOpType() & UDO::OPROW) {
1514  try {
1515  CHKERR op.opRhs(*op_data[0], false);
1516  }
1517  CATCH_OP_ERRORS(op);
1518  }
1519 
1520  if (op.getOpType() & UDO::OPCOL) {
1521  try {
1522  CHKERR op.opRhs(*op_data[1], false);
1523  }
1524  CATCH_OP_ERRORS(op);
1525  }
1526 
1527  if (op.getOpType() & UDO::OPROWCOL) {
1528  try {
1529  CHKERR op.opLhs(*op_data[0], *op_data[1]);
1530  }
1531  CATCH_OP_ERRORS(op);
1532  }
1534  };
1535 
1536  // collect bit ref level on all entities, fields and spaces
1538 
1539  auto oit = opPtrVector.begin();
1540  auto hi_oit = opPtrVector.end();
1541 
1542  // iterate over all operators on element
1543  for (; oit != hi_oit; oit++) {
1544 
1545  try {
1546 
1547  CHKERR oit->setPtrFE(this);
1548 
1549  if ((oit->opType & UDO::OPSPACE) == UDO::OPSPACE) {
1550 
1551  // run operator for space but specific field
1552  CHKERR evaluate_op_space(*oit);
1553 
1554  } else if (
1555 
1556  (oit->opType & (UDO::OPROW | UDO::OPCOL | UDO::OPROWCOL)) ==
1557  oit->opType
1558 
1559  ) {
1560 
1561  field_name[0] = oit->rowFieldName;
1562  field_name[1] = oit->colFieldName;
1563 
1564  for (size_t ss = 0; ss != 2; ss++) {
1565 
1566 #ifndef NDEBUG
1567  if (field_name[ss].empty()) {
1568  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1569  "Not set Field name in operator %d (0-row, 1-column) in "
1570  "operator %s",
1571  ss,
1572  (boost::typeindex::type_id_runtime(*oit).pretty_name())
1573  .c_str());
1574  }
1575 #endif
1576 
1577  field_struture[ss] = mField.get_field_structure(field_name[ss]);
1578  field_id[ss] = field_struture[ss]->getBitNumber();
1579  space[ss] = field_struture[ss]->getSpace();
1580  base[ss] = field_struture[ss]->getApproxBase();
1581  }
1582 
1583  // not that if field name do not change between operators, entity field
1584  // data are nor rebuild
1585  for (size_t ss = 0; ss != 2; ss++) {
1586 
1587  if (oit->getOpType() & types[ss] ||
1588  oit->getOpType() & UDO::OPROWCOL) {
1589  if (last_eval_field_id[ss] != field_id[ss]) {
1590  CHKERR set_op_entities_data(ss, *oit);
1591  last_eval_field_id[ss] = field_id[ss];
1592  }
1593  }
1594  }
1595 
1596  CHKERR swap_bases(*oit);
1597 
1598  // run operator for given field or row, column or both
1599  CHKERR evaluate_op_for_fields(*oit);
1600 
1601  CHKERR swap_bases(*oit);
1602 
1603  } else {
1604 
1605  for (int i = 0; i != 5; ++i)
1606  if (oit->opType & (1 << i))
1607  MOFEM_LOG("SELF", Sev::error) << UDO::OpTypeNames[i];
1608  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1609  "Impossible operator type");
1610  }
1611  }
1612  CATCH_OP_ERRORS(*oit);
1613  }
1615 }
1616 
1618  "OPROW", " OPCOL", "OPROWCOL", "OPSPACE", "OPLAST"};
1619 
1621  const std::string field_name, const EntityType type, const int side,
1622  VectorInt &indices) const {
1624  if (ptrFE == NULL)
1625  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1626 
1627  switch (type) {
1628  case MBVERTEX:
1630  break;
1631  default:
1633  }
1635 }
1636 
1638  const std::string field_name, const EntityType type, const int side,
1639  VectorInt &indices) const {
1641  if (ptrFE == NULL)
1642  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1643 
1644  switch (type) {
1645  case MBVERTEX:
1646  CHKERR ptrFE->getProblemNodesColIndices(field_name, indices);
1647  break;
1648  default:
1649  CHKERR ptrFE->getProblemTypeColIndices(field_name, type, side, indices);
1650  }
1652 }
1653 
1657  sidePtrFE = const_cast<ForcesAndSourcesCore *>(side_fe_ptr);
1659 }
1660 
1662  const ForcesAndSourcesCore *refine_fe_ptr) {
1664  refinePtrFE = const_cast<ForcesAndSourcesCore *>(refine_fe_ptr);
1666 }
1667 
1669  const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t side_dim,
1670  const EntityHandle ent_for_side, const int verb,
1671  const LogManager::SeverityLevel sev, AdjCache *adj_cache) {
1673 
1674  const auto *problem_ptr = getFEMethod()->problemPtr;
1675  auto &numered_fe =
1676  problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1677 
1678  auto fe_miit = ptrFE->mField.get_finite_elements()
1679  ->get<FiniteElement_name_mi_tag>()
1680  .find(fe_name);
1681  if (fe_miit != ptrFE->mField.get_finite_elements()
1682  ->get<FiniteElement_name_mi_tag>()
1683  .end()) {
1684 
1685  const auto ent = ent_for_side ? ent_for_side : getFEEntityHandle();
1686 
1687  side_fe->feName = fe_name;
1688 
1689  CHKERR side_fe->setSideFEPtr(ptrFE);
1690  CHKERR side_fe->copyBasicMethod(*getFEMethod());
1691  CHKERR side_fe->copyPetscData(*getFEMethod());
1692  CHKERR side_fe->copyKsp(*getFEMethod());
1693  CHKERR side_fe->copySnes(*getFEMethod());
1694  CHKERR side_fe->copyTs(*getFEMethod());
1695 
1696  side_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1697 
1698  CHKERR side_fe->preProcess();
1699 
1700  std::vector<boost::weak_ptr<NumeredEntFiniteElement>> fe_vec;
1701  auto get_numered_fe_ptr = [&](auto &fe_uid, Range &&adjacent_ents)
1702  -> std::vector<boost::weak_ptr<NumeredEntFiniteElement>> & {
1703  fe_vec.reserve(adjacent_ents.size());
1704  for (auto fe_ent : adjacent_ents) {
1705  auto miit = numered_fe.find(
1707  if (miit != numered_fe.end()) {
1708  fe_vec.emplace_back(*miit);
1709  }
1710  }
1711  return fe_vec;
1712  };
1713 
1714  auto get_bit_entity_adjacency = [&]() {
1715  Range adjacent_ents;
1716  CHKERR ptrFE->mField.getInterface<BitRefManager>()->getAdjacenciesAny(
1717  ent, side_dim, adjacent_ents);
1718  return adjacent_ents;
1719  };
1720 
1721  auto get_adj = [&](auto &fe_uid)
1722  -> std::vector<boost::weak_ptr<NumeredEntFiniteElement>> & {
1723  if (adj_cache) {
1724  try {
1725  return (*adj_cache).at(ent);
1726  } catch (const std::out_of_range &) {
1727  return (*adj_cache)[ent] =
1728  get_numered_fe_ptr(fe_uid, get_bit_entity_adjacency());
1729  }
1730  } else {
1731  return get_numered_fe_ptr(fe_uid, get_bit_entity_adjacency());
1732  }
1733  };
1734 
1735  auto adj = get_adj((*fe_miit)->getFEUId());
1736 
1737  int nn = 0;
1738  side_fe->loopSize = adj.size();
1739  for (auto fe_weak_ptr : adj) {
1740  if (auto fe_ptr = fe_weak_ptr.lock()) {
1741  if (verb >= VERBOSE)
1742  MOFEM_LOG("SELF", sev) << "Side finite element "
1743  << "(" << nn << "): " << *fe_ptr;
1744  side_fe->nInTheLoop = nn;
1745  side_fe->numeredEntFiniteElementPtr = fe_ptr;
1746  CHKERR (*side_fe)();
1747  }
1748  ++nn;
1749  }
1750 
1751  CHKERR side_fe->postProcess();
1752  }
1753 
1755 }
1756 
1758  const string &fe_name, ForcesAndSourcesCore *this_fe, const int verb,
1759  const LogManager::SeverityLevel sev) {
1761 
1762  if (verb >= VERBOSE)
1763  MOFEM_LOG("SELF", sev) << "This finite element: "
1764  << *getNumeredEntFiniteElementPtr();
1765 
1766  this_fe->feName = fe_name;
1767 
1768  CHKERR this_fe->setRefineFEPtr(ptrFE);
1769  CHKERR this_fe->copyBasicMethod(*getFEMethod());
1770  CHKERR this_fe->copyPetscData(*getFEMethod());
1771  CHKERR this_fe->copyKsp(*getFEMethod());
1772  CHKERR this_fe->copySnes(*getFEMethod());
1773  CHKERR this_fe->copyTs(*getFEMethod());
1774 
1775  this_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1776 
1777  CHKERR this_fe->preProcess();
1778 
1779  this_fe->nInTheLoop = getNinTheLoop();
1780  this_fe->numeredEntFiniteElementPtr = getNumeredEntFiniteElementPtr();
1781 
1782  CHKERR (*this_fe)();
1783 
1784  CHKERR this_fe->postProcess();
1785 
1787 }
1788 
1790  const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb,
1791  const LogManager::SeverityLevel sev) {
1793 
1794  auto &fes =
1795  ptrFE->mField.get_finite_elements()->get<FiniteElement_name_mi_tag>();
1796  auto fe_miit = fes.find(fe_name);
1797  if (fe_miit != fes.end()) {
1798 
1799  const auto *problem_ptr = getFEMethod()->problemPtr;
1800  auto &numered_fe =
1801  problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1802 
1803 
1804  parent_fe->feName = fe_name;
1805  CHKERR parent_fe->setRefineFEPtr(ptrFE);
1806  CHKERR parent_fe->copyBasicMethod(*getFEMethod());
1807  CHKERR parent_fe->copyPetscData(*getFEMethod());
1808  CHKERR parent_fe->copyKsp(*getFEMethod());
1809  CHKERR parent_fe->copySnes(*getFEMethod());
1810  CHKERR parent_fe->copyTs(*getFEMethod());
1811 
1812  parent_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1813 
1814  const auto parent_ent = getNumeredEntFiniteElementPtr()->getParentEnt();
1815  auto miit = numered_fe.find(EntFiniteElement::getLocalUniqueIdCalculate(
1816  parent_ent, (*fe_miit)->getFEUId()));
1817  if (miit != numered_fe.end()) {
1818  if (verb >= VERBOSE)
1819  MOFEM_LOG("SELF", sev) << "Parent finite element: " << **miit;
1820  parent_fe->loopSize = 1;
1821  parent_fe->nInTheLoop = 0;
1822  parent_fe->numeredEntFiniteElementPtr = *miit;
1823  CHKERR parent_fe->preProcess();
1824  CHKERR (*parent_fe)();
1825  CHKERR parent_fe->postProcess();
1826  } else {
1827  if (verb >= VERBOSE)
1828  MOFEM_LOG("SELF", sev) << "Parent finite element: no parent";
1829  parent_fe->loopSize = 0;
1830  parent_fe->nInTheLoop = 0;
1831  CHKERR parent_fe->preProcess();
1832  CHKERR parent_fe->postProcess();
1833  }
1834  }
1835 
1837 }
1838 
1840  const string &fe_name, ForcesAndSourcesCore *child_fe, const int verb,
1841  const LogManager::SeverityLevel sev) {
1843 
1844  auto fe_miit = ptrFE->mField.get_finite_elements()
1845  ->get<FiniteElement_name_mi_tag>()
1846  .find(fe_name);
1847  if (fe_miit != ptrFE->mField.get_finite_elements()
1848  ->get<FiniteElement_name_mi_tag>()
1849  .end()) {
1850 
1851  const auto *problem_ptr = getFEMethod()->problemPtr;
1852  auto &ref_ents = *getPtrFE()->mField.get_ref_ents();
1853  auto &numered_fe =
1854  problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1855 
1856  const auto parent_ent = getNumeredEntFiniteElementPtr()->getEnt();
1857  const auto parent_type = getNumeredEntFiniteElementPtr()->getEntType();
1858  auto range =
1859  ref_ents.get<Composite_ParentEnt_And_EntType_mi_tag>().equal_range(
1860  boost::make_tuple(parent_type, parent_ent));
1861 
1862  if (auto size = std::distance(range.first, range.second)) {
1863 
1864  std::vector<EntityHandle> childs_vec;
1865  childs_vec.reserve(size);
1866  for (; range.first != range.second; ++range.first)
1867  childs_vec.emplace_back((*range.first)->getEnt());
1868 
1869  Range childs;
1870 
1871  if ((childs_vec.back() - childs_vec.front() + 1) == size)
1872  childs = Range(childs_vec.front(), childs_vec.back());
1873  else
1874  childs.insert_list(childs_vec.begin(), childs_vec.end());
1875 
1876  child_fe->feName = fe_name;
1877 
1878  CHKERR child_fe->setRefineFEPtr(ptrFE);
1879  CHKERR child_fe->copyBasicMethod(*getFEMethod());
1880  CHKERR child_fe->copyPetscData(*getFEMethod());
1881  CHKERR child_fe->copyKsp(*getFEMethod());
1882  CHKERR child_fe->copySnes(*getFEMethod());
1883  CHKERR child_fe->copyTs(*getFEMethod());
1884 
1885  child_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1886 
1887  CHKERR child_fe->preProcess();
1888 
1889  int nn = 0;
1890  child_fe->loopSize = size;
1891 
1892  for (auto p = childs.pair_begin(); p != childs.pair_end(); ++p) {
1893 
1894  auto miit =
1895  numered_fe.lower_bound(EntFiniteElement::getLocalUniqueIdCalculate(
1896  p->first, (*fe_miit)->getFEUId()));
1897  auto hi_miit =
1898  numered_fe.upper_bound(EntFiniteElement::getLocalUniqueIdCalculate(
1899  p->second, (*fe_miit)->getFEUId()));
1900 
1901  for (; miit != hi_miit; ++miit) {
1902 
1903  if (verb >= VERBOSE)
1904  MOFEM_LOG("SELF", sev) << "Child finite element "
1905  << "(" << nn << "): " << **miit;
1906 
1907  child_fe->nInTheLoop = nn++;
1908  child_fe->numeredEntFiniteElementPtr = *miit;
1909  CHKERR (*child_fe)();
1910  }
1911  }
1912  }
1913 
1914  CHKERR child_fe->postProcess();
1915  }
1916 
1918 }
1919 
1920 int ForcesAndSourcesCore::getRule(int order_row, int order_col,
1921  int order_data) {
1922  return getRuleHook ? getRuleHook(order_row, order_col, order_data)
1923  : getRule(order_data);
1924 }
1925 
1927  int order_data) {
1928  return setRuleHook ? setRuleHook(this, order_row, order_col, order_data)
1929  : setGaussPts(order_data);
1930 }
1931 
1933 
1934 /** \deprecated setGaussPts(int row_order, int col_order, int data order);
1935  */
1938  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Sorry, not implemented");
1940 }
1941 
1943  const char type,
1944  const bool symm)
1945  : DataOperator(symm), opType(type), sPace(space), ptrFE(nullptr) {}
1946 
1948  const std::string field_name, const char type, const bool symm)
1949  : DataOperator(symm), opType(type), rowFieldName(field_name),
1950  colFieldName(field_name), sPace(LASTSPACE), ptrFE(nullptr) {}
1951 
1953  const std::string row_field_name, const std::string col_field_name,
1954  const char type, const bool symm)
1955  : DataOperator(symm), opType(type), rowFieldName(row_field_name),
1956  colFieldName(col_field_name), sPace(LASTSPACE), ptrFE(nullptr) {}
1957 
1960  if (preProcessHook) {
1961  ierr = preProcessHook();
1962  CHKERRG(ierr);
1963  }
1965 }
1968  if (operatorHook) {
1969  ierr = operatorHook();
1970  CHKERRG(ierr);
1971  } else {
1972 #ifndef NDEBUG
1973  MOFEM_LOG("SELF", Sev::warning)
1974  << "No method operator() overloaded on element entity on finite "
1975  "element <"
1976  << boost::typeindex::type_id_runtime(*this).pretty_name() << ">";
1977 #endif
1978  }
1980 }
1983  if (postProcessHook) {
1984  ierr = postProcessHook();
1985  CHKERRG(ierr);
1986  }
1988 }
1989 
1993  if (!(ptrFE = dynamic_cast<ForcesAndSourcesCore *>(ptr)))
1994  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1995  "User operator and finite element do not work together");
1997 }
1998 
1999 } // namespace MoFEM
NOSPACE
@ NOSPACE
Definition: definitions.h:83
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
MoFEM::FEMethod::getRowFieldEnts
auto & getRowFieldEnts() const
Definition: LoopMethods.hpp:423
MoFEM::ForcesAndSourcesCore::getEntityRowIndices
MoFEMErrorCode getEntityRowIndices(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
Definition: ForcesAndSourcesCore.cpp:412
MoFEM::BasicMethod::getNinTheLoop
int getNinTheLoop() const
get number of evaluated element in the loop
Definition: LoopMethods.hpp:207
MoFEM::ForcesAndSourcesCore::UserDataOperator::AdjCache
std::map< EntityHandle, std::vector< boost::weak_ptr< NumeredEntFiniteElement > >> AdjCache
Definition: ForcesAndSourcesCore.hpp:904
MoFEM::Field::getNbOfCoeffs
FieldCoefficientsNumber getNbOfCoeffs() const
Get number of field coefficients.
Definition: FieldMultiIndices.hpp:188
CATCH_OP_ERRORS
#define CATCH_OP_ERRORS(OP)
Definition: ForcesAndSourcesCore.cpp:1323
MoFEM::EntitiesFieldData::facesNodesOrder
MatrixInt facesNodesOrder
order of face nodes on element
Definition: EntitiesFieldData.hpp:47
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::ForcesAndSourcesCore::getNoFieldFieldData
MoFEMErrorCode getNoFieldFieldData(const int bit_number, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.
Definition: ForcesAndSourcesCore.cpp:855
MoFEM::ForcesAndSourcesCore::calBernsteinBezierBaseFunctionsOnElement
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement()
Calculate Bernstein-Bezier base.
Definition: ForcesAndSourcesCore.cpp:1106
MoFEM::DataOperator
base operator to do operations at Gauss Pt. level
Definition: DataOperators.hpp:24
MoFEM::BasicMethod::nInTheLoop
int nInTheLoop
number currently of processed method
Definition: LoopMethods.hpp:198
LASTBASE
@ LASTBASE
Definition: definitions.h:69
MoFEM::ForcesAndSourcesCore::getProblemNodesIndices
MoFEMErrorCode getProblemNodesIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
get indices of nodal indices which are declared for problem but not this particular element
Definition: ForcesAndSourcesCore.cpp:490
MoFEM::ForcesAndSourcesCore::getMaxRowOrder
int getMaxRowOrder() const
Get max order of approximation for field in rows.
Definition: ForcesAndSourcesCore.cpp:131
MoFEM::ForcesAndSourcesCore::UserDataOperator::getProblemColIndices
MoFEMErrorCode getProblemColIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get col indices.
Definition: ForcesAndSourcesCore.cpp:1637
EntityHandle
MoFEM::FEMethod::getColDofsPtr
auto getColDofsPtr() const
Definition: LoopMethods.hpp:443
MoFEM::ForcesAndSourcesCore::derivedDataOnElement
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
Definition: ForcesAndSourcesCore.hpp:468
MoFEM::EntitiesFieldData::basesOnEntities
std::array< std::bitset< LASTBASE >, MBMAXTYPE > basesOnEntities
bases on entity types
Definition: EntitiesFieldData.hpp:52
MoFEM::ForcesAndSourcesCore::UserDataOperator::OpTypeNames
static const char *const OpTypeNames[]
Definition: ForcesAndSourcesCore.hpp:574
MoFEM::ForcesAndSourcesCore::getFaceNodes
MoFEMErrorCode getFaceNodes(EntitiesFieldData &data) const
Get nodes on faces.
Definition: ForcesAndSourcesCore.cpp:921
MoFEM::ForcesAndSourcesCore::getBitRefLevelOnData
MoFEMErrorCode getBitRefLevelOnData()
Definition: ForcesAndSourcesCore.cpp:595
MoFEM::ForcesAndSourcesCore::getRuleHook
RuleHookFun getRuleHook
Hook to get rule.
Definition: ForcesAndSourcesCore.hpp:42
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::FEMethod::feName
std::string feName
Name of finite element.
Definition: LoopMethods.hpp:380
MoFEM::FieldEntity::getLoLocalEntityBitNumber
static UId getLoLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
Definition: FieldEntsMultiIndices.hpp:247
MoFEM::cmp_uid_hi
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
Definition: ForcesAndSourcesCore.cpp:29
NOBASE
@ NOBASE
Definition: definitions.h:59
MoFEM::PetscData::copyPetscData
MoFEMErrorCode copyPetscData(const PetscData &petsc_data)
Definition: LoopMethods.cpp:31
ApproximationBaseNames
const static char *const ApproximationBaseNames[]
Definition: definitions.h:72
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
MoFEM::CoreInterface::get_field_structure
virtual const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const =0
get field structure
MoFEM::ForcesAndSourcesCore::setRefineFEPtr
MoFEMErrorCode setRefineFEPtr(const ForcesAndSourcesCore *refine_fe_ptr)
Set the pointer to face element refined.
Definition: ForcesAndSourcesCore.cpp:1661
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Field::getSpace
FieldSpace getSpace() const
Get field approximation space.
Definition: FieldMultiIndices.hpp:150
MoFEM::FieldEntity_vector_view
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
Definition: FieldEntsMultiIndices.hpp:478
MoFEM::ForcesAndSourcesCore::setGaussPts
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
Definition: ForcesAndSourcesCore.cpp:1926
MoFEM::CoreInterface::get_field_bit_number
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
MoFEM::FEMethod::getDataFieldEnts
const FieldEntity_vector_view & getDataFieldEnts() const
Definition: LoopMethods.hpp:412
MoFEM::ForcesAndSourcesCore::getMaxDataOrder
int getMaxDataOrder() const
Get max order of approximation for data fields.
Definition: ForcesAndSourcesCore.cpp:120
MoFEM::Types::BitFieldId
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:42
MoFEM::ForcesAndSourcesCore::getEntityDataOrder
MoFEMErrorCode getEntityDataOrder(const EntityType type, const FieldSpace space, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
Get the entity data order.
Definition: ForcesAndSourcesCore.cpp:139
MoFEM::BasicMethod::preProcessHook
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
Definition: LoopMethods.hpp:294
MoFEM::ForcesAndSourcesCore::sidePtrFE
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
Definition: ForcesAndSourcesCore.hpp:507
MoFEM::ForcesAndSourcesCore::getElementPolynomialBase
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
Definition: ForcesAndSourcesCore.hpp:90
MoFEM::ForcesAndSourcesCore::preProcess
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
Definition: ForcesAndSourcesCore.cpp:1958
MoFEM::ForcesAndSourcesCore::getMaxColOrder
int getMaxColOrder() const
Get max order of approximation for field in columns.
Definition: ForcesAndSourcesCore.cpp:135
MoFEM::ForcesAndSourcesCore::calHierarchicalBaseFunctionsOnElement
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
Definition: ForcesAndSourcesCore.cpp:1086
MoFEM::ForcesAndSourcesCore::getNoFieldRowIndices
MoFEMErrorCode getNoFieldRowIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
Definition: ForcesAndSourcesCore.cpp:459
MoFEM::ForcesAndSourcesCore::getNoFieldIndices
MoFEMErrorCode getNoFieldIndices(const int bit_number, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
get NoField indices
Definition: ForcesAndSourcesCore.cpp:442
MoFEM::EntitiesFieldData::spacesOnEntities
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
Definition: EntitiesFieldData.hpp:50
MoFEM::Types::BitFEId
std::bitset< BITFEID_SIZE > BitFEId
Finite element Id.
Definition: Types.hpp:43
MoFEM::ForcesAndSourcesCore::getProblemTypeRowIndices
MoFEMErrorCode getProblemTypeRowIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
Definition: ForcesAndSourcesCore.cpp:572
THROW_MESSAGE
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
MoFEM::ForcesAndSourcesCore::getProblemTypeColIndices
MoFEMErrorCode getProblemTypeColIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
Definition: ForcesAndSourcesCore.cpp:586
MoFEM::EntitiesFieldData::facesNodes
MatrixInt facesNodes
nodes on finite element faces
Definition: EntitiesFieldData.hpp:46
sdf.r
int r
Definition: sdf.py:8
MoFEM::Field
Provide data structure for (tensor) field approximation.
Definition: FieldMultiIndices.hpp:51
USER_BASE
@ USER_BASE
user implemented approximation base
Definition: definitions.h:68
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::ForcesAndSourcesCore::getEntityFieldData
MoFEMErrorCode getEntityFieldData(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
Definition: ForcesAndSourcesCore.cpp:770
MoFEM::EntitiesFieldData::basesOnSpaces
std::array< std::bitset< LASTBASE >, LASTSPACE > basesOnSpaces
base on spaces
Definition: EntitiesFieldData.hpp:54
MoFEM::BasicMethod::fieldsPtr
const Field_multiIndex * fieldsPtr
raw pointer to fields container
Definition: LoopMethods.hpp:252
NumeredDofEntity_multiIndex
multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, UId, &NumeredDofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Part_mi_tag >, member< NumeredDofEntity, unsigned int, &NumeredDofEntity::pArt > >, ordered_non_unique< tag< Idx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::dofIdx > >, ordered_non_unique< tag< PetscGlobalIdx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::petscGloablDofIdx > >, ordered_non_unique< tag< PetscLocalIdx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::petscLocalDofIdx > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, EntityHandle, &NumeredDofEntity::getEnt > > > > NumeredDofEntity_multiIndex
MultiIndex container keeps NumeredDofEntity.
Definition: DofsMultiIndices.hpp:469
FieldSpace
FieldSpace
approximation spaces
Definition: definitions.h:82
VERBOSE
@ VERBOSE
Definition: definitions.h:209
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEM::ForcesAndSourcesCore::getNoFieldColIndices
MoFEMErrorCode getNoFieldColIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
Definition: ForcesAndSourcesCore.cpp:474
MoFEM::FEMethod::getRowDofsPtr
auto getRowDofsPtr() const
Definition: LoopMethods.hpp:439
MoFEM::CoreInterface::get_field_id
virtual BitFieldId get_field_id(const std::string &name) const =0
Get field Id.
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM::ForcesAndSourcesCore::refinePtrFE
ForcesAndSourcesCore * refinePtrFE
Element to integrate parent or child.
Definition: ForcesAndSourcesCore.hpp:524
MoFEM::ForcesAndSourcesCore::getUserPolynomialBase
auto & getUserPolynomialBase()
Get the User Polynomial Base object.
Definition: ForcesAndSourcesCore.hpp:97
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
a
constexpr double a
Definition: approx_sphere.cpp:30
MoFEM::ForcesAndSourcesCore::UserDataOperator::getProblemRowIndices
MoFEMErrorCode getProblemRowIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get row indices.
Definition: ForcesAndSourcesCore.cpp:1620
MoFEM::ForcesAndSourcesCore::getEntityIndices
MoFEMErrorCode getEntityIndices(EntitiesFieldData &data, const int bit_number, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
Definition: ForcesAndSourcesCore.cpp:345
MoFEM::ForcesAndSourcesCore::getColNodesIndices
MoFEMErrorCode getColNodesIndices(EntitiesFieldData &data, const int bit_number) const
get col node indices from FENumeredDofEntity_multiIndex
Definition: ForcesAndSourcesCore.cpp:328
MoFEM::Problem::numeredColDofsPtr
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this problem
Definition: ProblemsMultiIndices.hpp:75
MoFEM::Types::UId
uint128_t UId
Unique Id.
Definition: Types.hpp:31
MoFEM::FEMethod::getNumberOfNodes
auto getNumberOfNodes() const
Definition: LoopMethods.hpp:456
convert.type
type
Definition: convert.py:64
MoFEM::FieldEntity::getLocalUniqueIdCalculate
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
Definition: FieldEntsMultiIndices.hpp:136
SideNumber_multiIndex
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, signed char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
Definition: RefEntsMultiIndices.hpp:101
MoFEM::EntFiniteElement::getLocalUniqueIdCalculate
UId getLocalUniqueIdCalculate() const
Generate UId for finite element entity.
Definition: FEMultiIndices.hpp:528
MoFEM::ForcesAndSourcesCore::opPtrVector
boost::ptr_deque< UserDataOperator > opPtrVector
Vector of finite element users data operators.
Definition: ForcesAndSourcesCore.hpp:480
MoFEM::ForcesAndSourcesCore::dataOnElement
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
Definition: ForcesAndSourcesCore.hpp:461
MoFEM::SnesMethod::copySnes
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:75
MoFEM::ForcesAndSourcesCore::createDataOnElement
MoFEMErrorCode createDataOnElement(EntityType type)
Create a entity data on element object.
Definition: ForcesAndSourcesCore.cpp:1305
MoFEM::FieldEntity::getLoBitNumberUId
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
Definition: FieldEntsMultiIndices.hpp:222
MoFEM::cmp_uid_lo
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
Definition: ForcesAndSourcesCore.cpp:18
MoFEM::ForcesAndSourcesCore::UserDataOperator::setPtrFE
virtual MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr)
Definition: ForcesAndSourcesCore.cpp:1991
MoFEM::ForcesAndSourcesCore::ForcesAndSourcesCore
ForcesAndSourcesCore(Interface &m_field)
Definition: ForcesAndSourcesCore.cpp:40
MoFEM::VectorDofs
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
Definition: EntitiesFieldData.hpp:23
lapack_wrap.h
MoFEM::EntitiesFieldData::sPace
std::bitset< LASTSPACE > sPace
spaces on element
Definition: EntitiesFieldData.hpp:42
MoFEM::BasicMethod::postProcessHook
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
Definition: LoopMethods.hpp:299
LASTSPACE
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:89
MoFEM::ForcesAndSourcesCore::getRule
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
Definition: ForcesAndSourcesCore.cpp:1920
MoFEM::get_id_for_max_type
EntityHandle get_id_for_max_type()
Definition: RefEntsMultiIndices.hpp:13
MoFEM::ForcesAndSourcesCore::getEntityColIndices
MoFEMErrorCode getEntityColIndices(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
Definition: ForcesAndSourcesCore.cpp:427
MoFEM::ForcesAndSourcesCore::operator()
virtual MoFEMErrorCode operator()()
function is run for every finite element
Definition: ForcesAndSourcesCore.cpp:1966
AINSWORTH_LOBATTO_BASE
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
MoFEM::FEMethod::getFEName
auto getFEName() const
get finite element name
Definition: LoopMethods.hpp:391
MoFEM::LogManager::SeverityLevel
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
MoFEM::BasicMethod::loopSize
int loopSize
local number oe methods to process
Definition: LoopMethods.hpp:203
t
constexpr double t
plate stiffness
Definition: plate.cpp:59
MoFEM::EntitiesFieldData::bAse
std::bitset< LASTBASE > bAse
bases on element
Definition: EntitiesFieldData.hpp:45
MoFEM::ForcesAndSourcesCore::UserDataOperator
friend class UserDataOperator
Definition: ForcesAndSourcesCore.hpp:482
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::ForcesAndSourcesCore::getProblemTypeIndices
MoFEMErrorCode getProblemTypeIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, EntityType type, int side_number, VectorInt &indices) const
get indices by type (generic function) which are declared for problem but not this particular element
Definition: ForcesAndSourcesCore.cpp:531
MoFEM::ForcesAndSourcesCore::getSpacesAndBaseOnEntities
MoFEMErrorCode getSpacesAndBaseOnEntities(EntitiesFieldData &data) const
Get field approximation space and base on entities.
Definition: ForcesAndSourcesCore.cpp:990
EntData
EntitiesFieldData::EntData EntData
Definition: child_and_parent.cpp:37
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
AINSWORTH_BERNSTEIN_BEZIER_BASE
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
MoFEM::ForcesAndSourcesCore::setSideFEPtr
MoFEMErrorCode setSideFEPtr(const ForcesAndSourcesCore *side_fe_ptr)
Set the pointer to face element on the side.
Definition: ForcesAndSourcesCore.cpp:1655
convert.n
n
Definition: convert.py:82
MoFEM::ForcesAndSourcesCore::getRowNodesIndices
MoFEMErrorCode getRowNodesIndices(EntitiesFieldData &data, const int bit_number) const
get row node indices from FENumeredDofEntity_multiIndex
Definition: ForcesAndSourcesCore.cpp:311
MoFEM::interface_RefEntity::getSideNumberPtr
boost::shared_ptr< SideNumber > getSideNumberPtr() const
Get the Side number.
Definition: RefEntsMultiIndices.hpp:589
MoFEM::BasicMethod::copyBasicMethod
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
Definition: LoopMethods.cpp:117
MoFEM::Problem::numeredRowDofsPtr
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
Definition: ProblemsMultiIndices.hpp:73
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
Range
FieldSpaceNames
const static char *const FieldSpaceNames[]
Definition: definitions.h:92
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::BasicMethod::operatorHook
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
Definition: LoopMethods.hpp:304
MoFEM::Types::ApproximationOrder
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:26
DEMKOWICZ_JACOBI_BASE
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
MoFEM::Field::getBitNumber
FieldBitNumber getBitNumber() const
Get number of set bit in Field ID. Each field has uid, get getBitNumber get number of bit set for giv...
Definition: FieldMultiIndices.hpp:197
MoFEM::BitFieldId_mi_tag
Definition: TagMultiIndices.hpp:17
MoFEM::ForcesAndSourcesCore::UserDataOperator::loopSide
MoFEMErrorCode loopSide(const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy, AdjCache *adj_cache=nullptr)
User calls this function to loop over elements on the side of face. This function calls finite elemen...
Definition: ForcesAndSourcesCore.cpp:1668
MoFEM::ForcesAndSourcesCore::UserDataOperator::loopChildren
MoFEMErrorCode loopChildren(const string &fe_name, ForcesAndSourcesCore *child_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User calls this function to loop over parent elements. This function calls finite element with its op...
Definition: ForcesAndSourcesCore.cpp:1839
MoFEM::FiniteElement_name_mi_tag
Definition: TagMultiIndices.hpp:26
MoFEM::Types::VectorInt
UBlasVector< int > VectorInt
Definition: Types.hpp:67
MoFEM::ForcesAndSourcesCore::getProblemNodesColIndices
MoFEMErrorCode getProblemNodesColIndices(const std::string &field_name, VectorInt &nodes_indices) const
Definition: ForcesAndSourcesCore.cpp:579
MoFEM::TSMethod::copyTs
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data.
Definition: LoopMethods.cpp:96
UBlasVector< double >
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
MoFEM::BasicMethod::cacheWeakPtr
boost::weak_ptr< CacheTuple > cacheWeakPtr
Definition: LoopMethods.hpp:351
MoFEM::getMaxOrder
static int getMaxOrder(const ENTMULTIINDEX &multi_index)
Definition: ForcesAndSourcesCore.cpp:110
MoFEM::KspMethod::copyKsp
MoFEMErrorCode copyKsp(const KspMethod &ksp)
copy data form another method
Definition: LoopMethods.cpp:55
HCURL
@ HCURL
field with continuous tangents
Definition: definitions.h:86
MoFEM::BitFieldId_space_mi_tag
Definition: TagMultiIndices.hpp:68
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::BasicMethod::problemPtr
const Problem * problemPtr
raw pointer to problem
Definition: LoopMethods.hpp:250
MoFEM::ForcesAndSourcesCore::gaussPts
MatrixDouble gaussPts
Matrix of integration points.
Definition: ForcesAndSourcesCore.hpp:109
MoFEM::Composite_ParentEnt_And_EntType_mi_tag
Definition: TagMultiIndices.hpp:78
MoFEM::FEMethod::getColFieldEnts
auto & getColFieldEnts() const
Definition: LoopMethods.hpp:431
FieldApproximationBase
FieldApproximationBase
approximation base
Definition: definitions.h:58
MoFEM::FieldEntity::getHiBitNumberUId
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
Definition: FieldEntsMultiIndices.hpp:228
MoFEM::EntitiesFieldData::dataOnEntities
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
Definition: EntitiesFieldData.hpp:56
MoFEM::FieldEntity::getHiLocalEntityBitNumber
static UId getHiLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
Definition: FieldEntsMultiIndices.hpp:258
MoFEM::VectorFieldEntities
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
Definition: EntitiesFieldData.hpp:31
MoFEM::ForcesAndSourcesCore::loopOverOperators
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
Definition: ForcesAndSourcesCore.cpp:1347
MoFEM::BitRefManager
Managing BitRefLevels.
Definition: BitRefManager.hpp:21
MoFEM::FEMethod::getFEEntityHandle
EntityHandle getFEEntityHandle() const
Definition: LoopMethods.hpp:460
MoFEM::ForcesAndSourcesCore::setRuleHook
GaussHookFun setRuleHook
Set function to calculate integration rule.
Definition: ForcesAndSourcesCore.hpp:48
MoFEM::FEMethod::numeredEntFiniteElementPtr
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
Definition: LoopMethods.hpp:383
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
MoFEM::Unique_mi_tag
Definition: TagMultiIndices.hpp:18
MoFEM::ForcesAndSourcesCore::getNodesIndices
MoFEMErrorCode getNodesIndices(const int bit_number, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
get node indices
Definition: ForcesAndSourcesCore.cpp:219
MoFEM::field_it
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
Definition: Projection10NodeCoordsOnField.cpp:123
MoFEM::ForcesAndSourcesCore::UserDataOperator::loopParent
MoFEMErrorCode loopParent(const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User calls this function to loop over parent elements. This function calls finite element with its op...
Definition: ForcesAndSourcesCore.cpp:1789
MoFEM::ForcesAndSourcesCore::getNodesFieldData
MoFEMErrorCode getNodesFieldData(EntitiesFieldData &data, const int bit_number) const
Get data on nodes.
Definition: ForcesAndSourcesCore.cpp:642
MoFEM::FEDofEntity
keeps information about indexed dofs for the finite element
Definition: DofsMultiIndices.hpp:267
MoFEM::EntitiesFieldData
data structure for finite element entity
Definition: EntitiesFieldData.hpp:40
MoFEM::FEMethod::getDataFieldEntsPtr
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const
Definition: LoopMethods.hpp:417
MoFEM::ForcesAndSourcesCore::postProcess
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
Definition: ForcesAndSourcesCore.cpp:1981
MoFEM::get_id_for_min_type
EntityHandle get_id_for_min_type()
Definition: RefEntsMultiIndices.hpp:18
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEM::ForcesAndSourcesCore::UserDataOperator::loopThis
MoFEMErrorCode loopThis(const string &fe_name, ForcesAndSourcesCore *this_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User calls this function to loop over the same element using a different set of integration points....
Definition: ForcesAndSourcesCore.cpp:1757
MoFEM::ForcesAndSourcesCore::getProblemNodesRowIndices
MoFEMErrorCode getProblemNodesRowIndices(const std::string &field_name, VectorInt &nodes_indices) const
Definition: ForcesAndSourcesCore.cpp:565
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
MoFEM::ForcesAndSourcesCore::UserDataOperator::ptrFE
ForcesAndSourcesCore * ptrFE
Definition: ForcesAndSourcesCore.hpp:984
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
MoFEM::ForcesAndSourcesCore::mField
Interface & mField
Definition: ForcesAndSourcesCore.hpp:24
MoFEM::ForcesAndSourcesCore::getEntitySense
MoFEMErrorCode getEntitySense(const EntityType type, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
get sense (orientation) of entity
Definition: ForcesAndSourcesCore.cpp:84
NOFIELD
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84
MoFEM::ForcesAndSourcesCore::UserDataOperator::UserDataOperator
UserDataOperator(const FieldSpace space, const char type=OPSPACE, const bool symm=true)
Definition: ForcesAndSourcesCore.cpp:1942
MoFEM::DefaultElementAdjacency::getDefTypeMap
static bool getDefTypeMap(const EntityType fe_type, const EntityType ent_type)
Definition: FEMultiIndices.hpp:304