v0.9.0
MatrixManager.cpp
Go to the documentation of this file.
1 /**
2  * \brief Create adjacent matrices using different indices
3 
4  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
5  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
6  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
7  * License for more details.
8  *
9  * You should have received a copy of the GNU Lesser General Public
10  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
11 
12 */
13 
14 namespace MoFEM {
15 
16 /** \brief Create compressed matrix
17 
18  \note Only function class members are allowed in this class. NO VARIABLES.
19 
20  \todo It is obsolete implementation, code should be moved to interface
21  MatrixManager.
22 
23  \todo While matrix is created is assumed that all entities on element are
24  adjacent to each other, in some cases, this is created denser matrices than it
25  should be. Some kind of filtering can be added.
26 
27  \todo Creation of the block matrices
28 
29  \todo Some efficiency improvemnt are possible
30 
31 
32  */
34 
36  MPI_Comm comm = PETSC_COMM_WORLD, int verbose = 1)
37  : Core(moab, comm, verbose) {}
38 
39  typedef FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
41  typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
42  typedef NumeredDofEntity_multiIndex::index<PetscGlobalIdx_mi_tag>::type
44 
45  /** \brief Create matrix adjacencies
46 
47  Depending on TAG type, which some structure used to number dofs, matrix is
48  partitioned using part number stored in multi-index, or is partitioned on
49  parts based only on index number.
50 
51  See: Idx_mi_tag PetscGlobalIdx_mi_tag and PetscLocalIdx_mi_tag
52 
53  */
54  template <typename TAG>
56  createMatArrays(ProblemsByName::iterator p_miit, const MatType type,
57  std::vector<PetscInt> &i, std::vector<PetscInt> &j,
58  const bool no_diagonals = true, int verb = QUIET) const;
59 
60  /** \brief Get element adjacencies
61  */
62  template <typename TAG>
64  ProblemsByName::iterator p_miit,
65  typename boost::multi_index::index<NumeredDofEntity_multiIndex,
66  TAG>::type::iterator mit_row,
67  boost::shared_ptr<FieldEntity> mofem_ent_ptr,
68  std::vector<int> &dofs_col_view, int verb) const;
69 
71  buildFECol(ProblemsByName::iterator p_miit,
72  boost::shared_ptr<EntFiniteElement> ent_fe_ptr, bool do_cols_prob,
73  boost::shared_ptr<NumeredEntFiniteElement> &fe_ptr) const;
74 };
75 
77  ProblemsByName::iterator p_miit,
78  boost::shared_ptr<EntFiniteElement> ent_fe_ptr, bool do_cols_prob,
79  boost::shared_ptr<NumeredEntFiniteElement> &fe_ptr) const {
81 
82  if (!ent_fe_ptr)
83  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
84  "Pointer to EntFiniteElement not given");
85 
86  // if element is not part of problem
87  if ((ent_fe_ptr->getId() & p_miit->getBitFEId()).none())
89 
90  BitRefLevel prb_bit = p_miit->getBitRefLevel();
91  BitRefLevel prb_mask = p_miit->getMaskBitRefLevel();
92  BitRefLevel fe_bit = ent_fe_ptr->getBitRefLevel();
93  // if entity is not problem refinement level
94  if ((fe_bit & prb_mask) != fe_bit)
96  if ((fe_bit & prb_bit) != prb_bit)
98 
99  auto fe_it =
100  p_miit->numeredFiniteElements->find(ent_fe_ptr->getGlobalUniqueId());
101 
102  // Create element if is not there
103  if (fe_it == p_miit->numeredFiniteElements->end()) {
104  auto p = p_miit->numeredFiniteElements->insert(
105  boost::make_shared<NumeredEntFiniteElement>(ent_fe_ptr));
106  fe_it = p.first;
107  }
108  fe_ptr = *fe_it;
109 
110  if (fe_ptr) {
111 
112  // Build DOFs on columns
113  if (fe_ptr->cols_dofs->empty()) {
114 
115  // Get dofs on columns
117  CHKERR fe_ptr->getEntFiniteElement()->getColDofView(
118  *(p_miit->numeredDofsCols), cols_view, moab::Interface::UNION);
119 
120  // Reserve memory for field dofs
121  boost::shared_ptr<std::vector<FENumeredDofEntity>> dofs_array(
122  new std::vector<FENumeredDofEntity>());
123  fe_ptr->getColDofsSequence() = dofs_array;
124  dofs_array->reserve(cols_view.size());
125 
126  // Create dofs objects
127  for (auto &dof : cols_view) {
128  auto side_number_ptr = fe_ptr->getSideNumberPtr(dof->getEnt());
129  dofs_array->emplace_back(side_number_ptr, dof);
130  }
131 
132  // Finally add DoFS to multi-indices
133  auto hint = fe_ptr->cols_dofs->end();
134  for (auto &dof : *dofs_array)
135  hint = fe_ptr->cols_dofs->emplace_hint(hint, dofs_array, &dof);
136  }
137 
138  } else
139  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
140  "At that point ptr to finite element should be well known");
141 
143 }
144 
145 template <typename TAG>
147  ProblemsByName::iterator p_miit,
148  typename boost::multi_index::index<NumeredDofEntity_multiIndex,
149  TAG>::type::iterator mit_row,
150  boost::shared_ptr<FieldEntity> mofem_ent_ptr,
151  std::vector<int> &dofs_col_view, int verb) const {
153 
154  BitRefLevel prb_bit = p_miit->getBitRefLevel();
155  BitRefLevel prb_mask = p_miit->getMaskBitRefLevel();
156 
157  // check if dofs and columns are the same, i.e. structurally symmetric problem
158  bool do_cols_prob;
159  if (p_miit->numeredDofsRows == p_miit->numeredDofsCols)
160  do_cols_prob = false;
161  else
162  do_cols_prob = true;
163 
164  dofs_col_view.clear();
165  for (auto r = entFEAdjacencies.get<Unique_mi_tag>().equal_range(
166  mofem_ent_ptr->getGlobalUniqueId());
167  r.first != r.second; ++r.first) {
168 
169  if (r.first->byWhat & BYROW) {
170 
171  if ((r.first->entFePtr->getId() & p_miit->getBitFEId()).none()) {
172  // if element is not part of problem
173  continue;
174  }
175 
176  BitRefLevel fe_bit = r.first->entFePtr->getBitRefLevel();
177  // if entity is not problem refinement level
178  if ((fe_bit & prb_mask) != fe_bit)
179  continue;
180  if ((fe_bit & prb_bit) != prb_bit)
181  continue;
182  BitRefLevel dof_bit = mit_row->get()->getBitRefLevel();
183  // if entity is not problem refinement level
184  if ((fe_bit & dof_bit).none())
185  continue;
186 
187  boost::shared_ptr<NumeredEntFiniteElement> fe_ptr;
188  // get element, if element is not in database build columns dofs
189  CHKERR buildFECol(p_miit, r.first->entFePtr, do_cols_prob, fe_ptr);
190 
191  if (fe_ptr) {
192  for (FENumeredDofEntity_multiIndex::iterator vit =
193  fe_ptr.get()->cols_dofs->begin();
194  vit != fe_ptr.get()->cols_dofs->end(); vit++) {
195  const int idx = TAG::get_index(vit);
196  if (idx >= 0)
197  dofs_col_view.push_back(idx);
198  if (idx >= p_miit->getNbDofsCol()) {
199  std::ostringstream zz;
200  zz << "rank " << rAnk << " ";
201  zz << *(*vit) << std::endl;
202  SETERRQ(cOmm, PETSC_ERR_ARG_SIZ, zz.str().c_str());
203  }
204  }
205  if (verb >= NOISY) {
206  std::stringstream ss;
207  ss << "rank " << rAnk << ": numeredDofsCols" << std::endl;
208  FENumeredDofEntity_multiIndex::iterator dit, hi_dit;
209  dit = fe_ptr.get()->cols_dofs->begin();
210  hi_dit = fe_ptr.get()->cols_dofs->end();
211  for (; dit != hi_dit; dit++) {
212  ss << "\t" << **dit << std::endl;
213  }
214  PetscSynchronizedPrintf(cOmm, "%s", ss.str().c_str());
215  }
216  } else
217  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
218  "Element should be here, otherwise matrix will have missing "
219  "elements");
220  }
221  }
222 
224 }
225 
226 template <typename TAG>
228  ProblemsByName::iterator p_miit, const MatType type,
229  std::vector<PetscInt> &i, std::vector<PetscInt> &j, const bool no_diagonals,
230  int verb) const {
232  PetscLogEventBegin(MOFEM_EVENT_createMat, 0, 0, 0, 0);
233 
234  typedef
235  typename boost::multi_index::index<NumeredDofEntity_multiIndex, TAG>::type
236  NumeredDofEntitysByIdx;
237 
238  // Get multi-indices for rows and columns
239  const NumeredDofEntitysByIdx &dofs_row_by_idx =
240  p_miit->numeredDofsRows->get<TAG>();
241  int nb_dofs_row = p_miit->getNbDofsRow();
242  if (nb_dofs_row == 0) {
243  SETERRQ1(cOmm, MOFEM_DATA_INCONSISTENCY, "problem <%s> has zero rows",
244  p_miit->getName().c_str());
245  }
246 
247  // Get adjacencies form other processors
248  std::map<int, std::vector<int>> adjacent_dofs_on_other_parts;
249 
250  // If not partitioned set petsc layout for matrix. If partitioned need to get
251  // adjacencies form other parts. Note if algebra is only partitioned no need
252  // to collect adjacencies form other entities. Those are already on mesh
253  // which is assumed that is on each processor the same.
254  typename boost::multi_index::index<NumeredDofEntity_multiIndex,
255  TAG>::type::iterator miit_row,
256  hi_miit_row;
257 
258  if (TAG::IamNotPartitioned) {
259 
260  // Get range of local indices
261  PetscLayout layout;
262  CHKERR PetscLayoutCreate(cOmm, &layout);
263  CHKERR PetscLayoutSetBlockSize(layout, 1);
264  CHKERR PetscLayoutSetSize(layout, nb_dofs_row);
265  CHKERR PetscLayoutSetUp(layout);
266  PetscInt rstart, rend;
267  CHKERR PetscLayoutGetRange(layout, &rstart, &rend);
268  CHKERR PetscLayoutDestroy(&layout);
269  if (verb >= VERBOSE) {
270  PetscSynchronizedPrintf(cOmm, "\tcreate_Mat: row lower %d row upper %d\n",
271  rstart, rend);
272  PetscSynchronizedFlush(cOmm, PETSC_STDOUT);
273  }
274  miit_row = dofs_row_by_idx.lower_bound(rstart);
275  hi_miit_row = dofs_row_by_idx.lower_bound(rend);
276  if (std::distance(miit_row, hi_miit_row) != rend - rstart) {
277  SETERRQ4(
278  cOmm, PETSC_ERR_ARG_SIZ,
279  "data inconsistency, std::distance(miit_row,hi_miit_row) != rend - "
280  "rstart (%d != %d - %d = %d) ",
281  std::distance(miit_row, hi_miit_row), rend, rstart, rend - rstart);
282  }
283 
284  } else {
285 
286  // Make sure it is a PETSc cOmm
287  CHKERR PetscCommDuplicate(cOmm, &cOmm, NULL);
288 
289  // get adjacent nodes on other partitions
290  std::vector<std::vector<int>> dofs_vec(sIze);
291 
292  boost::shared_ptr<FieldEntity> mofem_ent_ptr;
293  std::vector<int> dofs_col_view;
294 
295  typename boost::multi_index::index<NumeredDofEntity_multiIndex,
296  TAG>::type::iterator mit_row,
297  hi_mit_row;
298 
299  mit_row = dofs_row_by_idx.begin();
300  hi_mit_row = dofs_row_by_idx.end();
301  for (; mit_row != hi_mit_row; mit_row++) {
302 
303  // Shared or multishared row and not owned. Those data should be send to
304  // other side.
305 
306  // Get entity adjacencies, no need to repeat that operation for dofs when
307  // are on the same entity. For simplicity is assumed that those sheered
308  // the same adjacencies.
309  unsigned char pstatus = (*mit_row)->getPStatus();
310  if ((pstatus & PSTATUS_NOT_OWNED) &&
311  (pstatus & (PSTATUS_SHARED | PSTATUS_MULTISHARED))) {
312 
313  bool get_adj_col = true;
314  if (mofem_ent_ptr) {
315  if (mofem_ent_ptr->getGlobalUniqueId() ==
316  (*mit_row)->getFieldEntityPtr()->getGlobalUniqueId()) {
317  get_adj_col = false;
318  }
319  }
320 
321  if (get_adj_col) {
322  // Get entity adjacencies
323  mofem_ent_ptr = (*mit_row)->getFieldEntityPtr();
324  CHKERR getEntityAdjacenies<TAG>(p_miit, mit_row, mofem_ent_ptr,
325  dofs_col_view, verb);
326  // Sort, unique and resize dofs_col_view
327  {
328  std::sort(dofs_col_view.begin(), dofs_col_view.end());
329  std::vector<int>::iterator new_end =
330  std::unique(dofs_col_view.begin(), dofs_col_view.end());
331  int new_size = std::distance(dofs_col_view.begin(), new_end);
332  dofs_col_view.resize(new_size);
333  }
334  // Add that row. Patterns is that first index is row index, second is
335  // size of adjacencies after that follows column adjacencies.
336  int owner = (*mit_row)->getOwnerProc();
337  dofs_vec[owner].emplace_back(TAG::get_index(mit_row)); // row index
338  dofs_vec[owner].emplace_back(
339  dofs_col_view.size()); // nb. of column adjacencies
340  // add adjacent cools
341  dofs_vec[owner].insert(dofs_vec[owner].end(), dofs_col_view.begin(),
342  dofs_col_view.end());
343  }
344  }
345  }
346 
347  int nsends = 0; // number of messages to send
348  std::vector<int> dofs_vec_length(sIze); // length of the message to proc
349  for (int proc = 0; proc < sIze; proc++) {
350 
351  if (!dofs_vec[proc].empty()) {
352 
353  dofs_vec_length[proc] = dofs_vec[proc].size();
354  nsends++;
355 
356  } else {
357 
358  dofs_vec_length[proc] = 0;
359  }
360  }
361 
362  std::vector<MPI_Status> status(sIze);
363 
364  // Computes the number of messages a node expects to receive
365  int nrecvs; // number of messages received
366  CHKERR PetscGatherNumberOfMessages(cOmm, NULL, &dofs_vec_length[0],
367  &nrecvs);
368 
369  // Computes info about messages that a MPI-node will receive, including
370  // (from-id,length) pairs for each message.
371  int *onodes; // list of node-ids from which messages are expected
372  int *olengths; // corresponding message lengths
373  CHKERR PetscGatherMessageLengths(cOmm, nsends, nrecvs, &dofs_vec_length[0],
374  &onodes, &olengths);
375 
376  // Gets a unique new tag from a PETSc communicator.
377  int tag;
378  CHKERR PetscCommGetNewTag(cOmm, &tag);
379 
380  // Allocate a buffer sufficient to hold messages of size specified in
381  // olengths. And post Irecvs on these buffers using node info from onodes
382  int **rbuf; // must bee freed by user
383  MPI_Request *r_waits; // must bee freed by user
384 
385  // rbuf has a pointers to messages. It has size of of nrecvs (number of
386  // messages) +1. In the first index a block is allocated,
387  // such that rbuf[i] = rbuf[i-1]+olengths[i-1].
388  CHKERR PetscPostIrecvInt(cOmm, tag, nrecvs, onodes, olengths, &rbuf,
389  &r_waits);
390 
391  MPI_Request *s_waits; // status of sens messages
392  CHKERR PetscMalloc1(nsends, &s_waits);
393 
394  // Send messages
395  for (int proc = 0, kk = 0; proc < sIze; proc++) {
396  if (!dofs_vec_length[proc])
397  continue; // no message to send to this proc
398  CHKERR MPI_Isend(&(dofs_vec[proc])[0], // buffer to send
399  dofs_vec_length[proc], // message length
400  MPIU_INT, proc, // to proc
401  tag, cOmm, s_waits + kk);
402  kk++;
403  }
404 
405  // Wait for received
406  if (nrecvs) {
407  CHKERR MPI_Waitall(nrecvs, r_waits, &status[0]);
408  }
409  // Wait for send messages
410  if (nsends) {
411  CHKERR MPI_Waitall(nsends, s_waits, &status[0]);
412  }
413 
414  for (int kk = 0; kk < nrecvs; kk++) {
415 
416  int len = olengths[kk];
417  int *data_from_proc = rbuf[kk];
418 
419  for (int ii = 0; ii < len;) {
420 
421  int row_idx = data_from_proc[ii++]; // get row number
422  int nb_adj_dofs = data_from_proc[ii++]; // get nb. of adjacent dofs
423 
424  if (debug) {
425 
426  DofByGlobalPetscIndex::iterator dit;
427  dit = p_miit->numeredDofsRows->get<PetscGlobalIdx_mi_tag>().find(
428  row_idx);
429  if (dit ==
430  p_miit->numeredDofsRows->get<PetscGlobalIdx_mi_tag>().end()) {
431  SETERRQ1(cOmm, MOFEM_DATA_INCONSISTENCY,
432  "dof %d can not be found in problem", row_idx);
433  }
434  }
435 
436  for (int jj = 0; jj < nb_adj_dofs; jj++) {
437  adjacent_dofs_on_other_parts[row_idx].push_back(data_from_proc[ii++]);
438  }
439  }
440  }
441 
442  // Cleaning
443  CHKERR PetscFree(s_waits);
444  CHKERR PetscFree(rbuf[0]);
445  CHKERR PetscFree(rbuf);
446  CHKERR PetscFree(r_waits);
447  CHKERR PetscFree(onodes);
448  CHKERR PetscFree(olengths);
449 
450  miit_row = dofs_row_by_idx.begin();
451  hi_miit_row = dofs_row_by_idx.end();
452  }
453 
454  boost::shared_ptr<FieldEntity> mofem_ent_ptr;
455  int row_last_evaluated_idx = -1;
456 
457  std::vector<int> dofs_vec;
458  std::vector<int> dofs_col_view;
459 
460  // loop local rows
461  int nb_loc_row_from_iterators = 0;
462  unsigned int rows_to_fill = p_miit->getNbLocalDofsRow();
463  i.reserve(rows_to_fill + 1);
464  for (; miit_row != hi_miit_row; miit_row++) {
465 
466  if (!TAG::IamNotPartitioned) {
467  if (static_cast<int>((*miit_row)->getPart()) != rAnk)
468  continue;
469  }
470  // This is only for cross-check if everything is ok
471  nb_loc_row_from_iterators++;
472 
473  // add next row to compressed matrix
474  i.push_back(j.size());
475 
476  // Get entity adjacencies, no need to repeat that operation for dofs when
477  // are on the same entity. For simplicity is assumed that those share the
478  // same adjacencies.
479  if ((!mofem_ent_ptr)
480  ? 1
481  : (mofem_ent_ptr->getGlobalUniqueId() !=
482  (*miit_row)->getFieldEntityPtr()->getGlobalUniqueId())) {
483 
484  if (verb >= NOISY) {
485  std::stringstream ss;
486  ss << "rank " << rAnk << ": row " << **miit_row << std::endl;
487  PetscSynchronizedPrintf(cOmm, "%s", ss.str().c_str());
488  }
489 
490  // get entity adjacencies
491  mofem_ent_ptr = (*miit_row)->getFieldEntityPtr();
492  CHKERR getEntityAdjacenies<TAG>(p_miit, miit_row, mofem_ent_ptr,
493  dofs_col_view, verb);
494  row_last_evaluated_idx = TAG::get_index(miit_row);
495 
496  dofs_vec.resize(0);
497  // insert dofs_col_view
498  dofs_vec.insert(dofs_vec.end(), dofs_col_view.begin(),
499  dofs_col_view.end());
500 
501  unsigned char pstatus = (*miit_row)->getPStatus();
502  if (pstatus > 0) {
503  std::map<int, std::vector<int>>::iterator mit;
504  mit = adjacent_dofs_on_other_parts.find(row_last_evaluated_idx);
505  if (mit == adjacent_dofs_on_other_parts.end()) {
506  // NOTE: Dof can adjacent to other part but no elements are there
507  // which use that dof std::cerr << *miit_row << std::endl; SETERRQ1(
508  // cOmm,MOFEM_DATA_INCONSISTENCY,
509  // "data inconsistency row_last_evaluated_idx = %d",
510  // row_last_evaluated_idx
511  // );
512  } else {
513  dofs_vec.insert(dofs_vec.end(), mit->second.begin(),
514  mit->second.end());
515  }
516  }
517 
518  // sort and make unique
519  sort(dofs_vec.begin(), dofs_vec.end());
520  std::vector<int>::iterator new_end =
521  unique(dofs_vec.begin(), dofs_vec.end());
522  int new_size = std::distance(dofs_vec.begin(), new_end);
523  dofs_vec.resize(new_size);
524  if (verb >= NOISY) {
525  std::stringstream ss;
526  ss << "rank " << rAnk << ": dofs_vec for " << *mofem_ent_ptr
527  << std::endl;
528  PetscSynchronizedPrintf(cOmm, "%s", ss.str().c_str());
529  }
530  }
531 
532  // Try to be smart reserving memory
533  if (j.capacity() < j.size() + dofs_vec.size()) {
534 
535  unsigned int nb_nonzero = j.size() + dofs_vec.size();
536  unsigned int average_row_fill =
537  nb_nonzero / i.size() + nb_nonzero % i.size();
538  if (j.capacity() < rows_to_fill * average_row_fill) {
539  j.reserve(rows_to_fill * average_row_fill);
540  }
541  }
542 
543  // add indices to compressed matrix
544  if (verb >= VERY_VERBOSE) {
545  PetscSynchronizedPrintf(cOmm, "rank %d: ", rAnk);
546  }
547  std::vector<int>::iterator diit, hi_diit;
548  diit = dofs_vec.begin();
549  hi_diit = dofs_vec.end();
550  for (; diit != hi_diit; diit++) {
551 
552  if (no_diagonals) {
553  if (*diit == TAG::get_index(miit_row)) {
554  continue;
555  }
556  }
557  j.push_back(*diit);
558 
559  if (verb >= VERY_VERBOSE) {
560  PetscSynchronizedPrintf(cOmm, "%d ", *diit);
561  }
562  }
563  if (verb >= VERY_VERBOSE) {
564  PetscSynchronizedPrintf(cOmm, "\n", *diit);
565  }
566  }
567 
568  if (verb >= VERY_VERBOSE) {
569  PetscSynchronizedFlush(cOmm, PETSC_STDOUT);
570  }
571 
572  // build adj matrix
573  i.push_back(j.size());
574 
575  if (strcmp(type, MATMPIADJ) == 0) {
576 
577  // Adjacency matrix used to partition problems, f.e. METIS
578  if (i.size() - 1 != (unsigned int)nb_loc_row_from_iterators) {
579  SETERRQ(cOmm, PETSC_ERR_ARG_SIZ, "data inconsistency");
580  }
581 
582  } else if (strcmp(type, MATMPIAIJ) == 0) {
583 
584  // Compressed MPIADJ matrix
585  if (i.size() - 1 != (unsigned int)nb_loc_row_from_iterators) {
586  SETERRQ(cOmm, PETSC_ERR_ARG_SIZ, "data inconsistency");
587  }
588  PetscInt nb_local_dofs_row = p_miit->getNbLocalDofsRow();
589  if ((unsigned int)nb_local_dofs_row != i.size() - 1) {
590  SETERRQ(cOmm, PETSC_ERR_ARG_SIZ, "data inconsistency");
591  }
592 
593  } else if (strcmp(type, MATAIJ) == 0) {
594 
595  // Sequential compressed ADJ matrix
596  if (i.size() - 1 != (unsigned int)nb_loc_row_from_iterators) {
597  SETERRQ(cOmm, PETSC_ERR_ARG_SIZ, "data inconsistency");
598  }
599  PetscInt nb_local_dofs_row = p_miit->getNbLocalDofsRow();
600  if ((unsigned int)nb_local_dofs_row != i.size() - 1) {
601  SETERRQ(cOmm, PETSC_ERR_ARG_SIZ, "data inconsistency");
602  }
603 
604  } else {
605 
606  SETERRQ(cOmm, PETSC_ERR_ARG_NULL, "not implemented");
607  }
608 
609  PetscLogEventEnd(MOFEM_EVENT_createMat, 0, 0, 0, 0);
611 }
612 
614  UnknownInterface **iface) const {
616  *iface = NULL;
617  if (uuid == IDD_MOFEMMatrixManager) {
618  *iface = const_cast<MatrixManager *>(this);
620  }
621  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
623 }
624 
626  : cOre(const_cast<MoFEM::Core &>(core)) {
627  PetscLogEventRegister("MatrixManagerCreateMPIAIJWithArrays", 0,
629  PetscLogEventRegister("MatrixManagerCreateMPIAdjWithArrays", 0,
631  PetscLogEventRegister("MatrixManagerCreateSeqAIJWithArrays", 0,
633  PetscLogEventRegister("MatrixManagerCheckMPIAIJWithArraysMatrixFillIn", 0,
635 }
637 
638 template <>
639 MoFEMErrorCode MatrixManager::createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(
640  const std::string name, Mat *Aij, int verb) {
641  MoFEM::CoreInterface &m_field = cOre;
642  CreateRowComressedADJMatrix *core_ptr =
643  static_cast<CreateRowComressedADJMatrix *>(&cOre);
645  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
646 
647  const Problem_multiIndex *problems_ptr;
648  CHKERR m_field.get_problems(&problems_ptr);
649  auto &prb = problems_ptr->get<Problem_mi_tag>();
650  auto p_miit = prb.find(name);
651  if (p_miit == prb.end()) {
652  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
653  "problem < %s > is not found (top tip: check spelling)",
654  name.c_str());
655  }
656 
657  std::vector<int> i_vec, j_vec;
659  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
660 
661  int nb_row_dofs = p_miit->getNbDofsRow();
662  int nb_col_dofs = p_miit->getNbDofsCol();
663  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
664  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
665 
666  CHKERR ::MatCreateMPIAIJWithArrays(
667  m_field.get_comm(), nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
668  nb_col_dofs, &*i_vec.begin(), &*j_vec.begin(), PETSC_NULL, Aij);
669 
670  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
672 }
673 
674 template <>
676 MatrixManager::createMPIAdjWithArrays<Idx_mi_tag>(const std::string name,
677  Mat *Adj, int verb) {
678  MoFEM::CoreInterface &m_field = cOre;
679  CreateRowComressedADJMatrix *core_ptr =
680  static_cast<CreateRowComressedADJMatrix *>(&cOre);
682  PetscLogEventBegin(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
683 
684  const Problem_multiIndex *problems_ptr;
685  CHKERR m_field.get_problems(&problems_ptr);
686  auto &prb = problems_ptr->get<Problem_mi_tag>();
687  auto p_miit = prb.find(name);
688  if (p_miit == prb.end()) {
689  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
690  "problem < %s > is not found (top tip: check spelling)",
691  name.c_str());
692  }
693 
694  std::vector<int> i_vec, j_vec;
695  CHKERR core_ptr->createMatArrays<Idx_mi_tag>(p_miit, MATMPIADJ, i_vec, j_vec,
696  true, verb);
697  int *_i, *_j;
698  CHKERR PetscMalloc(i_vec.size() * sizeof(int), &_i);
699  CHKERR PetscMalloc(j_vec.size() * sizeof(int), &_j);
700  copy(i_vec.begin(), i_vec.end(), _i);
701  copy(j_vec.begin(), j_vec.end(), _j);
702 
703  int nb_col_dofs = p_miit->getNbDofsCol();
704  CHKERR MatCreateMPIAdj(m_field.get_comm(), i_vec.size() - 1, nb_col_dofs, _i,
705  _j, PETSC_NULL, Adj);
706  CHKERR MatSetOption(*Adj, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
707 
708  PetscLogEventEnd(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
710 }
711 
712 template <>
713 MoFEMErrorCode MatrixManager::createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(
714  const std::string name, Mat *Aij, int verb) {
715  MoFEM::CoreInterface &m_field = cOre;
716  CreateRowComressedADJMatrix *core_ptr =
717  static_cast<CreateRowComressedADJMatrix *>(&cOre);
719  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
720 
721  const Problem_multiIndex *problems_ptr;
722  CHKERR m_field.get_problems(&problems_ptr);
723  auto &prb = problems_ptr->get<Problem_mi_tag>();
724  auto p_miit = prb.find(name);
725  if (p_miit == prb.end()) {
726  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
727  "problem < %s > is not found (top tip: check spelling)",
728  name.c_str());
729  }
730 
731  std::vector<int> i_vec, j_vec;
732  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(p_miit, MATAIJ, i_vec,
733  j_vec, false, verb);
734 
735  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
736  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
737 
738  double *_a;
739  CHKERR PetscMalloc(j_vec.size() * sizeof(double), &_a);
740 
741  Mat tmpMat;
742  CHKERR ::MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, nb_local_dofs_row,
743  nb_local_dofs_col, &*i_vec.begin(),
744  &*j_vec.begin(), _a, &tmpMat);
745  CHKERR MatDuplicate(tmpMat, MAT_SHARE_NONZERO_PATTERN, Aij);
746  CHKERR MatDestroy(&tmpMat);
747 
748  CHKERR PetscFree(_a);
749 
750  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
752 }
753 
754 template <>
756 MatrixManager::checkMPIAIJWithArraysMatrixFillIn<PetscGlobalIdx_mi_tag>(
757  const std::string problem_name, int row_print, int col_print, int verb) {
758  MoFEM::CoreInterface &m_field = cOre;
760  PetscLogEventBegin(MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn, 0, 0, 0, 0);
761 
762  struct TestMatrixFillIn : public FEMethod {
763  CoreInterface *mFieldPtr;
764 
765  Mat A;
766 
767  int rowPrint, colPrint;
768 
769  TestMatrixFillIn(CoreInterface *m_field_ptr, Mat a, int row_print,
770  int col_print)
771  : mFieldPtr(m_field_ptr), A(a), rowPrint(row_print),
772  colPrint(col_print){};
773 
774  MoFEMErrorCode preProcess() {
777  }
778 
779  MoFEMErrorCode operator()() {
781 
782  if (refinedFiniteElementsPtr->find(
783  numeredEntFiniteElementPtr->getEnt()) ==
784  refinedFiniteElementsPtr->end()) {
785  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
786  "data inconsistency");
787  }
788 
789  for (FENumeredDofEntity_multiIndex::iterator cit = colPtr->begin();
790  cit != colPtr->end(); cit++) {
791 
792  if (refinedEntitiesPtr->find((*cit)->getEnt()) ==
793  refinedEntitiesPtr->end()) {
794  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
795  "data inconsistency");
796  }
797  if (!(*cit)->getActive()) {
798  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
799  "data inconsistency");
800  }
801 
802  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
803  Composite_Unique_mi_tag>::type::iterator ait;
804  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
805  boost::make_tuple((*cit)->getFieldEntityPtr()->getGlobalUniqueId(),
806  numeredEntFiniteElementPtr->getGlobalUniqueId()));
807  if (ait == adjacenciesPtr->end()) {
808  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
809  "adjacencies data inconsistency");
810  } else {
811  UId uid = ait->getEntUniqueId();
812  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
813  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
814  "data inconsistency");
815  }
816  if (dofsPtr->find((*cit)->getGlobalUniqueId()) == dofsPtr->end()) {
817  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
818  "data inconsistency");
819  }
820  }
821 
822  if ((*cit)->getEntType() != MBVERTEX) {
823 
824  FENumeredDofEntity_multiIndex::index<
825  Composite_Name_Type_And_Side_Number_mi_tag>::type::iterator dit,
826  hi_dit;
827  dit = colPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
828  .lower_bound(boost::make_tuple(
829  (*cit)->getName(), (*cit)->getEntType(),
830  (*cit)->sideNumberPtr->side_number));
831  hi_dit = colPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
832  .upper_bound(boost::make_tuple(
833  (*cit)->getName(), (*cit)->getEntType(),
834  (*cit)->sideNumberPtr->side_number));
835  int nb_dofs_on_ent = std::distance(dit, hi_dit);
836 
837  int max_order = (*cit)->getMaxOrder();
838  if ((*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order) !=
839  nb_dofs_on_ent) {
840  std::cerr << "Warning: Number of Dofs in Col diffrent than number "
841  "of dofs for given entity order "
842  << (*cit)->getNbOfCoeffs() *
843  (*cit)->getOrderNbDofs(max_order)
844  << " " << nb_dofs_on_ent << std::endl;
845  }
846  }
847  }
848 
849  FENumeredDofEntity_multiIndex::iterator rit = rowPtr->begin();
850  for (; rit != rowPtr->end(); rit++) {
851 
852  if (refinedEntitiesPtr->find((*rit)->getEnt()) ==
853  refinedEntitiesPtr->end()) {
854  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
855  "data inconsistency");
856  }
857  if (!(*rit)->getActive()) {
858  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
859  "data inconsistency");
860  }
861 
862  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
863  Composite_Unique_mi_tag>::type::iterator ait;
864  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
865  boost::make_tuple((*rit)->getFieldEntityPtr()->getGlobalUniqueId(),
866  numeredEntFiniteElementPtr->getGlobalUniqueId()));
867  if (ait == adjacenciesPtr->end()) {
868  std::ostringstream ss;
869  ss << *(*rit) << std::endl;
870  ss << *numeredEntFiniteElementPtr << std::endl;
871  ss << "dof: " << (*rit)->getBitRefLevel() << std::endl;
872  ss << "fe: " << numeredEntFiniteElementPtr->getBitRefLevel()
873  << std::endl;
874  ss << "problem: " << problemPtr->getBitRefLevel() << std::endl;
875  ss << "problem mask: " << problemPtr->getMaskBitRefLevel()
876  << std::endl;
877  PetscPrintf(mFieldPtr->get_comm(), "%s", ss.str().c_str());
878  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
879  "adjacencies data inconsistency");
880  } else {
881  UId uid = ait->getEntUniqueId();
882  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
883  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
884  "data inconsistency");
885  }
886  if (dofsPtr->find((*rit)->getGlobalUniqueId()) == dofsPtr->end()) {
887  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
888  "data inconsistency");
889  }
890  }
891  int row = (*rit)->getPetscGlobalDofIdx();
892 
893  FENumeredDofEntity_multiIndex::iterator cit = colPtr->begin();
894  for (; cit != colPtr->end(); cit++) {
895 
896  int col = (*cit)->getPetscGlobalDofIdx();
897 
898  if (row == rowPrint && col == colPrint) {
899 
900  std::ostringstream ss;
901  ss << "fe:\n" << *numeredEntFiniteElementPtr << std::endl;
902  ss << "row:\n" << *(*rit) << std::endl;
903  ss << "col:\n" << *(*cit) << std::endl;
904 
905  ss << "fe:\n"
906  << numeredEntFiniteElementPtr->getBitRefLevel() << std::endl;
907  ss << "row:\n" << (*rit)->getBitRefLevel() << std::endl;
908  ss << "col:\n" << (*cit)->getBitRefLevel() << std::endl;
909 
910  std::cerr << ss.str() << std::endl;
911 
912  // PetscPrintf(mFieldPtr->get_comm(),"%s\n",ss.str().c_str());
913  }
914 
915  CHKERR MatSetValue(A, row, col, 1, INSERT_VALUES);
916  }
917 
918  if ((*rit)->getEntType() != MBVERTEX) {
919 
920  FENumeredDofEntity_multiIndex::index<
921  Composite_Name_Type_And_Side_Number_mi_tag>::type::iterator dit,
922  hi_dit;
923  dit = rowPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
924  .lower_bound(boost::make_tuple(
925  (*rit)->getName(), (*rit)->getEntType(),
926  (*rit)->sideNumberPtr->side_number));
927  hi_dit = rowPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
928  .upper_bound(boost::make_tuple(
929  (*rit)->getName(), (*rit)->getEntType(),
930  (*rit)->sideNumberPtr->side_number));
931  int nb_dofs_on_ent = std::distance(dit, hi_dit);
932 
933  int max_order = (*rit)->getMaxOrder();
934  if ((*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order) !=
935  nb_dofs_on_ent) {
936  std::cerr << "Warning: Number of Dofs in Row diffrent than number "
937  "of dofs for given entity order "
938  << (*rit)->getNbOfCoeffs() *
939  (*rit)->getOrderNbDofs(max_order)
940  << " " << nb_dofs_on_ent << std::endl;
941  }
942  }
943  }
944 
946  }
947 
948  MoFEMErrorCode postProcess() {
950 
951  // cerr << mFieldPtr->get_comm_rank() << endl;
952  CHKERR MatAssemblyBegin(A, MAT_FLUSH_ASSEMBLY);
953  CHKERR MatAssemblyEnd(A, MAT_FLUSH_ASSEMBLY);
954 
956  }
957  };
958 
959  // create matrix
960  Mat A;
961  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, &A, verb);
962  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
963 
964  if (verb >= VERY_VERBOSE) {
965  MatView(A, PETSC_VIEWER_STDOUT_WORLD);
966  }
967 
968  if (verb >= NOISY) {
969  MatView(A, PETSC_VIEWER_DRAW_WORLD);
970  std::string wait;
971  std::cin >> wait;
972  }
973 
974  TestMatrixFillIn method(&m_field, A, row_print, col_print);
975 
976  // get problem
977  const Problem_multiIndex *problems_ptr;
978  CHKERR m_field.get_problems(&problems_ptr);
979  auto &prb_set = problems_ptr->get<Problem_mi_tag>();
980  auto p_miit = prb_set.find(problem_name);
981  if (p_miit == prb_set.end())
982  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
983  "problem < %s > not found (top tip: check spelling)",
984  problem_name.c_str());
985 
986  if (verb >= VERBOSE)
987  PetscPrintf(m_field.get_comm(), "check problem < %s >\n",
988  problem_name.c_str());
989 
990  // loop all elements in problem and check if assemble is without error
991  const FiniteElement_multiIndex *fe_ptr;
992  CHKERR m_field.get_finite_elements(&fe_ptr);
993  for (auto &fe : *fe_ptr) {
994  if (verb >= VERBOSE)
995  PetscPrintf(m_field.get_comm(), "\tcheck element %s\n",
996  fe->getName().c_str());
997  CHKERR m_field.loop_finite_elements(problem_name, fe->getName(), method,
998  nullptr, MF_EXIST, verb);
999  }
1000 
1001  CHKERR MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
1002  CHKERR MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
1003  CHKERR MatDestroy(&A);
1004 
1005  PetscLogEventEnd(MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn, 0, 0, 0, 0);
1006 
1008 }
1009 
1010 } // namespace MoFEM
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
MoFEM interface unique ID.
int sIze
MoFEM communicator size.
Definition: Core.hpp:849
MPI_Comm cOmm
MoFEM communicator.
Definition: Core.hpp:846
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscLogEvent MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
PetscLogEvent MOFEM_EVENT_createMat
Definition: Core.hpp:838
virtual MoFEMErrorCode loop_finite_elements(const std::string &problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
Core (interface) class.
Definition: Core.hpp:50
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
std::reference_wrapper< moab::Interface > moab
moab database
Definition: Core.hpp:265
Create compressed matrix.
FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index< Unique_mi_tag >::type AdjByEnt
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< const_mem_fun< NumeredDofEntity::interface_type_DofEntity, const UId &, &NumeredDofEntity::getGlobalUniqueId > > > > NumeredDofEntity_multiIndex_uid_view_ordered
virtual MoFEMErrorCode get_finite_elements(const FiniteElement_multiIndex **fe_ptr) const =0
Get finite elements multi-index.
int rAnk
MOFEM communicator rank.
Definition: Core.hpp:850
~MatrixManager()
Destructor.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
static const bool debug
FieldEntityEntFiniteElementAdjacencyMap_multiIndex entFEAdjacencies
adjacencies of elements to dofs
Definition: Core.hpp:255
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
virtual MoFEMErrorCode get_problems(const Problem_multiIndex **problems_ptr) const =0
Get pointer to problems multi-index.
CreateRowComressedADJMatrix(moab::Interface &moab, MPI_Comm comm=PETSC_COMM_WORLD, int verbose=1)
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
MatrixManager(const MoFEM::Core &core)
#define CHKERR
Inline error check.
Definition: definitions.h:596
static const MOFEMuuid IDD_MOFEMMatrixManager
MoFEMErrorCode buildFECol(ProblemsByName::iterator p_miit, boost::shared_ptr< EntFiniteElement > ent_fe_ptr, bool do_cols_prob, boost::shared_ptr< NumeredEntFiniteElement > &fe_ptr) const
NumeredDofEntity_multiIndex::index< PetscGlobalIdx_mi_tag >::type DofByGlobalPetscIndex
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Problem_multiIndex::index< Problem_mi_tag >::type ProblemsByName
int verbose
Verbosity level.
Definition: Core.hpp:885
MoFEMErrorCode createMatArrays(ProblemsByName::iterator p_miit, const MatType type, std::vector< PetscInt > &i, std::vector< PetscInt > &j, const bool no_diagonals=true, int verb=QUIET) const
Create matrix adjacencies.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
virtual MPI_Comm & get_comm() const =0
MoFEMErrorCode getEntityAdjacenies(ProblemsByName::iterator p_miit, typename boost::multi_index::index< NumeredDofEntity_multiIndex, TAG >::type::iterator mit_row, boost::shared_ptr< FieldEntity > mofem_ent_ptr, std::vector< int > &dofs_col_view, int verb) const
Get element adjacencies.
uint128_t UId
Unique Id.
Definition: Types.hpp:41
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.