Reference documentation for deal.II version GIT 6da2e5d553 2022-07-01 18:55:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
connectivity.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2020 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_tria_connectivity_h
17 #define dealii_tria_connectivity_h
18 
19 #include <deal.II/base/config.h>
20 
22 #include <deal.II/base/ndarray.h>
23 
26 
27 
29 
30 
31 namespace internal
32 {
33  namespace TriangulationImplementation
34  {
39  struct CellTypeBase
40  {
44  virtual ~CellTypeBase() = default;
45 
49  virtual unsigned int
50  n_entities(const unsigned int d) const
51  {
52  Assert(false, ExcNotImplemented());
53  (void)d;
54 
55  return 0;
56  }
57 
61  virtual ::ArrayView<const unsigned int>
62  vertices_of_entity(const unsigned int d, const unsigned int e) const
63  {
64  Assert(false, ExcNotImplemented());
65  (void)d;
66  (void)e;
67 
68  return {};
69  }
70 
74  virtual ::ReferenceCell
75  type_of_entity(const unsigned int d, const unsigned int e) const
76  {
77  Assert(false, ExcNotImplemented());
78  (void)d;
79  (void)e;
80 
82  }
83 
87  virtual unsigned int
88  n_lines_of_surface(const unsigned int face) const
89  {
90  Assert(false, ExcNotImplemented());
91  (void)face;
92 
93  return 0;
94  }
95 
99  virtual unsigned int
100  nth_line_of_surface(const unsigned int line,
101  const unsigned int face) const
102  {
103  Assert(false, ExcNotImplemented());
104  (void)line;
105  (void)face;
106 
107  return 0;
108  }
109 
113  virtual const std::array<unsigned int, 2> &
114  vertices_of_nth_line_of_surface(const unsigned int line,
115  const unsigned int face) const
116  {
117  Assert(false, ExcNotImplemented());
118  (void)line;
119  (void)face;
120 
121  const static std::array<unsigned int, 2> table = {};
122 
123  return table;
124  }
125  };
126 
127 
128 
132  struct CellTypeLine : public CellTypeBase
133  {
135  vertices_of_entity(const unsigned int d,
136  const unsigned int e) const override
137  {
138  (void)e;
139 
140  if (d == 1)
141  {
142  static const std::array<unsigned int, 2> table = {{0, 1}};
143 
144  AssertDimension(e, 0);
145 
146  return {table};
147  }
148 
149  Assert(false, ExcNotImplemented());
150 
151  return {};
152  }
153 
155  type_of_entity(const unsigned int d, const unsigned int e) const override
156  {
157  (void)e;
158 
159  if (d == 1)
161 
162  Assert(false, ExcNotImplemented());
163 
165  }
166 
167  unsigned int
168  n_entities(const unsigned int d) const override
169  {
170  static std::array<unsigned int, 2> table = {{2, 1}};
171  return table[d];
172  }
173  };
174 
175 
176 
180  struct CellTypeTri : public CellTypeBase
181  {
183  vertices_of_entity(const unsigned int d,
184  const unsigned int e) const override
185  {
186  if (d == 2)
187  {
188  static const std::array<unsigned int, 3> table = {{0, 1, 2}};
189 
190  AssertDimension(e, 0);
191 
192  return {table};
193  }
194 
195  if (d == 1)
196  {
197  static const ::ndarray<unsigned int, 3, 2> table = {
198  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
199 
200  return {table[e]};
201  }
202 
203  Assert(false, ExcNotImplemented());
204 
205  return {};
206  }
207 
208  virtual ::ReferenceCell
209  type_of_entity(const unsigned int d, const unsigned int e) const override
210  {
211  (void)e;
212 
213  if (d == 2)
215 
216  if (d == 1)
218 
219  Assert(false, ExcNotImplemented());
220 
222  }
223 
224  unsigned int
225  n_entities(const unsigned int d) const override
226  {
227  static std::array<unsigned int, 3> table = {{3, 3, 1}};
228  return table[d];
229  }
230  };
231 
232 
233 
237  struct CellTypeQuad : public CellTypeBase
238  {
240  vertices_of_entity(const unsigned int d,
241  const unsigned int e) const override
242  {
243  if (d == 2)
244  {
245  static const std::array<unsigned int, 4> table = {{0, 1, 2, 3}};
246 
247  AssertDimension(e, 0);
248 
249  return {table};
250  }
251 
252  if (d == 1)
253  {
254  static const ::ndarray<unsigned int, 4, 2> table = {
255  {{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}};
256 
257  return {table[e]};
258  }
259 
260  Assert(false, ExcNotImplemented());
261 
262  return {};
263  }
264 
265  virtual ::ReferenceCell
266  type_of_entity(const unsigned int d, const unsigned int e) const override
267  {
268  (void)e;
269 
270  if (d == 2)
272 
273  if (d == 1)
275 
276  Assert(false, ExcNotImplemented());
277 
279  }
280 
281  unsigned int
282  n_entities(const unsigned int d) const override
283  {
284  static std::array<unsigned int, 3> table = {{4, 4, 1}};
285  return table[d];
286  }
287  };
288 
289 
290 
294  struct CellTypeTet : public CellTypeBase
295  {
297  vertices_of_entity(const unsigned int d,
298  const unsigned int e) const override
299  {
300  if (d == 3)
301  {
302  static const std::array<unsigned int, 4> table = {{0, 1, 2, 3}};
303 
304  AssertDimension(e, 0);
305 
306  return {table};
307  }
308 
309  if (d == 2)
310  {
311  static const ::ndarray<unsigned int, 4, 3> table = {
312  {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
313 
314  return {table[e]};
315  }
316 
317  if (d == 1)
318  {
319  static const ::ndarray<unsigned int, 6, 2> table = {
320  {{{0, 1}}, {{1, 2}}, {{2, 0}}, {{0, 3}}, {{1, 3}}, {{2, 3}}}};
321 
322  return {table[e]};
323  }
324 
325  Assert(false, ExcNotImplemented());
326 
327  return {};
328  }
329 
330  virtual ::ReferenceCell
331  type_of_entity(const unsigned int d, const unsigned int e) const override
332  {
333  (void)e;
334 
335  if (d == 3)
337 
338  if (d == 2)
340 
341  if (d == 1)
343 
344  Assert(false, ExcNotImplemented());
345 
347  }
348 
349  unsigned int
350  n_entities(const unsigned int d) const override
351  {
352  static std::array<unsigned int, 4> table = {{4, 6, 4, 1}};
353  return table[d];
354  }
355 
356  unsigned int
357  n_lines_of_surface(const unsigned int line) const override
358  {
359  (void)line;
360  return 3;
361  }
362 
363  unsigned int
364  nth_line_of_surface(const unsigned int line,
365  const unsigned int face) const override
366  {
367  const static ::ndarray<unsigned int, 4, 3> table = {
368  {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
369 
370  return table[face][line];
371  }
372 
373  const std::array<unsigned int, 2> &
374  vertices_of_nth_line_of_surface(const unsigned int line,
375  const unsigned int face) const override
376  {
377  const static ::ndarray<unsigned int, 4, 3, 2> table = {
378  {{{{{0, 1}}, {{1, 2}}, {{2, 0}}}},
379  {{{{1, 0}}, {{0, 3}}, {{3, 1}}}},
380  {{{{0, 2}}, {{2, 3}}, {{3, 0}}}},
381  {{{{2, 1}}, {{1, 3}}, {{3, 2}}}}}};
382 
383  return table[face][line];
384  }
385  };
386 
387 
393  {
395  vertices_of_entity(const unsigned int d,
396  const unsigned int e) const override
397  {
398  if (d == 3)
399  {
400  static const std::array<unsigned int, 5> table = {{0, 1, 2, 3, 4}};
401 
402  AssertDimension(e, 0);
403 
404  return {table};
405  }
406 
407  if (d == 2)
408  {
409  if (e == 0)
410  {
411  static const std::array<unsigned int, 4> table = {{0, 1, 2, 3}};
412  return {table};
413  }
414 
415  static const ::ndarray<unsigned int, 4, 3> table = {
416  {{{0, 2, 4}}, {{3, 1, 4}}, {{1, 0, 4}}, {{2, 3, 4}}}};
417 
418  return {table[e - 1]};
419  }
420 
421  if (d == 1)
422  {
423  static const ::ndarray<unsigned int, 8, 2> table = {
424  {{{0, 2}},
425  {{1, 3}},
426  {{0, 1}},
427  {{2, 3}},
428  {{0, 4}},
429  {{1, 4}},
430  {{2, 4}},
431  {{3, 4}}}};
432 
433  return {table[e]};
434  }
435 
436  Assert(false, ExcNotImplemented());
437 
438  return {};
439  }
440 
441  virtual ::ReferenceCell
442  type_of_entity(const unsigned int d, const unsigned int e) const override
443  {
444  (void)e;
445 
446  if (d == 3)
448 
449  if (d == 2 && e == 0)
451  else if (d == 2)
453 
454  if (d == 1)
456 
457  Assert(false, ExcNotImplemented());
458 
460  }
461 
462  unsigned int
463  n_entities(const unsigned int d) const override
464  {
465  static std::array<unsigned int, 4> table = {{5, 8, 5, 1}};
466  return table[d];
467  }
468 
469  unsigned int
470  n_lines_of_surface(const unsigned int surface) const override
471  {
472  if (surface == 0)
473  return 4;
474 
475  return 3;
476  }
477 
478  unsigned int
479  nth_line_of_surface(const unsigned int line,
480  const unsigned int face) const override
481  {
482  const static ::ndarray<unsigned int, 5, 4> table = {
483  {{{0, 1, 2, 3}},
484  {{0, 6, 4, numbers::invalid_unsigned_int}},
485  {{1, 5, 7, numbers::invalid_unsigned_int}},
486  {{2, 4, 5, numbers::invalid_unsigned_int}},
487  {{3, 7, 6, numbers::invalid_unsigned_int}}}};
488 
489  return table[face][line];
490  }
491 
492  const std::array<unsigned int, 2> &
493  vertices_of_nth_line_of_surface(const unsigned int line,
494  const unsigned int face) const override
495  {
496  static const unsigned int X = static_cast<unsigned int>(-1);
497 
498  const static ::ndarray<unsigned int, 5, 4, 2> table = {
499  {{{{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}},
500  {{{{0, 2}}, {{2, 4}}, {{4, 0}}, {{X, X}}}},
501  {{{{3, 1}}, {{1, 4}}, {{4, 3}}, {{X, X}}}},
502  {{{{1, 0}}, {{0, 4}}, {{4, 1}}, {{X, X}}}},
503  {{{{2, 3}}, {{3, 4}}, {{4, 2}}, {{X, X}}}}}};
504 
505  return table[face][line];
506  }
507  };
508 
509 
510 
514  struct CellTypeWedge : public CellTypeBase
515  {
517  vertices_of_entity(const unsigned int d,
518  const unsigned int e) const override
519  {
520  if (d == 3)
521  {
522  static const std::array<unsigned int, 6> table = {
523  {0, 1, 2, 3, 4, 5}};
524 
525  AssertDimension(e, 0);
526 
527  return {table};
528  }
529 
530  if (d == 2)
531  {
532  if (e == 0 || e == 1)
533  {
534  static const ::ndarray<unsigned int, 2, 3> table = {
535  {{{1, 0, 2}}, {{3, 4, 5}}}};
536 
537  return {table[e]};
538  }
539 
540  static const ::ndarray<unsigned int, 3, 4> table = {
541  {{{0, 1, 3, 4}}, {{1, 2, 4, 5}}, {{2, 0, 5, 3}}}};
542 
543  return {table[e - 2]};
544  }
545 
546  if (d == 1)
547  {
548  static const ::ndarray<unsigned int, 9, 2> table = {
549  {{{0, 1}},
550  {{1, 2}},
551  {{2, 0}},
552  {{3, 4}},
553  {{4, 5}},
554  {{5, 3}},
555  {{0, 3}},
556  {{1, 4}},
557  {{2, 5}}}};
558 
559  return {table[e]};
560  }
561 
562  Assert(false, ExcNotImplemented());
563 
564  return {};
565  }
566 
567  virtual ::ReferenceCell
568  type_of_entity(const unsigned int d, const unsigned int e) const override
569  {
570  (void)e;
571 
572  if (d == 3)
574 
575  if (d == 2 && e > 1)
577  else if (d == 2)
579 
580  if (d == 1)
582 
583  Assert(false, ExcNotImplemented());
584 
586  }
587 
588  unsigned int
589  n_entities(const unsigned int d) const override
590  {
591  static std::array<unsigned int, 4> table = {{6, 9, 5, 1}};
592  return table[d];
593  }
594 
595  unsigned int
596  n_lines_of_surface(const unsigned int surface) const override
597  {
598  if (surface > 1)
599  return 4;
600 
601  return 3;
602  }
603 
604  unsigned int
605  nth_line_of_surface(const unsigned int line,
606  const unsigned int face) const override
607  {
608  static const unsigned int X = static_cast<unsigned int>(-1);
609 
610  const static ::ndarray<unsigned int, 5, 4> table = {
611  {{{0, 2, 1, X}},
612  {{3, 4, 5, X}},
613  {{6, 7, 0, 3}},
614  {{7, 8, 1, 4}},
615  {{8, 6, 5, 2}}}};
616 
617  return table[face][line];
618  }
619 
620  const std::array<unsigned int, 2> &
621  vertices_of_nth_line_of_surface(const unsigned int line,
622  const unsigned int face) const override
623  {
624  static const unsigned int X = static_cast<unsigned int>(-1);
625 
626  const static ::ndarray<unsigned int, 5, 4, 2> table = {
627  {{{{{1, 0}}, {{0, 2}}, {{2, 1}}, {{X, X}}}},
628  {{{{3, 4}}, {{4, 5}}, {{5, 3}}, {{X, X}}}},
629  {{{{0, 3}}, {{1, 4}}, {{0, 1}}, {{3, 4}}}},
630  {{{{1, 4}}, {{2, 5}}, {{1, 2}}, {{4, 5}}}},
631  {{{{2, 5}}, {{0, 3}}, {{2, 0}}, {{5, 3}}}}}};
632 
633  return table[face][line];
634  }
635  };
636 
637 
638 
642  struct CellTypeHex : public CellTypeBase
643  {
645  vertices_of_entity(const unsigned int d,
646  const unsigned int e) const override
647  {
648  if (d == 3)
649  {
650  static const std::array<unsigned int, 8> table = {
651  {0, 1, 2, 3, 4, 5, 6, 7}};
652 
653  AssertDimension(e, 0);
654 
655  return {table};
656  }
657 
658  if (d == 2)
659  {
660  static const ::ndarray<unsigned int, 6, 4> table = {
661  {{{0, 2, 4, 6}},
662  {{1, 3, 5, 7}},
663  {{0, 4, 1, 5}},
664  {{2, 6, 3, 7}},
665  {{0, 1, 2, 3}},
666  {{4, 5, 6, 7}}}};
667 
668  return {table[e]};
669  }
670 
671  if (d == 1)
672  {
673  static const ::ndarray<unsigned int, 12, 2> table = {
674  {{{0, 2}},
675  {{1, 3}},
676  {{0, 1}},
677  {{2, 3}},
678  {{4, 6}},
679  {{5, 7}},
680  {{4, 5}},
681  {{6, 7}},
682  {{0, 4}},
683  {{1, 5}},
684  {{2, 6}},
685  {{3, 7}}}};
686 
687  return {table[e]};
688  }
689 
690  Assert(false, ExcNotImplemented());
691 
692  return {};
693  }
694 
695  virtual ::ReferenceCell
696  type_of_entity(const unsigned int d, const unsigned int e) const override
697  {
698  (void)e;
699 
700  if (d == 3)
702 
703  if (d == 2)
705 
706  if (d == 1)
708 
709  Assert(false, ExcNotImplemented());
710 
712  }
713 
714  unsigned int
715  n_entities(const unsigned int d) const override
716  {
717  static std::array<unsigned int, 4> table = {{8, 12, 6, 1}};
718  return table[d];
719  }
720 
721  unsigned int
722  n_lines_of_surface(const unsigned int surface) const override
723  {
724  (void)surface;
725  return 4;
726  }
727 
728  unsigned int
729  nth_line_of_surface(const unsigned int line,
730  const unsigned int face) const override
731  {
732  const static ::ndarray<unsigned int, 6, 4> table = {
733  {{{8, 10, 0, 4}},
734  {{9, 11, 1, 5}},
735  {{2, 6, 8, 9}},
736  {{3, 7, 10, 11}},
737  {{0, 1, 2, 3}},
738  {{4, 5, 6, 7}}}};
739 
740  return table[face][line];
741  }
742 
743  const std::array<unsigned int, 2> &
744  vertices_of_nth_line_of_surface(const unsigned int line,
745  const unsigned int face) const override
746  {
747  const static ::ndarray<unsigned int, 6, 4, 2> table = {
748  {{{{{0, 4}}, {{2, 6}}, {{0, 2}}, {{4, 6}}}},
749  {{{{1, 5}}, {{3, 7}}, {{1, 3}}, {{5, 7}}}},
750  {{{{0, 1}}, {{4, 5}}, {{0, 4}}, {{1, 5}}}},
751  {{{{2, 3}}, {{6, 7}}, {{2, 6}}, {{3, 7}}}},
752  {{{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}},
753  {{{{4, 6}}, {{5, 7}}, {{4, 5}}, {{6, 7}}}}}};
754 
755  return table[face][line];
756  }
757  };
758 
759 
760 
767  template <typename T = unsigned int>
768  struct CRS
769  {
773  CRS()
774  : ptr{0} {};
775 
779  CRS(const std::vector<std::size_t> &ptr, const std::vector<T> &col)
780  : ptr(ptr)
781  , col(col)
782  {}
783 
784  // row index
785  std::vector<std::size_t> ptr;
786 
787  // column index
788  std::vector<T> col;
789  };
790 
791 
792 
821  template <typename T = unsigned int>
823  {
824  Connectivity(const unsigned int dim,
825  const std::vector<::ReferenceCell> &cell_types)
826  : dim(dim)
828  {}
829 
830  inline std::vector<unsigned char> &
831  entity_orientations(const unsigned int structdim)
832  {
833  if (structdim == 1)
834  return line_orientation;
835 
836  AssertDimension(structdim, 2);
837 
838  return quad_orientation;
839  }
840 
841  inline const std::vector<unsigned char> &
842  entity_orientations(const unsigned int structdim) const
843  {
844  if (structdim == 1)
845  return line_orientation;
846 
847  AssertDimension(structdim, 2);
848 
849  return quad_orientation;
850  }
851 
852  inline std::vector<::ReferenceCell> &
853  entity_types(const unsigned int structdim)
854  {
855  if (structdim == dim)
856  return cell_types;
857 
858  // for vertices/lines the entity types are clear (0/1)
859  AssertDimension(structdim, 2);
860  AssertDimension(dim, 3);
861 
862  return quad_types;
863  }
864 
865  inline const std::vector<::ReferenceCell> &
866  entity_types(const unsigned int structdim) const
867  {
868  if (structdim == dim)
869  return cell_types;
870 
871  // for vertices/lines the entity types are clear (0/1)
872  AssertDimension(structdim, 2);
873  AssertDimension(dim, 3);
874 
875  return quad_types;
876  }
877 
878  inline CRS<T> &
879  entity_to_entities(const unsigned int from, const unsigned int to)
880  {
881  if (from == dim && to == dim)
882  return neighbors;
883  else if (from == dim && to == dim - 1)
884  return cell_entities;
885  else if (dim == 3 && from == 2 && to == 0)
886  return quad_vertices;
887  else if (dim == 3 && from == 2 && to == 1)
888  return quad_lines;
889  else if (from == 1 && to == 0)
890  return line_vertices;
891 
892  Assert(false, ExcNotImplemented());
893 
894  return cell_entities;
895  }
896 
897  inline const CRS<T> &
898  entity_to_entities(const unsigned int from, const unsigned int to) const
899  {
900  if (from == dim && to == dim)
901  return neighbors;
902  else if (from == dim && to == dim - 1)
903  return cell_entities;
904  else if (dim == 3 && from == 2 && to == 0)
905  return quad_vertices;
906  else if (dim == 3 && from == 2 && to == 1)
907  return quad_lines;
908  else if (from == 1 && to == 0)
909  return line_vertices;
910 
911  Assert(false, ExcNotImplemented());
912 
913  return cell_entities;
914  }
915 
916  private:
917  const unsigned int dim;
918  std::vector<::ReferenceCell> cell_types;
919 
921 
922  std::vector<unsigned char> line_orientation;
923 
926 
927  std::vector<unsigned char> quad_orientation;
928 
931 
932  std::vector<::ReferenceCell> quad_types;
933  };
934 
935 
936 
944  template <typename T>
945  void
946  determine_neighbors(const CRS<T> &con_cf, CRS<T> &con_cc)
947  {
948  const auto &col_cf = con_cf.col;
949  const auto &ptr_cf = con_cf.ptr;
950 
951  auto &col_cc = con_cc.col;
952  auto &ptr_cc = con_cc.ptr;
953 
954  const unsigned int n_faces =
955  *std::max_element(col_cf.begin(), col_cf.end()) + 1;
956 
957  // clear and initialize with -1 (assume that all faces are at the
958  // boundary)
959  col_cc = std::vector<T>(col_cf.size(), -1);
960  ptr_cc = ptr_cf;
961 
962  std::vector<std::pair<T, unsigned int>> neighbors(n_faces, {-1, -1});
963 
964  // loop over all cells
965  for (unsigned int i_0 = 0; i_0 < ptr_cf.size() - 1; ++i_0)
966  {
967  // ... and all its faces
968  for (std::size_t j_0 = ptr_cf[i_0]; j_0 < ptr_cf[i_0 + 1]; ++j_0)
969  {
970  if (neighbors[col_cf[j_0]].first == static_cast<unsigned int>(-1))
971  {
972  // face is visited the first time -> save the visiting cell
973  // and the face pointer
974  neighbors[col_cf[j_0]] = std::pair<T, unsigned int>(i_0, j_0);
975  }
976  else
977  {
978  // face is visited the second time -> now we know the cells
979  // on both sides of the face and we can determine for both
980  // cells the neighbor
981  col_cc[j_0] = neighbors[col_cf[j_0]].first;
982  col_cc[neighbors[col_cf[j_0]].second] = i_0;
983  }
984  }
985  }
986  }
987 
988 
989 
997  template <int key_length, typename FU>
998  void
1000  const unsigned int d,
1001  const std::vector<std::shared_ptr<CellTypeBase>> &cell_types,
1002  const std::vector<::ReferenceCell> & cell_types_index,
1003  const CRS<unsigned int> & crs,
1004  CRS<unsigned int> & crs_d, // result
1005  CRS<unsigned int> & crs_0, // result
1006  std::vector<unsigned char> & orientations, // result
1007  const FU & second_key_function)
1008  {
1009  const bool compatibility_mode = true;
1010 
1011  const std::vector<std::size_t> & cell_ptr = crs.ptr;
1012  const std::vector<unsigned int> &cell_vertices = crs.col;
1013  std::vector<std::size_t> & ptr_d = crs_d.ptr;
1014  std::vector<unsigned int> & col_d = crs_d.col;
1015 
1016  // note: we do not pre-allocate memory for these arrays because it turned
1017  // out that counting unique entities is more expensive than push_back().
1018  std::vector<std::size_t> & ptr_0 = crs_0.ptr;
1019  std::vector<unsigned int> &col_0 = crs_0.col;
1020 
1021  // clear
1022  ptr_0 = {};
1023  col_0 = {};
1024 
1025  unsigned int n_entities = 0;
1026 
1027  for (const auto &c : cell_types_index)
1028  n_entities +=
1029  cell_types[static_cast<types::geometric_entity_type>(c)]->n_entities(
1030  d);
1031 
1032  // step 1: store each d-dimensional entity of a cell (described by their
1033  // vertices) into a vector and create a key for them
1034  //
1035  // note: it turned out to be more efficient to have a vector of tuples
1036  // than to have two vectors (sorting becomes inefficient)
1037  std::vector<
1038  std::tuple<std::array<unsigned int, key_length>, unsigned int>>
1039  keys; // key (sorted vertices), cell-entity index
1040 
1041  std::vector<std::array<unsigned int, key_length>> ad_entity_vertices;
1042  std::vector<::ReferenceCell> ad_entity_types;
1043  std::vector<std::array<unsigned int, key_length>> ad_compatibility;
1044 
1045  keys.reserve(n_entities);
1046  ad_entity_vertices.reserve(n_entities);
1047  ad_entity_types.reserve(n_entities);
1048  ad_compatibility.reserve(n_entities);
1049 
1050  ptr_d.resize(cell_types_index.size() + 1);
1051  ptr_d[0] = 0;
1052 
1053  static const unsigned int offset = 1;
1054 
1055  // loop over all cells
1056  for (unsigned int c = 0, counter = 0; c < cell_types_index.size(); ++c)
1057  {
1058  const auto &cell_type =
1059  cell_types[static_cast<types::geometric_entity_type>(
1060  cell_types_index[c])];
1061  ptr_d[c + 1] = ptr_d[c] + cell_type->n_entities(d);
1062 
1063  // ... collect vertices of cell
1064  const ::ArrayView<const unsigned int> cell_vertice(
1065  cell_vertices.data() + cell_ptr[c], cell_ptr[c + 1] - cell_ptr[c]);
1066 
1067  // ... loop over all its entities
1068  for (unsigned int e = 0; e < cell_type->n_entities(d); ++e)
1069  {
1070  // ... determine global entity vertices
1071  const auto &local_entity_vertices =
1072  cell_type->vertices_of_entity(d, e);
1073 
1074  std::array<unsigned int, key_length> entity_vertices;
1075  std::fill(entity_vertices.begin(), entity_vertices.end(), 0);
1076 
1077  for (unsigned int i = 0; i < local_entity_vertices.size(); ++i)
1078  entity_vertices[i] =
1079  cell_vertice[local_entity_vertices[i]] + offset;
1080 
1081  // ... create key
1082  std::array<unsigned int, key_length> key = entity_vertices;
1083  std::sort(key.begin(), key.end());
1084  keys.emplace_back(key, counter++);
1085 
1086  ad_entity_vertices.emplace_back(entity_vertices);
1087 
1088  ad_entity_types.emplace_back(cell_type->type_of_entity(d, e));
1089 
1090  if (compatibility_mode)
1091  ad_compatibility.emplace_back(
1092  second_key_function(entity_vertices, cell_type, c, e));
1093  }
1094  }
1095 
1096  col_d.resize(keys.size());
1097  orientations.resize(keys.size());
1098 
1099  // step 2: sort according to key so that entities with same key can be
1100  // merged
1101  std::sort(keys.begin(), keys.end());
1102 
1103 
1104  if (compatibility_mode)
1105  {
1106  unsigned int n_unique_entities = 0;
1107  unsigned int n_unique_entity_vertices = 0;
1108 
1109  std::array<unsigned int, key_length> ref_key, new_key;
1110  std::fill(ref_key.begin(), ref_key.end(), 0);
1111  for (unsigned int i = 0; i < keys.size(); ++i)
1112  {
1113  const auto offset_i = std::get<1>(keys[i]);
1114 
1115  if (ref_key != std::get<0>(keys[i]))
1116  {
1117  ref_key = std::get<0>(keys[i]);
1118 
1119  n_unique_entities++;
1120  n_unique_entity_vertices +=
1121  cell_types[static_cast<types::geometric_entity_type>(
1122  ad_entity_types[offset_i])]
1123  ->n_entities(0);
1124 
1125  new_key = ad_compatibility[offset_i];
1126  }
1127 
1128  std::get<0>(keys[i]) = new_key;
1129  }
1130 
1131  std::sort(keys.begin(), keys.end());
1132 
1133  ptr_0.reserve(n_unique_entities);
1134  col_0.reserve(n_unique_entity_vertices);
1135  }
1136 
1137 
1138  std::array<unsigned int, key_length> ref_key;
1139  std::array<unsigned int, key_length> ref_indices;
1140  std::fill(ref_key.begin(), ref_key.end(), 0);
1141 
1142  for (unsigned int i = 0, counter = ::numbers::invalid_unsigned_int;
1143  i < keys.size();
1144  i++)
1145  {
1146  const auto offset_i = std::get<1>(keys[i]);
1147 
1148  if (ref_key != std::get<0>(keys[i]))
1149  {
1150  // new key
1151  counter++;
1152  ref_key = std::get<0>(keys[i]);
1153  ref_indices = ad_entity_vertices[offset_i];
1154 
1155  ptr_0.push_back(col_0.size());
1156  for (const auto j : ad_entity_vertices[offset_i])
1157  if (j != 0)
1158  col_0.push_back(j - offset);
1159 
1160  // take its orientation as default
1161  col_d[offset_i] = counter;
1162  orientations[offset_i] = 1;
1163  }
1164  else
1165  {
1166  col_d[offset_i] = counter;
1167  orientations[offset_i] =
1168  ad_entity_types[offset_i].compute_orientation(
1169  ad_entity_vertices[offset_i], ref_indices);
1170  }
1171  }
1172  ptr_0.push_back(col_0.size());
1173  }
1174 
1175 
1176 
1181  template <typename FU>
1182  void
1183  build_entity(const unsigned int d,
1184  const std::vector<std::shared_ptr<CellTypeBase>> &cell_types,
1185  const std::vector<::ReferenceCell> &cell_types_index,
1186  const CRS<unsigned int> & crs,
1187  CRS<unsigned int> & crs_d,
1188  CRS<unsigned int> & crs_0,
1189  std::vector<unsigned char> & orientations,
1190  const FU & second_key_function)
1191  {
1192  std::size_t key_length = 0;
1193 
1194  for (const auto &c : cell_types_index)
1195  {
1196  const auto &cell_type =
1197  cell_types[static_cast<types::geometric_entity_type>(c)];
1198  for (unsigned int e = 0; e < cell_type->n_entities(d); ++e)
1199  key_length =
1200  std::max(key_length, cell_type->vertices_of_entity(d, e).size());
1201  }
1202 
1203  if (key_length == 2)
1204  build_entity_templated<2>(d,
1205  cell_types,
1206  cell_types_index,
1207  crs,
1208  crs_d,
1209  crs_0,
1210  orientations,
1211  second_key_function);
1212  else if (key_length == 3)
1213  build_entity_templated<3>(d,
1214  cell_types,
1215  cell_types_index,
1216  crs,
1217  crs_d,
1218  crs_0,
1219  orientations,
1220  second_key_function);
1221  else if (key_length == 4)
1222  build_entity_templated<4>(d,
1223  cell_types,
1224  cell_types_index,
1225  crs,
1226  crs_d,
1227  crs_0,
1228  orientations,
1229  second_key_function);
1230  else
1232  }
1233 
1234 
1235 
1243  inline void
1245  const std::vector<std::shared_ptr<CellTypeBase>> &cell_types,
1246  const std::vector<::ReferenceCell> & cell_types_index,
1247  const CRS<unsigned int> & con_cv,
1248  const CRS<unsigned int> & con_cl,
1249  const CRS<unsigned int> & con_lv,
1250  const CRS<unsigned int> & con_cq,
1251  const CRS<unsigned int> & con_qv,
1252  const std::vector<unsigned char> & ori_cq,
1253  CRS<unsigned int> & con_ql, // result
1254  std::vector<unsigned char> & ori_ql, // result
1255  std::vector<::ReferenceCell> & quad_t_id // result
1256  )
1257  {
1258  // reset output
1259  ori_ql = {};
1260  con_ql.ptr = {};
1261  con_ql.col = {};
1262 
1263  con_ql.ptr.resize(con_qv.ptr.size());
1264  con_ql.ptr[0] = 0;
1265 
1266  quad_t_id.resize(con_qv.ptr.size() - 1);
1267 
1268  // count the number of lines of each face
1269  for (unsigned int c = 0; c < con_cq.ptr.size() - 1; ++c)
1270  {
1271  const auto &cell_type =
1272  cell_types[static_cast<types::geometric_entity_type>(
1273  cell_types_index[c])];
1274 
1275  // loop over faces
1276  for (unsigned int f_ = con_cq.ptr[c], f_index = 0;
1277  f_ < con_cq.ptr[c + 1];
1278  ++f_, ++f_index)
1279  {
1280  const unsigned int f = con_cq.col[f_];
1281 
1282  con_ql.ptr[f + 1] = cell_type->n_lines_of_surface(f_index);
1283  }
1284  }
1285 
1286  // use the counts to determine the offsets -> prefix sum
1287  for (unsigned int i = 0; i < con_ql.ptr.size() - 1; ++i)
1288  con_ql.ptr[i + 1] += con_ql.ptr[i];
1289 
1290  // allocate memory
1291  con_ql.col.resize(con_ql.ptr.back());
1292  ori_ql.resize(con_ql.ptr.back());
1293 
1294  // loop over cells
1295  for (unsigned int c = 0; c < con_cq.ptr.size() - 1; ++c)
1296  {
1297  const auto &cell_type =
1298  cell_types[static_cast<types::geometric_entity_type>(
1299  cell_types_index[c])];
1300 
1301  // loop over faces
1302  for (unsigned int f_ = con_cq.ptr[c], f_index = 0;
1303  f_ < con_cq.ptr[c + 1];
1304  ++f_, ++f_index)
1305  {
1306  const unsigned int f = con_cq.col[f_];
1307 
1308  // only faces with default orientation have to do something
1309  if (ori_cq[f_] != 1)
1310  continue;
1311 
1312  // determine entity type of face
1313  quad_t_id[f] = cell_type->type_of_entity(2, f_index);
1314 
1315  // loop over lines
1316  for (unsigned int l = 0;
1317  l < cell_type->n_lines_of_surface(f_index);
1318  ++l)
1319  {
1320  // determine global index of line
1321  const unsigned int local_line_index =
1322  cell_type->nth_line_of_surface(l, f_index);
1323  const unsigned int global_line_index =
1324  con_cl.col[con_cl.ptr[c] + local_line_index];
1325  con_ql.col[con_ql.ptr[f] + l] = global_line_index;
1326 
1327  // determine orientation of line
1328  const auto line_vertices_1_ref =
1329  cell_type->vertices_of_nth_line_of_surface(l, f_index);
1330 
1331  bool same = true;
1332  for (unsigned int v = 0; v < line_vertices_1_ref.size(); ++v)
1333  if (con_cv.col[con_cv.ptr[c] + line_vertices_1_ref[v]] !=
1334  con_lv.col[con_lv.ptr[global_line_index] + v])
1335  {
1336  same = false;
1337  break;
1338  }
1339 
1340  // ... comparison gives orientation
1341  ori_ql[con_ql.ptr[f] + l] = (same ? 1 : 0);
1342  }
1343  }
1344  }
1345  }
1346 
1347 
1348 
1357  template <typename T>
1358  Connectivity<T>
1359  build_connectivity(const unsigned int dim,
1360  const std::vector<std::shared_ptr<CellTypeBase>> &cell_t,
1361  const std::vector<::ReferenceCell> &cell_t_id,
1362  const CRS<T> & con_cv)
1363  {
1364  Connectivity<T> connectivity(dim, cell_t_id);
1365 
1366  CRS<T> temp1; // needed for 3D
1367 
1368  if (dim == 1)
1369  connectivity.entity_to_entities(1, 0) = con_cv;
1370 
1371  if (dim == 2 || dim == 3) // build lines
1372  {
1373  std::vector<unsigned char> dummy;
1374 
1375  build_entity(1,
1376  cell_t,
1377  connectivity.entity_types(dim),
1378  con_cv,
1379  dim == 2 ? connectivity.entity_to_entities(2, 1) : temp1,
1380  connectivity.entity_to_entities(1, 0),
1381  dim == 2 ? connectivity.entity_orientations(1) : dummy,
1382  [](auto key, const auto &, const auto &, const auto &) {
1383  // to ensure same enumeration as in deal.II
1384  return key;
1385  });
1386  }
1387 
1388  if (dim == 3) // build quads
1389  {
1390  build_entity(
1391  2,
1392  cell_t,
1393  connectivity.entity_types(3),
1394  con_cv,
1395  connectivity.entity_to_entities(3, 2),
1396  connectivity.entity_to_entities(2, 0),
1397  connectivity.entity_orientations(2),
1398  [&](auto key, const auto &cell_type, const auto &c, const auto &f) {
1399  // to ensure same enumeration as in deal.II
1400  AssertIndexRange(cell_type->n_lines_of_surface(f),
1401  key.size() + 1);
1402 
1403  unsigned int l = 0;
1404 
1405  for (; l < cell_type->n_lines_of_surface(f); ++l)
1406  key[l] =
1407  temp1
1408  .col[temp1.ptr[c] + cell_type->nth_line_of_surface(l, f)] +
1409  1 /*offset!*/;
1410 
1411  for (; l < key.size(); ++l)
1412  key[l] = 0;
1413 
1414  return key;
1415  });
1416 
1417  // create connectivity: quad -> line
1418  build_intersection(cell_t,
1419  connectivity.entity_types(3),
1420  con_cv,
1421  temp1,
1422  connectivity.entity_to_entities(1, 0),
1423  connectivity.entity_to_entities(3, 2),
1424  connectivity.entity_to_entities(2, 0),
1425  connectivity.entity_orientations(2),
1426  connectivity.entity_to_entities(2, 1),
1427  connectivity.entity_orientations(1),
1428  connectivity.entity_types(2));
1429  }
1430 
1431  // determine neighbors
1432  determine_neighbors(connectivity.entity_to_entities(dim, dim - 1),
1433  connectivity.entity_to_entities(dim, dim));
1434 
1435  return connectivity;
1436  }
1437 
1438 
1439 
1443  template <typename T, int dim>
1444  Connectivity<T>
1445  build_connectivity(const std::vector<CellData<dim>> &cells)
1446  {
1447  AssertThrow(cells.size() > 0, ExcMessage("No cells have been provided!"));
1448 
1449  // vector of possible cell entity types
1450  std::vector<std::shared_ptr<CellTypeBase>> cell_types_impl(8);
1451 
1452  cell_types_impl[static_cast<types::geometric_entity_type>(
1454  .reset(new CellTypeLine());
1455  cell_types_impl[static_cast<types::geometric_entity_type>(
1457  .reset(new CellTypeTri());
1458  cell_types_impl[static_cast<types::geometric_entity_type>(
1460  .reset(new CellTypeQuad());
1461  cell_types_impl[static_cast<types::geometric_entity_type>(
1463  .reset(new CellTypeTet());
1464  cell_types_impl[static_cast<types::geometric_entity_type>(
1466  .reset(new CellTypePyramid());
1467  cell_types_impl[static_cast<types::geometric_entity_type>(
1469  .reset(new CellTypeWedge());
1470  cell_types_impl[static_cast<types::geometric_entity_type>(
1472  .reset(new CellTypeHex());
1473 
1474  // determine cell types and process vertices
1475  std::vector<T> cell_vertices;
1476  cell_vertices.reserve(
1477  std::accumulate(cells.begin(),
1478  cells.end(),
1479  0,
1480  [](const auto &result, const auto &cell) {
1481  return result + cell.vertices.size();
1482  }));
1483 
1484  std::vector<std::size_t> cell_vertices_ptr;
1485  cell_vertices_ptr.reserve(cells.size() + 1);
1486  cell_vertices_ptr.push_back(0);
1487 
1488  std::vector<::ReferenceCell> cell_types_indices;
1489  cell_types_indices.reserve(cells.size());
1490 
1491  // loop over cells and create CRS
1492  for (const auto &cell : cells)
1493  {
1494 #ifdef DEBUG
1495  auto vertices_unique = cell.vertices;
1496  std::sort(vertices_unique.begin(), vertices_unique.end());
1497  vertices_unique.erase(std::unique(vertices_unique.begin(),
1498  vertices_unique.end()),
1499  vertices_unique.end());
1500 
1501  Assert(vertices_unique.size() == cell.vertices.size(),
1502  ExcMessage(
1503  "The definition of a cell refers to the same vertex several "
1504  "times. This is not possible. A common reason is that "
1505  "CellData::vertices has a size that does not match the "
1506  "size expected from the reference cell. Please resize "
1507  "CellData::vertices or use the appropriate constructor of "
1508  "CellData."));
1509 #endif
1510 
1511  const ::ReferenceCell reference_cell =
1513  cell.vertices.size());
1514 
1516  ExcNotImplemented());
1518  reference_cell),
1519  cell_types_impl.size());
1520  Assert(cell_types_impl[static_cast<types::geometric_entity_type>(
1521  reference_cell)]
1522  .get() != nullptr,
1523  ExcNotImplemented());
1524 
1525  cell_types_indices.push_back(reference_cell);
1526 
1527  // create CRS of vertices (to remove template argument dim)
1528  for (const auto &vertex : cell.vertices)
1529  cell_vertices.push_back(vertex);
1530 
1531  cell_vertices_ptr.push_back(cell_vertices.size());
1532  }
1533 
1534  // do the actual work
1535  return build_connectivity<T>(dim,
1536  cell_types_impl,
1537  cell_types_indices,
1538  {cell_vertices_ptr, cell_vertices});
1539  }
1540  } // namespace TriangulationImplementation
1541 } // namespace internal
1542 
1543 
1545 
1546 #endif
static ReferenceCell n_vertices_to_type(const int dim, const unsigned int n_vertices)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
Point< 2 > first
Definition: grid_out.cc:4604
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1667
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
graph_traits< Graph >::vertex_descriptor Vertex
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
constexpr const ReferenceCell Tetrahedron
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Pyramid
constexpr const ReferenceCell Invalid
constexpr const ReferenceCell Triangle
constexpr const ReferenceCell Hexahedron
constexpr const ReferenceCell Line
void determine_neighbors(const CRS< T > &con_cf, CRS< T > &con_cc)
Definition: connectivity.h:946
void build_intersection(const std::vector< std::shared_ptr< CellTypeBase >> &cell_types, const std::vector<::ReferenceCell > &cell_types_index, const CRS< unsigned int > &con_cv, const CRS< unsigned int > &con_cl, const CRS< unsigned int > &con_lv, const CRS< unsigned int > &con_cq, const CRS< unsigned int > &con_qv, const std::vector< unsigned char > &ori_cq, CRS< unsigned int > &con_ql, std::vector< unsigned char > &ori_ql, std::vector<::ReferenceCell > &quad_t_id)
void build_entity_templated(const unsigned int d, const std::vector< std::shared_ptr< CellTypeBase >> &cell_types, const std::vector<::ReferenceCell > &cell_types_index, const CRS< unsigned int > &crs, CRS< unsigned int > &crs_d, CRS< unsigned int > &crs_0, std::vector< unsigned char > &orientations, const FU &second_key_function)
Definition: connectivity.h:999
void build_entity(const unsigned int d, const std::vector< std::shared_ptr< CellTypeBase >> &cell_types, const std::vector<::ReferenceCell > &cell_types_index, const CRS< unsigned int > &crs, CRS< unsigned int > &crs_d, CRS< unsigned int > &crs_0, std::vector< unsigned char > &orientations, const FU &second_key_function)
Connectivity< T > build_connectivity(const unsigned int dim, const std::vector< std::shared_ptr< CellTypeBase >> &cell_t, const std::vector<::ReferenceCell > &cell_t_id, const CRS< T > &con_cv)
static const unsigned int invalid_unsigned_int
Definition: types.h:201
std::uint8_t geometric_entity_type
Definition: types.h:157
CRS(const std::vector< std::size_t > &ptr, const std::vector< T > &col)
Definition: connectivity.h:779
virtual unsigned int n_entities(const unsigned int d) const
Definition: connectivity.h:50
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const
Definition: connectivity.h:75
virtual ::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const
Definition: connectivity.h:62
virtual const std::array< unsigned int, 2 > & vertices_of_nth_line_of_surface(const unsigned int line, const unsigned int face) const
Definition: connectivity.h:114
virtual unsigned int n_lines_of_surface(const unsigned int face) const
Definition: connectivity.h:88
virtual unsigned int nth_line_of_surface(const unsigned int line, const unsigned int face) const
Definition: connectivity.h:100
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:696
const std::array< unsigned int, 2 > & vertices_of_nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:744
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:715
unsigned int nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:729
unsigned int n_lines_of_surface(const unsigned int surface) const override
Definition: connectivity.h:722
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:645
::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:155
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:168
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:135
unsigned int n_lines_of_surface(const unsigned int surface) const override
Definition: connectivity.h:470
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:442
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:395
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:463
unsigned int nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:479
const std::array< unsigned int, 2 > & vertices_of_nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:493
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:240
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:282
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:266
unsigned int nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:364
unsigned int n_lines_of_surface(const unsigned int line) const override
Definition: connectivity.h:357
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:331
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:350
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:297
const std::array< unsigned int, 2 > & vertices_of_nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:374
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:225
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:183
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:209
unsigned int n_entities(const unsigned int d) const override
Definition: connectivity.h:589
unsigned int n_lines_of_surface(const unsigned int surface) const override
Definition: connectivity.h:596
::ArrayView< const unsigned int > vertices_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:517
virtual ::ReferenceCell type_of_entity(const unsigned int d, const unsigned int e) const override
Definition: connectivity.h:568
unsigned int nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:605
const std::array< unsigned int, 2 > & vertices_of_nth_line_of_surface(const unsigned int line, const unsigned int face) const override
Definition: connectivity.h:621
std::vector<::ReferenceCell > & entity_types(const unsigned int structdim)
Definition: connectivity.h:853
Connectivity(const unsigned int dim, const std::vector<::ReferenceCell > &cell_types)
Definition: connectivity.h:824
const std::vector< unsigned char > & entity_orientations(const unsigned int structdim) const
Definition: connectivity.h:842
const CRS< T > & entity_to_entities(const unsigned int from, const unsigned int to) const
Definition: connectivity.h:898
std::vector< unsigned char > & entity_orientations(const unsigned int structdim)
Definition: connectivity.h:831
const std::vector<::ReferenceCell > & entity_types(const unsigned int structdim) const
Definition: connectivity.h:866
CRS< T > & entity_to_entities(const unsigned int from, const unsigned int to)
Definition: connectivity.h:879