v0.9.2
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 /* This file is part of MoFEM.
7  * MoFEM is free software: you can redistribute it and/or modify it under
8  * the terms of the GNU Lesser General Public License as published by the
9  * Free Software Foundation, either version 3 of the License, or (at your
10  * option) any later version.
11  *
12  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
13  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15  * License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
19 
20 #ifdef __cplusplus
21 extern "C" {
22 #endif
23 #include <cblas.h>
24 #include <lapack_wrap.h>
25 // #include <gm_rule.h>
26 #ifdef __cplusplus
27 }
28 #endif
29 
30 namespace MoFEM {
31 
33  :
34 
35  mField(m_field), getRuleHook(0), setRuleHook(0),
36  dataOnElement{
37 
38  nullptr,
39  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // NOFIELD
40  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // H1
41  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HCURL
42  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HDIV
43  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET) // L2
44 
45  },
46  derivedDataOnElement{
47 
48  nullptr,
49  boost::make_shared<DerivedDataForcesAndSourcesCore>(
50  dataOnElement[NOFIELD]), // NOFIELD
51  boost::make_shared<DerivedDataForcesAndSourcesCore>(
52  dataOnElement[H1]), // H1
53  boost::make_shared<DerivedDataForcesAndSourcesCore>(
54  dataOnElement[HCURL]), // HCURL
55  boost::make_shared<DerivedDataForcesAndSourcesCore>(
56  dataOnElement[HDIV]), // HDIV
57  boost::make_shared<DerivedDataForcesAndSourcesCore>(
58  dataOnElement[L2]) // L2
59 
60  },
61  dataNoField(*dataOnElement[NOFIELD].get()),
62  dataH1(*dataOnElement[H1].get()), dataHcurl(*dataOnElement[HCURL].get()),
63  dataHdiv(*dataOnElement[HDIV].get()), dataL2(*dataOnElement[L2].get()),
64  lastEvaluatedElementEntityType(MBMAXTYPE), sidePtrFE(nullptr) {}
65 
66 // ** Sense **
67 
69  const EntityType type,
70  boost::ptr_vector<DataForcesAndSourcesCore::EntData> &data) const {
72 
73  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable().get<2>();
74  for (auto r = side_table.equal_range(type); r.first != r.second; ++r.first) {
75  const int side_number = (*r.first)->side_number;
76  if (side_number >= 0) {
77  const int brother_side_number = (*r.first)->brother_side_number;
78  const int sense = (*r.first)->sense;
79 
80  data[side_number].getSense() = sense;
81  if (brother_side_number != -1)
82  data[brother_side_number].getSense() = sense;
83  }
84  }
86 }
87 
88 // ** Order **
89 
90 template <typename ENTMULTIINDEX>
91 static inline int getMaxOrder(const ENTMULTIINDEX &multi_index) {
92  int max_order = 0;
93  for (auto ent_field_weak_ptr : multi_index)
94  if (auto e = ent_field_weak_ptr.lock()) {
95  const int order = e->getMaxOrder();
96  max_order = (max_order < order) ? order : max_order;
97  }
98  return max_order;
99 }
100 
102  int max_order = 0;
103  for (auto e : *dataFieldEntsPtr) {
104  const int order = e->getMaxOrder();
105  max_order = (max_order < order) ? order : max_order;
106  }
107  return max_order;
108 }
109 
111  return getMaxOrder(*rowFieldEntsPtr);
112 }
113 
115  return getMaxOrder(*colFieldEntsPtr);
116 }
117 
119  const EntityType type, const FieldSpace space,
120  boost::ptr_vector<DataForcesAndSourcesCore::EntData> &data) const {
122 
123  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable();
124 
125  for (unsigned int s = 0; s != data.size(); ++s)
126  data[s].getDataOrder() = 0;
127 
128  auto &fields_ents =
130 
131  for (auto r = fields_ents.equal_range(boost::make_tuple(type, space));
132  r.first != r.second; ++r.first) {
133 
134  auto &e = **r.first;
135 
136  auto sit = side_table.find(e.getEnt());
137  if (sit != side_table.end()) {
138 
139  auto &side = *sit;
140  const int side_number = side->side_number;
141  if (side_number >= 0) {
142  ApproximationOrder ent_order = e.getMaxOrder();
143  auto &dat = data[side_number];
144  dat.getDataOrder() =
145  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
146  }
147  } else
148  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
149  "Entity on side of the element not found");
150  }
151 
152  for (auto r = side_table.get<2>().equal_range(type); r.first != r.second;
153  ++r.first) {
154  const int brother_side_number = (*r.first)->brother_side_number;
155  if (brother_side_number != -1) {
156  const int side_number = (*r.first)->side_number;
157  data[brother_side_number].getDataOrder() =
158  data[side_number].getDataOrder();
159  }
160  }
161 
163 }
164 
165 // ** Indices **
166 
168  const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs,
169  VectorInt &nodes_indices, VectorInt &local_nodes_indices) const {
171  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
172  auto tuple = boost::make_tuple(field_name, MBVERTEX);
173  auto dit = dofs_by_type.lower_bound(tuple);
174  auto hi_dit = dofs_by_type.upper_bound(tuple);
175 
176  if (dit != hi_dit) {
177 
178  int num_nodes;
179  CHKERR getNumberOfNodes(num_nodes);
180  int max_nb_dofs = 0;
181  const int nb_dofs_on_vert = (*dit)->getNbOfCoeffs();
182  max_nb_dofs = nb_dofs_on_vert * num_nodes;
183  nodes_indices.resize(max_nb_dofs, false);
184  local_nodes_indices.resize(max_nb_dofs, false);
185  if (std::distance(dit, hi_dit) != max_nb_dofs) {
186  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
187  std::fill(local_nodes_indices.begin(), local_nodes_indices.end(), -1);
188  }
189 
190  for (; dit != hi_dit; dit++) {
191  auto &dof = **dit;
192  const int idx = dof.getPetscGlobalDofIdx();
193  const int local_idx = dof.getPetscLocalDofIdx();
194  const int side_number = dof.sideNumberPtr->side_number;
195  const int pos = side_number * nb_dofs_on_vert + dof.getDofCoeffIdx();
196  nodes_indices[pos] = idx;
197  local_nodes_indices[pos] = local_idx;
198  const int brother_side_number =
199  (*dit)->sideNumberPtr->brother_side_number;
200  if (brother_side_number != -1) {
201  const int elem_idx =
202  brother_side_number * nb_dofs_on_vert + (*dit)->getDofCoeffIdx();
203  nodes_indices[elem_idx] = idx;
204  local_nodes_indices[elem_idx] = local_idx;
205  }
206  }
207 
208  } else {
209  nodes_indices.resize(0, false);
210  local_nodes_indices.resize(0, false);
211  }
212 
214 }
215 
218  const std::string &field_name) const {
219  return getNodesIndices(field_name,
220  const_cast<FENumeredDofEntity_multiIndex &>(
221  numeredEntFiniteElementPtr->getRowsDofs()),
222  data.dataOnEntities[MBVERTEX][0].getIndices(),
223  data.dataOnEntities[MBVERTEX][0].getLocalIndices());
224 }
225 
228  const std::string &field_name) const {
229  return getNodesIndices(field_name,
230  const_cast<FENumeredDofEntity_multiIndex &>(
231  numeredEntFiniteElementPtr->getColsDofs()),
232  data.dataOnEntities[MBVERTEX][0].getIndices(),
233  data.dataOnEntities[MBVERTEX][0].getLocalIndices());
234 }
235 
237  DataForcesAndSourcesCore &data, const std::string &field_name,
238  FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo,
239  const EntityType type_hi) const {
241 
242  for (EntityType t = type_lo; t != type_hi; ++t) {
243  for (auto &dat : data.dataOnEntities[t]) {
244  dat.getIndices().resize(0, false);
245  dat.getLocalIndices().resize(0, false);
246  }
247  }
248 
249  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
250  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
251  if (dit == dofs_by_type.end())
253  auto hi_dit =
254  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
255  for (; dit != hi_dit; ++dit) {
256  auto &dof = **dit;
257  const EntityType type = dof.getEntType();
258  const int side = dof.sideNumberPtr->side_number;
259 
260  if (side >= 0) {
261 
262  auto &dat = data.dataOnEntities[type][side];
263  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
264  if (nb_dofs_on_ent) {
265  const int brother_side = dof.sideNumberPtr->brother_side_number;
266  auto &ent_field_indices = dat.getIndices();
267  auto &ent_field_local_indices = dat.getLocalIndices();
268  if (ent_field_indices.empty()) {
269  ent_field_indices.resize(nb_dofs_on_ent, false);
270  ent_field_local_indices.resize(nb_dofs_on_ent, false);
271  std::fill(ent_field_indices.data().begin(),
272  ent_field_indices.data().end(), -1);
273  std::fill(ent_field_local_indices.data().begin(),
274  ent_field_local_indices.data().end(), -1);
275  }
276  const int idx = dof.getEntDofIdx();
277  ent_field_indices[idx] = dof.getPetscGlobalDofIdx();
278  ent_field_local_indices[idx] = dof.getPetscLocalDofIdx();
279  if (brother_side != -1) {
280  auto &dat_brother = data.dataOnEntities[type][brother_side];
281  dat_brother.getIndices().resize(nb_dofs_on_ent, false);
282  dat_brother.getLocalIndices().resize(nb_dofs_on_ent, false);
283  dat_brother.getIndices()[idx] = dat.getIndices()[idx];
284  dat_brother.getLocalIndices()[idx] = dat.getLocalIndices()[idx];
285  }
286  }
287  }
288  }
289 
291 }
292 
294 ForcesAndSourcesCore::getNoFieldIndices(const std::string &field_name,
296  VectorInt &indices) const {
298  auto dit = dofs.get<FieldName_mi_tag>().lower_bound(field_name);
299  auto hi_dit = dofs.get<FieldName_mi_tag>().upper_bound(field_name);
300  indices.resize(std::distance(dit, hi_dit));
301  for (; dit != hi_dit; dit++) {
302  int idx = (*dit)->getPetscGlobalDofIdx();
303  indices[(*dit)->getDofCoeffIdx()] = idx;
304  }
306 }
307 
309  DataForcesAndSourcesCore &data, const std::string &field_name) const {
311  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
312  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
313  }
314  CHKERR getNoFieldIndices(field_name,
315  const_cast<FENumeredDofEntity_multiIndex &>(
316  numeredEntFiniteElementPtr->getRowsDofs()),
317  data.dataOnEntities[MBENTITYSET][0].getIndices());
319 }
320 
322  DataForcesAndSourcesCore &data, const std::string &field_name) const {
324  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
325  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
326  }
327  CHKERR getNoFieldIndices(field_name,
328  const_cast<FENumeredDofEntity_multiIndex &>(
329  numeredEntFiniteElementPtr->getColsDofs()),
330  data.dataOnEntities[MBENTITYSET][0].getIndices());
332 }
333 
334 // ** Indices from problem **
335 
337  const std::string &field_name, const NumeredDofEntity_multiIndex &dofs,
338  VectorInt &nodes_indices) const {
340 
341  const Field *field_struture = mField.get_field_structure(field_name);
342  if (field_struture->getSpace() == H1) {
343 
344  int num_nodes;
345  CHKERR getNumberOfNodes(num_nodes);
346  nodes_indices.resize(field_struture->getNbOfCoeffs() * num_nodes, false);
347  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
348 
349  auto &side_table = const_cast<SideNumber_multiIndex &>(
350  numeredEntFiniteElementPtr->getSideNumberTable());
351  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 0));
352  auto hi_siit =
353  side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 10000));
354 
355  int nn = 0;
356  for (; siit != hi_siit; siit++, nn++) {
357 
358  if (siit->get()->side_number == -1)
359  continue;
360 
361  const EntityHandle ent = siit->get()->ent;
362  auto dit =
363  dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().lower_bound(
364  boost::make_tuple(field_name, ent, 0));
365  auto hi_dit =
366  dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().upper_bound(
367  boost::make_tuple(field_name, ent, 10000)); /// very large number
368  for (; dit != hi_dit; dit++) {
369  nodes_indices[siit->get()->side_number * (*dit)->getNbOfCoeffs() +
370  (*dit)->getDofCoeffIdx()] =
371  (*dit)->getPetscGlobalDofIdx();
372  }
373  }
374  } else {
375  nodes_indices.resize(0, false);
376  }
377 
379 }
380 
382  const std::string &field_name, const NumeredDofEntity_multiIndex &dofs,
383  EntityType type, int side_number, VectorInt &indices) const {
385 
386  indices.resize(0);
387 
388  auto &side_table = const_cast<SideNumber_multiIndex &>(
389  numeredEntFiniteElementPtr->getSideNumberTable());
390  auto siit =
391  side_table.get<1>().lower_bound(boost::make_tuple(type, side_number));
392  auto hi_siit =
393  side_table.get<1>().upper_bound(boost::make_tuple(type, side_number));
394 
395  for (; siit != hi_siit; siit++) {
396 
397  if (siit->get()->side_number == -1)
398  continue;
399 
400  const EntityHandle ent = siit->get()->ent;
401  NumeredDofEntity_multiIndex::index<
402  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator dit,
403  hi_dit;
404  dit = dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().lower_bound(
405  boost::make_tuple(field_name, ent, 0));
406  hi_dit =
407  dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().upper_bound(
408  boost::make_tuple(field_name, ent, 10000)); /// very large number
409 
410  indices.resize(std::distance(dit, hi_dit));
411  for (; dit != hi_dit; dit++) {
412 
413  indices[(*dit)->getEntDofIdx()] = (*dit)->getPetscGlobalDofIdx();
414  }
415  }
416 
418 }
419 
421  const std::string &field_name, VectorInt &nodes_indices) const {
422  return getProblemNodesIndices(field_name, *(problemPtr->numeredDofsRows),
423  nodes_indices);
424 }
425 
428  EntityType type, int side_number,
429  VectorInt &indices) const {
430  return getProblemTypeIndices(field_name, *(problemPtr->numeredDofsRows), type,
431  side_number, indices);
432 }
433 
435  const std::string &field_name, VectorInt &nodes_indices) const {
436  return getProblemNodesIndices(field_name, *(problemPtr->numeredDofsCols),
437  nodes_indices);
438 }
439 
442  EntityType type, int side_number,
443  VectorInt &indices) const {
444  return getProblemTypeIndices(field_name, *(problemPtr->numeredDofsCols), type,
445  side_number, indices);
446 }
447 
448 // ** Data **
449 
452  const std::string &field_name) const {
453 
454  auto get_nodes_field_data = [&](FEDofEntity_multiIndex &dofs,
455  VectorDouble &nodes_data,
456  VectorDofs &nodes_dofs, FieldSpace &space,
458  VectorInt &bb_node_order) {
460 
461  auto &dofs_by_name_and_type = dofs.get<Composite_Name_And_Type_mi_tag>();
462  auto tuple = boost::make_tuple(field_name, MBVERTEX);
463  auto dit = dofs_by_name_and_type.lower_bound(tuple);
464  if (dit == dofs_by_name_and_type.end())
465  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
466  "No nodal dofs on element");
467  auto hi_dit = dofs.get<Composite_Name_And_Type_mi_tag>().upper_bound(tuple);
468 
469  if (dit != hi_dit) {
470  auto &first_dof = **dit;
471  space = first_dof.getSpace();
472  base = first_dof.getApproxBase();
473  int num_nodes;
474  CHKERR getNumberOfNodes(num_nodes);
475  bb_node_order.resize(num_nodes, false);
476  bb_node_order.clear();
477  const int nb_dof_idx = first_dof.getNbOfCoeffs();
478  const int max_nb_dofs = nb_dof_idx * num_nodes;
479  nodes_data.resize(max_nb_dofs, false);
480  nodes_dofs.resize(max_nb_dofs, false);
481  nodes_data.clear();
482 
483  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
484  for (; dit != hi_dit;) {
485  const auto &dof_ptr = *dit;
486  const auto &dof = *dof_ptr;
487  const auto &sn = *dof.sideNumberPtr;
488  const int side_number = sn.side_number;
489  const int brother_side_number = sn.brother_side_number;
490  if (brother_side_number != -1)
491  brother_dofs_vec.emplace_back(dof_ptr);
492 
493  bb_node_order[side_number] = dof.getMaxOrder();
494  int pos = side_number * nb_dof_idx;
495  auto ent_filed_data_vec = dof.getEntFieldData();
496  for (int ii = 0; ii != nb_dof_idx; ++ii) {
497  nodes_data[pos] = ent_filed_data_vec[ii];
498  nodes_dofs[pos] = (*dit).get();
499  ++pos;
500  ++dit;
501  }
502  }
503 
504  for (auto &dof_ptr : brother_dofs_vec) {
505  if (const auto d = dof_ptr.lock()) {
506  const auto &sn = d->sideNumberPtr;
507  const int side_number = sn->side_number;
508  const int brother_side_number = sn->brother_side_number;
509  bb_node_order[brother_side_number] = bb_node_order[side_number];
510  int pos = side_number * nb_dof_idx;
511  int brother_pos = brother_side_number * nb_dof_idx;
512  for (int ii = 0; ii != nb_dof_idx; ++ii) {
513  nodes_data[brother_pos] = nodes_data[pos];
514  nodes_dofs[brother_pos] = nodes_dofs[pos];
515  ++pos;
516  ++brother_pos;
517  }
518  }
519  }
520 
521  } else {
522  nodes_data.resize(0, false);
523  nodes_dofs.resize(0, false);
524  }
525 
527  };
528 
529  return get_nodes_field_data(
530  const_cast<FEDofEntity_multiIndex &>(
531  numeredEntFiniteElementPtr->getDataDofs()),
532  data.dataOnEntities[MBVERTEX][0].getFieldData(),
533  data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
534  data.dataOnEntities[MBVERTEX][0].getSpace(),
535  data.dataOnEntities[MBVERTEX][0].getBase(),
536  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder());
537 }
538 
540  DataForcesAndSourcesCore &data, const std::string &field_name,
541  const EntityType type_lo, const EntityType type_hi) const {
543  for (EntityType t = type_lo; t != type_hi; ++t) {
544  for (auto &dat : data.dataOnEntities[t]) {
545  dat.getDataOrder() = 0;
546  dat.getBase() = NOBASE;
547  dat.getSpace() = NOSPACE;
548  dat.getFieldData().resize(0, false);
549  dat.getFieldDofs().resize(0, false);
550  }
551  }
552 
553  auto &dofs = const_cast<FEDofEntity_multiIndex &>(
554  numeredEntFiniteElementPtr->getDataDofs());
555  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
556  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
557  if (dit == dofs_by_type.end())
559  auto hi_dit =
560  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
561  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
562  for (; dit != hi_dit;) {
563 
564  auto &dof = **dit;
565  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
566  if (nb_dofs_on_ent) {
567 
568  const EntityType type = dof.getEntType();
569  const int side = dof.sideNumberPtr->side_number;
570  if (side >= 0) {
571 
572  auto &dat = data.dataOnEntities[type][side];
573  auto &ent_field_dofs = dat.getFieldDofs();
574  auto &ent_field_data = dat.getFieldData();
575  const int brother_side = dof.sideNumberPtr->brother_side_number;
576  if (brother_side != -1)
577  brother_dofs_vec.emplace_back(*dit);
578 
579  if (ent_field_data.empty()) {
580 
581  dat.getBase() = dof.getApproxBase();
582  dat.getSpace() = dof.getSpace();
583  const int ent_order = dof.getMaxOrder();
584  dat.getDataOrder() =
585  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
586  ent_field_data.resize(nb_dofs_on_ent, false);
587  noalias(ent_field_data) = dof.getEntFieldData();
588  ent_field_dofs.resize(nb_dofs_on_ent, false);
589  for (int ii = 0; ii != nb_dofs_on_ent; ++ii) {
590  ent_field_dofs[ii] = (*dit).get();
591  ++dit;
592  }
593  }
594 
595  } else {
596 
597  for (int ii = 0; ii != nb_dofs_on_ent; ++ii)
598  ++dit;
599  }
600  }
601  }
602 
603  for (auto &dof_ptr : brother_dofs_vec) {
604  if (auto d = dof_ptr.lock()) {
605  const EntityType type = d->getEntType();
606  const int side = d->sideNumberPtr->side_number;
607  const int brother_side = d->sideNumberPtr->brother_side_number;
608  auto &dat = data.dataOnEntities[type][side];
609  auto &dat_brother = data.dataOnEntities[type][brother_side];
610  dat_brother.getBase() = dat.getBase();
611  dat_brother.getSpace() = dat.getSpace();
612  dat_brother.getDataOrder() = dat.getDataOrder();
613  dat_brother.getFieldData() = dat.getFieldData();
614  dat_brother.getFieldDofs() = dat.getFieldDofs();
615  }
616  }
617 
619 }
620 
622  const boost::string_ref field_name, FEDofEntity_multiIndex &dofs,
623  VectorDouble &ent_field_data, VectorDofs &ent_field_dofs) const {
625  auto dit = dofs.get<FieldName_mi_tag>().lower_bound(field_name);
626  auto hi_dit = dofs.get<FieldName_mi_tag>().upper_bound(field_name);
627  int size = std::distance(dit, hi_dit);
628  ent_field_data.resize(size, false);
629  ent_field_dofs.resize(size, false);
630  for (; dit != hi_dit; dit++) {
631  int idx = (*dit)->getDofCoeffIdx();
632  ent_field_data[idx] = (*dit)->getFieldData();
633  ent_field_dofs[idx] = (*dit).get();
634  }
636 }
637 
639  DataForcesAndSourcesCore &data, const boost::string_ref field_name) const {
641  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
642  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
643  }
645  field_name,
646  const_cast<FEDofEntity_multiIndex &>(
647  numeredEntFiniteElementPtr->getDataDofs()),
648  data.dataOnEntities[MBENTITYSET][0].getFieldData(),
649  data.dataOnEntities[MBENTITYSET][0].getFieldDofs());
651 }
652 
653 // ** Face **
654 
658  // PetscAttachDebugger();
659  data.facesNodes.resize(4, 3, false);
660  auto &side_table = const_cast<SideNumber_multiIndex &>(
661  numeredEntFiniteElementPtr->getSideNumberTable());
662  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBTRI, 0));
663  auto hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBTRI, 4));
664  if (std::distance(siit, hi_siit) != 4) {
665  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
666  "Should be 4 triangles on tet, side_table not initialized");
667  }
668  const int canonical_face_sense_p1[4][3] = {
669  {0, 1, 3},
670  {1, 2, 3},
671  {0, 3, 2} /**/,
672  {0, 2, 1} /**/}; // second index is offset (positive sense)
673  const int canonical_face_sense_m1[4][3] = {
674  {0, 3, 1},
675  {1, 3, 2},
676  {0, 2, 3},
677  {0, 1, 2}}; // second index is offset (negative sense
678  for (; siit != hi_siit; siit++) {
679  const boost::shared_ptr<SideNumber> side = *siit;
680  int face_conn[3] = {-1, -1, -1};
681  if (side->offset == 0) {
682  face_conn[0] = side->sense == 1
683  ? canonical_face_sense_p1[(int)side->side_number][0]
684  : canonical_face_sense_m1[(int)side->side_number][0];
685  face_conn[1] = side->sense == 1
686  ? canonical_face_sense_p1[(int)side->side_number][1]
687  : canonical_face_sense_m1[(int)side->side_number][1];
688  face_conn[2] = side->sense == 1
689  ? canonical_face_sense_p1[(int)side->side_number][2]
690  : canonical_face_sense_m1[(int)side->side_number][2];
691  }
692  if (side->offset == 1) {
693  face_conn[0] =
694  side->sense == 1
695  ? canonical_face_sense_p1[(int)side->side_number][1]
696  : canonical_face_sense_m1[(int)side->side_number][2] /**/;
697  face_conn[1] = side->sense == 1
698  ? canonical_face_sense_p1[(int)side->side_number][2]
699  : canonical_face_sense_m1[(int)side->side_number][0];
700  face_conn[2] = side->sense == 1
701  ? canonical_face_sense_p1[(int)side->side_number][0]
702  : canonical_face_sense_m1[(int)side->side_number][1];
703  }
704  if (side->offset == 2) {
705  face_conn[0] =
706  side->sense == 1
707  ? canonical_face_sense_p1[(int)side->side_number][2]
708  : canonical_face_sense_m1[(int)side->side_number][1] /**/;
709  face_conn[1] = side->sense == 1
710  ? canonical_face_sense_p1[(int)side->side_number][0]
711  : canonical_face_sense_m1[(int)side->side_number][2];
712  face_conn[2] = side->sense == 1
713  ? canonical_face_sense_p1[(int)side->side_number][1]
714  : canonical_face_sense_m1[(int)side->side_number][0];
715  }
716  for (int nn = 0; nn < 3; nn++)
717  data.facesNodes(side->side_number, nn) = face_conn[nn];
718  {
719  const EntityHandle *conn_tet;
720  int num_nodes_tet;
722  CHKERR mField.get_moab().get_connectivity(ent, conn_tet, num_nodes_tet,
723  true);
724  if (num_nodes_tet != 4)
725  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
726  "data inconsistency");
727  int num_nodes_face;
728  const EntityHandle *conn_face;
729  CHKERR mField.get_moab().get_connectivity(side->ent, conn_face,
730  num_nodes_face, true);
731  if (num_nodes_face != 3)
732  SETERRQ(PETSC_COMM_SELF, 1, "data inconsistency");
733  if (conn_face[0] != conn_tet[data.facesNodes(side->side_number, 0)])
734  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
735  "data inconsistency");
736  if (conn_face[1] != conn_tet[data.facesNodes(side->side_number, 1)])
737  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
738  "data inconsistency");
739  if (conn_face[2] != conn_tet[data.facesNodes(side->side_number, 2)])
740  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
741  "data inconsistency");
742  }
743  }
745 }
746 
747 // ** Space and Base **
748 
750  DataForcesAndSourcesCore &data) const {
752 
753  if (nInTheLoop == 0) {
754  data.sPace.reset();
755  data.bAse.reset();
756  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
757  data.spacesOnEntities[t].reset();
758  data.basesOnEntities[t].reset();
759  }
760  for (int s = 0; s != LASTSPACE; ++s) {
761  data.basesOnSpaces[s].reset();
762  }
763  }
764 
765  if (dataFieldEntsPtr)
766  for (auto e : *dataFieldEntsPtr) {
767  // get data from entity
768  const EntityType type = e->getEntType();
769  const FieldSpace space = e->getSpace();
770  const FieldApproximationBase approx = e->getApproxBase();
771 
772  // set data
773  data.sPace.set(space);
774  data.bAse.set(approx);
775  data.spacesOnEntities[type].set(space);
776  data.basesOnEntities[type].set(approx);
777  data.basesOnSpaces[space].set(approx);
778  }
779  else
780  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
781  "data fields ents not allocated on element");
782 
784 }
785 
787  const FieldApproximationBase b) {
789  if (dataOnElement[H1]->bAse.test(b)) {
790  switch (static_cast<FieldApproximationBase>(b)) {
791  case NOBASE:
792  break;
794  // BERNSTEIN_BEZIER_BASE is not hierarchical base
795  break;
800  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
801  "Functions genrating approximation base not defined");
802 
803  for (int space = H1; space != LASTSPACE; ++space) {
804  if (dataOnElement[H1]->sPace.test(space) &&
805  dataOnElement[H1]->bAse.test(b) &&
806  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
807  CHKERR getElementPolynomialBase()->getValue(
808  gaussPts,
809  boost::make_shared<EntPolynomialBaseCtx>(
810  *dataOnElement[space], static_cast<FieldSpace>(space),
811  static_cast<FieldApproximationBase>(b), NOBASE));
812  }
813  }
814  break;
815  case USER_BASE:
816  if (!getUserPolynomialBase())
817  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
818  "Functions genrating user approximation base not defined");
819 
820  for (int space = H1; space != LASTSPACE; ++space)
821  if (dataOnElement[H1]->sPace.test(space) &&
822  dataOnElement[H1]->bAse.test(b) &&
823  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
824  CHKERR getUserPolynomialBase()->getValue(
825  gaussPts,
826  boost::make_shared<EntPolynomialBaseCtx>(
827  *dataOnElement[space], static_cast<FieldSpace>(space),
828  static_cast<FieldApproximationBase>(b), NOBASE));
829  }
830  break;
831  default:
833  "Base <%s> not yet implemented",
834  ApproximationBaseNames[static_cast<FieldApproximationBase>(b)]);
835  }
836  }
838 }
839 
842  /// Use the some node base. Node base is usually used for construction other
843  /// bases.
844  for (int space = HCURL; space != LASTSPACE; ++space) {
845  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
846  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
847  }
848  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
850  static_cast<FieldApproximationBase>(b));
851  }
853 }
854 
858 
859  auto get_nodal_base_data = [&](DataForcesAndSourcesCore &data,
860  const std::string &field_name) {
862  auto &space = data.dataOnEntities[MBVERTEX][0].getSpace();
863  auto &base = data.dataOnEntities[MBVERTEX][0].getBase();
864  auto &bb_node_order = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder();
865 
866  auto &dofs_by_name_and_type =
868  auto tuple = boost::make_tuple(field_name, MBVERTEX);
869  auto dit = dofs_by_name_and_type.lower_bound(tuple);
870  if (dit == dofs_by_name_and_type.end())
871  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
872  "No nodal dofs on element");
873  auto hi_dit =
874  dataPtr->get<Composite_Name_And_Type_mi_tag>().upper_bound(tuple);
875 
876  if (dit != hi_dit) {
877  auto &first_dof = **dit;
878  space = first_dof.getSpace();
879  base = first_dof.getApproxBase();
880  int num_nodes;
881  CHKERR getNumberOfNodes(num_nodes);
882  bb_node_order.resize(num_nodes, false);
883  bb_node_order.clear();
884  const int nb_dof_idx = first_dof.getNbOfCoeffs();
885 
886  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
887  for (; dit != hi_dit;) {
888  const auto &dof_ptr = *dit;
889  const auto &dof = *dof_ptr;
890  const auto &sn = *dof.sideNumberPtr;
891  const int side_number = sn.side_number;
892  const int brother_side_number = sn.brother_side_number;
893  if (brother_side_number != -1)
894  brother_dofs_vec.emplace_back(dof_ptr);
895  bb_node_order[side_number] = dof.getMaxOrder();
896  for (int ii = 0; ii != nb_dof_idx; ++ii)
897  ++dit;
898  }
899 
900  for (auto &dof_ptr : brother_dofs_vec) {
901  if (const auto d = dof_ptr.lock()) {
902  const auto &sn = d->sideNumberPtr;
903  const int side_number = sn->side_number;
904  const int brother_side_number = sn->brother_side_number;
905  bb_node_order[brother_side_number] = bb_node_order[side_number];
906  }
907  }
908  }
910  };
911 
912  auto get_entity_base_data = [&](DataForcesAndSourcesCore &data,
913  const std::string &field_name,
914  const EntityType type_lo,
915  const EntityType type_hi) {
917  for (EntityType t = type_lo; t != type_hi; ++t) {
918  for (auto &dat : data.dataOnEntities[t]) {
919  dat.getDataOrder() = 0;
920  dat.getBase() = NOBASE;
921  dat.getSpace() = NOSPACE;
922  dat.getFieldData().resize(0, false);
923  dat.getFieldDofs().resize(0, false);
924  }
925  }
926 
927  auto &dofs = const_cast<FEDofEntity_multiIndex &>(
928  numeredEntFiniteElementPtr->getDataDofs());
929  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
930  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
931  if (dit == dofs_by_type.end())
933  auto hi_dit =
934  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
935  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
936  for (; dit != hi_dit;) {
937 
938  auto &dof = **dit;
939  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
940  if (nb_dofs_on_ent) {
941  const EntityType type = dof.getEntType();
942  const int side = dof.sideNumberPtr->side_number;
943  auto &dat = data.dataOnEntities[type][side];
944 
945  const int brother_side = dof.sideNumberPtr->brother_side_number;
946  if (brother_side != -1)
947  brother_dofs_vec.emplace_back(*dit);
948 
949  dat.getBase() = dof.getApproxBase();
950  dat.getSpace() = dof.getSpace();
951  const int ent_order = dof.getMaxOrder();
952  dat.getDataOrder() =
953  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
954  for (int ii = 0; ii != nb_dofs_on_ent; ++ii) {
955  ++dit;
956  }
957  }
958  }
959 
960  for (auto &dof_ptr : brother_dofs_vec) {
961  if (auto d = dof_ptr.lock()) {
962  const EntityType type = d->getEntType();
963  const int side = d->sideNumberPtr->side_number;
964  const int brother_side = d->sideNumberPtr->brother_side_number;
965  auto &dat = data.dataOnEntities[type][side];
966  auto &dat_brother = data.dataOnEntities[type][brother_side];
967  dat_brother.getBase() = dat.getBase();
968  dat_brother.getSpace() = dat.getSpace();
969  dat_brother.getDataOrder() = dat.getDataOrder();
970  }
971  }
973  };
974 
975  auto &ents_data = *dataFieldEntsPtr;
976  for (auto &e : ents_data) {
977  if (e->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
978  auto space = e->getSpace();
979  for (EntityType t = MBVERTEX; t != MBPOLYHEDRON; ++t) {
980  for (auto &dat : (*dataOnElement[space]).dataOnEntities[t]) {
981  for (auto &ptr : dat.getBBAlphaIndicesByOrderArray())
982  ptr.reset();
983  for (auto &ptr : dat.getBBNByOrderArray())
984  ptr.reset();
985  for (auto &ptr : dat.getBBDiffNByOrderArray())
986  ptr.reset();
987  }
988  }
989  }
990  }
991 
992  for (auto &e : ents_data) {
993  if (e->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
994  auto field_name = e->getName();
995  auto space = e->getSpace();
996  CHKERR get_nodal_base_data(*dataOnElement[space], field_name);
997  CHKERR get_entity_base_data(*dataOnElement[space], field_name, MBEDGE,
998  MBPOLYHEDRON);
999  CHKERR getElementPolynomialBase()->getValue(
1000  gaussPts,
1001  boost::make_shared<EntPolynomialBaseCtx>(
1002  *dataOnElement[space], field_name, static_cast<FieldSpace>(space),
1004  }
1005  }
1006 
1008 };
1009 
1012 
1013  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1016 
1017  // Data on elements for proper spaces
1018  for (int space = H1; space != LASTSPACE; ++space) {
1019  dataOnElement[space]->setElementType(type);
1020  derivedDataOnElement[space]->setElementType(type);
1021  }
1022 
1024 
1026 }
1027 
1028 #define FUNCTION_NAME_WITH_OP_NAME(OP) \
1029  std::ostringstream ss; \
1030  ss << "(Calling user data operator " \
1031  << boost::typeindex::type_id_runtime(OP).pretty_name() << " rowField " \
1032  << (OP).rowFieldName << " colField " << (OP).colFieldName << ") "
1033 
1034 #define CATCH_OP_ERRORS(OP) \
1035  catch (MoFEMExceptionInitial const &ex) { \
1036  FUNCTION_NAME_WITH_OP_NAME(OP) << PETSC_FUNCTION_NAME; \
1037  return PetscError(PETSC_COMM_SELF, ex.lINE, ss.str().c_str(), __FILE__, \
1038  ex.errorCode, PETSC_ERROR_INITIAL, ex.what()); \
1039  } \
1040  catch (MoFEMExceptionRepeat const &ex) { \
1041  FUNCTION_NAME_WITH_OP_NAME(OP) << PETSC_FUNCTION_NAME; \
1042  return PetscError(PETSC_COMM_SELF, ex.lINE, ss.str().c_str(), __FILE__, \
1043  ex.errorCode, PETSC_ERROR_REPEAT, " "); \
1044  } \
1045  catch (MoFEMException const &ex) { \
1046  FUNCTION_NAME_WITH_OP_NAME(OP) << ex.errorMessage; \
1047  SETERRQ(PETSC_COMM_SELF, ex.errorCode, ss.str().c_str()); \
1048  } \
1049  catch (std::exception const &ex) { \
1050  std::string message("Error: " + std::string(ex.what()) + " at " + \
1051  boost::lexical_cast<std::string>(__LINE__) + " : " + \
1052  std::string(__FILE__) + " in " + \
1053  std::string(PETSC_FUNCTION_NAME)); \
1054  FUNCTION_NAME_WITH_OP_NAME(OP) << message; \
1055  SETERRQ(PETSC_COMM_SELF, MOFEM_STD_EXCEPTION_THROW, ss.str().c_str()); \
1056  }
1057 
1060 
1061  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1064  std::vector<std::string> last_eval_field_name(2);
1065 
1066  boost::ptr_vector<UserDataOperator>::iterator oit, hi_oit;
1067  oit = opPtrVector.begin();
1068  hi_oit = opPtrVector.end();
1069 
1070  for (; oit != hi_oit; oit++) {
1071 
1072  try {
1073 
1074  CHKERR oit->setPtrFE(this);
1075 
1076  if (oit->opType == UserDataOperator::OPLAST) {
1077 
1078  // Set field
1079  switch (oit->sPace) {
1080  case NOSPACE:
1081  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Unknown space");
1082  case NOFIELD:
1083  case H1:
1084  case HCURL:
1085  case HDIV:
1086  case L2:
1087  break;
1088  default:
1089  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1090  "Not implemented for this space", oit->sPace);
1091  }
1092 
1093  // Reseat all data which all field dependent
1094  dataOnElement[oit->sPace]->resetFieldDependentData();
1095  last_eval_field_name[0] = "";
1096 
1097  // Run operator
1098  try {
1099  CHKERR oit->opRhs(*dataOnElement[oit->sPace], false);
1100  }
1101  CATCH_OP_ERRORS(*oit);
1102 
1103  } else {
1104 
1105  boost::shared_ptr<DataForcesAndSourcesCore> op_data[2];
1106  std::array<bool, 2> base_swap;
1107  std::array<std::pair<std::string, FieldApproximationBase>, 2>
1108  base_swap_data;
1109  auto swap_bases = [&]() {
1111  for (size_t ss = 0; ss != 2; ++ss)
1112  if (base_swap[ss])
1113  CHKERR op_data[ss]->baseSwap(base_swap_data[ss].first,
1114  base_swap_data[ss].second);
1116  };
1117 
1118  for (size_t ss = 0; ss != 2; ss++) {
1119 
1120  const std::string field_name =
1121  !ss ? oit->rowFieldName : oit->colFieldName;
1122  if (field_name.empty()) {
1123  SETERRQ2(
1124  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1125  "No field name in operator %d (0-row, 1-column) in operator %s",
1126  ss,
1127  (boost::typeindex::type_id_runtime(*oit).pretty_name())
1128  .c_str());
1129  }
1130  const Field *field_struture = mField.get_field_structure(field_name);
1131  const BitFieldId data_id = field_struture->getId();
1132  const FieldSpace space = field_struture->getSpace();
1133  const FieldApproximationBase base = field_struture->getApproxBase();
1134  op_data[ss] =
1135  !ss ? dataOnElement[space] : derivedDataOnElement[space];
1136 
1137  switch (base) {
1139  base_swap_data[ss] = std::pair<std::string, FieldApproximationBase>(
1140  field_name, AINSWORTH_BERNSTEIN_BEZIER_BASE);
1141  base_swap[ss] = true;
1142  break;
1143  default:
1144  base_swap[ss] = false;
1145  };
1146 
1147  if ((oit->getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
1148  data_id)
1149  .none()) {
1150  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1151  "no data field < %s > on finite element < %s >",
1152  field_name.c_str(), feName.c_str());
1153  }
1154 
1155  if (oit->getOpType() & types[ss] ||
1156  oit->getOpType() & UserDataOperator::OPROWCOL) {
1157 
1158  switch (space) {
1159  case NOSPACE:
1160  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1161  "unknown space");
1162  break;
1163  case NOFIELD:
1164  case H1:
1165  case HCURL:
1166  case HDIV:
1167  case L2:
1168  break;
1169  default:
1170  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1171  "Not implemented for this space", space);
1172  }
1173 
1174  if (last_eval_field_name[ss] != field_name) {
1175 
1176  CHKERR getEntityFieldData(*op_data[ss], field_name, MBEDGE);
1177  if (!ss)
1178  CHKERR getEntityRowIndices(*op_data[ss], field_name, MBEDGE);
1179  else
1180  CHKERR getEntityColIndices(*op_data[ss], field_name, MBEDGE);
1181 
1182  switch (space) {
1183  case NOSPACE:
1184  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1185  "unknown space");
1186  break;
1187  case H1:
1188  if (!ss)
1189  CHKERR getRowNodesIndices(*op_data[ss], field_name);
1190  else
1191  CHKERR getColNodesIndices(*op_data[ss], field_name);
1192  CHKERR getNodesFieldData(*op_data[ss], field_name);
1193  break;
1194  case HCURL:
1195  case HDIV:
1196  break;
1197  case L2:
1198  switch (type) {
1199  case MBVERTEX:
1200  CHKERR getNodesFieldData(*op_data[ss], field_name);
1201  break;
1202  default:
1203  break;
1204  }
1205  break;
1206  case NOFIELD:
1207  if (!getNinTheLoop()) {
1208  // NOFIELD data are the same for each element, can be
1209  // retrieved only once
1210  if (!ss) {
1211  CHKERR getNoFieldRowIndices(*op_data[ss], field_name);
1212  } else {
1213  CHKERR getNoFieldColIndices(*op_data[ss], field_name);
1214  }
1215  CHKERR getNoFieldFieldData(*op_data[ss], field_name);
1216  }
1217  break;
1218  default:
1219  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1220  "Not implemented for this space", space);
1221  }
1222  last_eval_field_name[ss] = field_name;
1223  }
1224  }
1225  }
1226 
1227  CHKERR swap_bases();
1228 
1229  if (oit->getOpType() & UserDataOperator::OPROW) {
1230  try {
1231  CHKERR oit->opRhs(*op_data[0], false);
1232  }
1233  CATCH_OP_ERRORS(*oit);
1234  }
1235 
1236  if (oit->getOpType() & UserDataOperator::OPCOL) {
1237  try {
1238  CHKERR oit->opRhs(*op_data[1], false);
1239  }
1240  CATCH_OP_ERRORS(*oit);
1241  }
1242 
1243  if (oit->getOpType() & UserDataOperator::OPROWCOL) {
1244  try {
1245  CHKERR oit->opLhs(*op_data[0], *op_data[1]);
1246  }
1247  CATCH_OP_ERRORS(*oit);
1248  }
1249 
1250  CHKERR swap_bases();
1251 
1252  }
1253  }
1254  CATCH_OP_ERRORS(*oit);
1255  }
1257 }
1258 
1260  const std::string field_name, const EntityType type, const int side,
1261  VectorInt &indices) const {
1263  if (ptrFE == NULL)
1264  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1265 
1266  switch (type) {
1267  case MBVERTEX:
1268  CHKERR ptrFE->getProblemNodesRowIndices(field_name, indices);
1269  break;
1270  default:
1271  CHKERR ptrFE->getProblemTypeRowIndices(field_name, type, side, indices);
1272  }
1274 }
1275 
1277  const std::string field_name, const EntityType type, const int side,
1278  VectorInt &indices) const {
1280  if (ptrFE == NULL)
1281  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1282 
1283  switch (type) {
1284  case MBVERTEX:
1285  CHKERR ptrFE->getProblemNodesColIndices(field_name, indices);
1286  break;
1287  default:
1288  CHKERR ptrFE->getProblemTypeColIndices(field_name, type, side, indices);
1289  }
1291 }
1292 
1296  sidePtrFE = const_cast<ForcesAndSourcesCore *>(side_fe_ptr);
1298 }
1299 
1301  const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t side_dim,
1302  const EntityHandle ent_for_side) {
1304  const EntityHandle ent = ent_for_side ? ent_for_side : getFEEntityHandle();
1305 
1306  const Problem *problem_ptr = getFEMethod()->problemPtr;
1307  Range adjacent_ents;
1308  CHKERR ptrFE->mField.getInterface<BitRefManager>()->getAdjacenciesAny(
1309  ent, side_dim, adjacent_ents);
1310  typedef NumeredEntFiniteElement_multiIndex::index<
1311  Composite_Name_And_Ent_mi_tag>::type FEByComposite;
1312  FEByComposite &numered_fe =
1314 
1315  side_fe->feName = fe_name;
1316 
1317  CHKERR side_fe->setSideFEPtr(ptrFE);
1318  CHKERR side_fe->copyBasicMethod(*getFEMethod());
1319  CHKERR side_fe->copyKsp(*getFEMethod());
1320  CHKERR side_fe->copySnes(*getFEMethod());
1321  CHKERR side_fe->copyTs(*getFEMethod());
1322 
1323  CHKERR side_fe->preProcess();
1324 
1325  int nn = 0;
1326  side_fe->loopSize = adjacent_ents.size();
1327  for (Range::iterator vit = adjacent_ents.begin(); vit != adjacent_ents.end();
1328  vit++) {
1329  FEByComposite::iterator miit =
1330  numered_fe.find(boost::make_tuple(fe_name, *vit));
1331  if (miit != numered_fe.end()) {
1332  side_fe->nInTheLoop = nn++;
1333  side_fe->numeredEntFiniteElementPtr = *miit;
1334  side_fe->dataFieldEntsPtr = (*miit)->sPtr->data_field_ents_view;
1335  side_fe->rowFieldEntsPtr = (*miit)->sPtr->row_field_ents_view;
1336  side_fe->colFieldEntsPtr = (*miit)->sPtr->col_field_ents_view;
1337  side_fe->dataPtr = (*miit)->sPtr->data_dofs;
1338  side_fe->rowPtr = (*miit)->rows_dofs;
1339  side_fe->colPtr = (*miit)->cols_dofs;
1340  CHKERR (*side_fe)();
1341  }
1342  }
1343 
1344  CHKERR side_fe->postProcess();
1345 
1347 }
1348 
1349 int ForcesAndSourcesCore::getRule(int order_row, int order_col,
1350  int order_data) {
1351  return getRuleHook ? getRuleHook(order_row, order_col, order_data)
1352  : getRule(order_data);
1353 }
1354 
1356  int order_data) {
1357  return setRuleHook ? setRuleHook(order_row, order_col, order_data)
1358  : setGaussPts(order_data);
1359 }
1360 
1362 
1363 /** \deprecated setGaussPts(int row_order, int col_order, int data order);
1364  */
1367  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Sorry, not implemented");
1369 }
1370 
1372  const char type,
1373  const bool symm)
1374  : DataOperator(symm), opType(type), sPace(space), ptrFE(nullptr) {}
1375 
1377  const std::string &field_name, const char type, const bool symm)
1378  : DataOperator(symm), opType(type), rowFieldName(field_name),
1379  colFieldName(field_name), sPace(LASTSPACE), ptrFE(nullptr) {}
1380 
1382  const std::string &row_field_name, const std::string &col_field_name,
1383  const char type, const bool symm)
1384  : DataOperator(symm), opType(type), rowFieldName(row_field_name),
1385  colFieldName(col_field_name), sPace(LASTSPACE), ptrFE(nullptr) {}
1386 
1389  if (preProcessHook) {
1390  ierr = preProcessHook();
1391  CHKERRG(ierr);
1392  }
1394 }
1397  if (operatorHook) {
1398  ierr = operatorHook();
1399  CHKERRG(ierr);
1400  }
1402 }
1405  if (postProcessHook) {
1406  ierr = postProcessHook();
1407  CHKERRG(ierr);
1408  }
1410 }
1411 
1415  if (!(ptrFE = dynamic_cast<ForcesAndSourcesCore *>(ptr)))
1416  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1417  "User operator and finite element do not work together");
1419 }
1420 
1421 } // namespace MoFEM
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from other base method to this base method.
ForcesAndSourcesCore(Interface &m_field)
int getMaxRowOrder() const
Get max order of approximation for field in rows.
Deprecated interface functions.
boost::ptr_vector< UserDataOperator > opPtrVector
Vector of finite element users data operators.
MoFEMErrorCode getProblemColIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get col indices.
MatrixDouble gaussPts
Matrix of integration points.
base operator to do operations at Gauss Pt. level
field with continuous normal traction
Definition: definitions.h:179
user implemented approximation base
Definition: definitions.h:160
MoFEMErrorCode getNoFieldColIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
virtual moab::Interface & get_moab()=0
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
MoFEMErrorCode getEntitySense(const EntityType type, boost::ptr_vector< DataForcesAndSourcesCore::EntData > &data) const
get sense (orientation) of entity
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
MoFEMErrorCode createDataOnElement()
Create a entity data on element object.
Provide data structure for (tensor) field approximation.The Field is intended to provide support for ...
MoFEMErrorCode getColNodesIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col node indices from FENumeredDofEntity_multiIndex
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement()
Calculate Bernstein-Bezier base.
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:37
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
int loopSize
local number oe methods to process
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
MoFEMErrorCode getEntityDataOrder(const EntityType type, const FieldSpace space, boost::ptr_vector< DataForcesAndSourcesCore::EntData > &data) const
Get the entity data order.
auto & getUserPolynomialBase()
Get the User Polynomial Base object.
multi_index_container< boost::shared_ptr< FENumeredDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< FENumeredDofEntity::interface_type_DofEntity, const UId &, &FENumeredDofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FENumeredDofEntity::interface_type_DofEntity, EntityHandle, &FENumeredDofEntity::getEnt > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef > >, ordered_non_unique< tag< Composite_Name_Type_And_Side_Number_mi_tag >, composite_key< FENumeredDofEntity, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef >, const_mem_fun< FENumeredDofEntity::interface_type_RefEntity, EntityType, &FENumeredDofEntity::getEntType >, KeyFromKey< member< SideNumber, char, &SideNumber::side_number >, member< FENumeredDofEntity::BaseFEEntity, boost::shared_ptr< SideNumber >, &FENumeredDofEntity::sideNumberPtr > > > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< FENumeredDofEntity, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef >, const_mem_fun< FENumeredDofEntity::interface_type_RefEntity, EntityType, &FENumeredDofEntity::getEntType > > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FENumeredDofEntity, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef >, const_mem_fun< FENumeredDofEntity::interface_type_DofEntity, EntityHandle, &FENumeredDofEntity::getEnt > > > > > FENumeredDofEntity_multiIndex
MultiIndex container keeps FENumeredDofEntity.
boost::shared_ptr< const FieldEntity_vector_view > colFieldEntsPtr
Pointer to finite element field entities column view.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
std::string feName
Name of finite element.
MoFEMErrorCode getFaceTriNodes(DataForcesAndSourcesCore &data) const
Get nodes on triangles.
FieldCoefficientsNumber getNbOfCoeffs() const
Get number of field coefficients.
int getNinTheLoop() const
get number of evaluated element in the loop
boost::shared_ptr< const FieldEntity_vector_view > rowFieldEntsPtr
Pointer to finite element field entities row view.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
keeps basic data about problemThis is low level structure with information about problem,...
MoFEMErrorCode getNoFieldRowIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
UserDataOperator(const FieldSpace space, const char type=OPLAST, const bool symm=true)
boost::shared_ptr< const FEDofEntity_multiIndex > dataPtr
Pointer to finite element data dofs.
std::array< std::bitset< LASTBASE >, LASTSPACE > basesOnSpaces
base on spaces
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
FieldSpace getSpace() const
Get field approximation space.
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
OpType
Controls loop over entities on element.
MoFEMErrorCode getNoFieldFieldData(const boost::string_ref field_name, FEDofEntity_multiIndex &dofs, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs) const
Get field data on nodes.
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
const Problem * problemPtr
raw pointer to problem
static const char *const ApproximationBaseNames[]
Definition: definitions.h:164
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode copyKsp(const KspMethod &ksp)
copy data form another method
Definition: LoopMethods.cpp:60
MoFEMErrorCode setSideFEPtr(const ForcesAndSourcesCore *side_fe_ptr)
Set the pointer to face element on the side.
MoFEMErrorCode getProblemNodesRowIndices(const std::string &field_name, VectorInt &nodes_indices) const
MoFEMErrorCode getNodesFieldData(DataForcesAndSourcesCore &data, const std::string &field_name) const
Get data on nodes.
RuleHookFun getRuleHook
Hook to get rule.
#define CATCH_OP_ERRORS(OP)
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data.
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
MoFEMErrorCode getProblemTypeRowIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
MatrixInt facesNodes
nodes on finite element faces
FieldApproximationBase
approximation base
Definition: definitions.h:150
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
const BitFieldId & getId() const
Get unique field id.
MoFEMErrorCode getNodesIndices(const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices, VectorInt &local_nodes_indices) const
get node indices
MoFEMErrorCode getProblemNodesColIndices(const std::string &field_name, VectorInt &nodes_indices) const
FieldApproximationBase getApproxBase() const
Get approximation base.
Managing BitRefLevels.
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
MoFEMErrorCode getProblemTypeColIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
EntityHandle getFEEntityHandle() const
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsRows
store DOFs on rows for this problem
field with continuous tangents
Definition: definitions.h:178
FieldSpace
approximation spaces
Definition: definitions.h:174
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:52
#define CHKERR
Inline error check.
Definition: definitions.h:602
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &data, const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
int nInTheLoop
number currently of processed method
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
MoFEMErrorCode getEntityFieldData(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
std::array< std::bitset< LASTBASE >, MBMAXTYPE > basesOnEntities
bases on entity types
constexpr int order
MoFEMErrorCode getNoFieldIndices(const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
get NoField indices
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
static int getMaxOrder(const ENTMULTIINDEX &multi_index)
multi_index_container< boost::shared_ptr< FEDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< FEDofEntity::interface_type_DofEntity, const UId &, &FEDofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FEDofEntity::interface_type_DofEntity, EntityHandle, &FEDofEntity::getEnt > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType > > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_DofEntity, EntityHandle, &FEDofEntity::getEnt > > >, ordered_non_unique< tag< Composite_Name_Type_And_Side_Number_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType >, KeyFromKey< member< SideNumber, char, &SideNumber::side_number >, member< FEDofEntity::BaseFEEntity, boost::shared_ptr< SideNumber >, &FEDofEntity::sideNumberPtr > > > > > > FEDofEntity_multiIndex
MultiIndex container keeps FEDofEntity.
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:72
MultiIndex Tag for field name.
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy snes data.
Definition: LoopMethods.cpp:84
virtual MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr)
std::bitset< LASTSPACE > sPace
spaces on element
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
virtual MoFEMErrorCode operator()()
function is run for every finite element
MoFEMErrorCode getEntityColIndices(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
structure to get information form mofem into DataForcesAndSourcesCore
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
virtual MPI_Comm & get_comm() const =0
continuous field
Definition: definitions.h:177
virtual const Field * get_field_structure(const std::string &name)=0
get field structure
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:73
MoFEMErrorCode getRowNodesIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get row node indices from FENumeredDofEntity_multiIndex
boost::shared_ptr< NumeredEntFiniteElement_multiIndex > numeredFiniteElements
store finite elements
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
MoFEMErrorCode loopSide(const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0)
User call this function to loop over elements on the side of face. This function calls finite element...
MoFEMErrorCode getSpacesAndBaseOnEntities(DataForcesAndSourcesCore &data) const
Get field approximation space and base on entities.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
boost::shared_ptr< const FENumeredDofEntity_multiIndex > colPtr
Pointer to finite element columns dofs view.
std::bitset< LASTBASE > bAse
bases on element
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsCols
store DOFs on columns for this problem
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.
int getMaxColOrder() const
Get max order of approximation for field in columns.
MoFEMErrorCode getProblemRowIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get row indices.
RuleHookFun setRuleHook
Set function to calculate integration rule.
int getMaxDataOrder() const
Get max order of approximation for data fields.
boost::shared_ptr< const FENumeredDofEntity_multiIndex > rowPtr
Pointer to finite element rows dofs view.
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.
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
MoFEMErrorCode getEntityRowIndices(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, const UId &, &NumeredDofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Unique_Ent_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, const UId &, &NumeredDofEntity::getEntGlobalUniqueId > >, 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< FieldName_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_Field, boost::string_ref, &NumeredDofEntity::getNameRef > >, 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 > >, ordered_non_unique< tag< Composite_Name_And_Ent_And_EntDofIdx_mi_tag >, composite_key< NumeredDofEntity, const_mem_fun< NumeredDofEntity::interface_type_Field, boost::string_ref, &NumeredDofEntity::getNameRef >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, EntityHandle, &NumeredDofEntity::getEnt >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, DofIdx, &NumeredDofEntity::getEntDofIdx > > >, ordered_non_unique< tag< Composite_Name_And_Part_mi_tag >, composite_key< NumeredDofEntity, const_mem_fun< NumeredDofEntity::interface_type_Field, boost::string_ref, &NumeredDofEntity::getNameRef >, member< NumeredDofEntity, unsigned int, &NumeredDofEntity::pArt > > >, ordered_non_unique< tag< Composite_Name_Ent_And_Part_mi_tag >, composite_key< NumeredDofEntity, const_mem_fun< NumeredDofEntity::interface_type_Field, boost::string_ref, &NumeredDofEntity::getNameRef >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, EntityHandle, &NumeredDofEntity::getEnt >, member< NumeredDofEntity, unsigned int, &NumeredDofEntity::pArt > > > > > NumeredDofEntity_multiIndex
MultiIndex container keeps NumeredDofEntity.
field with C-1 continuity
Definition: definitions.h:180