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