1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 2020 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_reference_cell_h 17 #define dealii_tria_reference_cell_h 18 19 #include <deal.II/base/config.h> 20 21 #include <deal.II/base/geometry_info.h> 22 23 24 DEAL_II_NAMESPACE_OPEN 25 26 27 /** 28 * A namespace for reference cells. 29 */ 30 namespace ReferenceCell 31 { 32 /** 33 * Supported reference cell types. 34 */ 35 enum class Type : std::uint8_t 36 { 37 Vertex = 0, 38 Line = 1, 39 Tri = 2, 40 Quad = 3, 41 Tet = 4, 42 Pyramid = 5, 43 Wedge = 6, 44 Hex = 7, 45 Invalid = static_cast<std::uint8_t>(-1) 46 }; 47 48 /** 49 * Return the correct simplex reference cell type for the given dimension 50 * @p dim. 51 */ 52 inline Type get_simplex(const unsigned int dim)53 get_simplex(const unsigned int dim) 54 { 55 switch (dim) 56 { 57 case 0: 58 return Type::Vertex; 59 case 1: 60 return Type::Line; 61 case 2: 62 return Type::Tri; 63 case 3: 64 return Type::Tet; 65 default: 66 Assert(false, ExcNotImplemented()); 67 return Type::Invalid; 68 } 69 } 70 71 /** 72 * Return the correct hypercube reference cell type for the given dimension 73 * @p dim. 74 */ 75 inline Type get_hypercube(const unsigned int dim)76 get_hypercube(const unsigned int dim) 77 { 78 switch (dim) 79 { 80 case 0: 81 return Type::Vertex; 82 case 1: 83 return Type::Line; 84 case 2: 85 return Type::Quad; 86 case 3: 87 return Type::Hex; 88 default: 89 Assert(false, ExcNotImplemented()); 90 return Type::Invalid; 91 } 92 } 93 94 namespace internal 95 { 96 /** 97 * Check if the bit at position @p n in @p number is set. 98 */ 99 inline static bool get_bit(const unsigned char number,const unsigned int n)100 get_bit(const unsigned char number, const unsigned int n) 101 { 102 AssertIndexRange(n, 8); 103 104 // source: 105 // https://stackoverflow.com/questions/47981/how-do-you-set-clear-and-toggle-a-single-bit 106 // "Checking a bit" 107 return (number >> n) & 1U; 108 } 109 110 111 112 /** 113 * Set the bit at position @p n in @p number to value @p x. 114 */ 115 inline static void set_bit(unsigned char & number,const unsigned int n,const bool x)116 set_bit(unsigned char &number, const unsigned int n, const bool x) 117 { 118 AssertIndexRange(n, 8); 119 120 // source: 121 // https://stackoverflow.com/questions/47981/how-do-you-set-clear-and-toggle-a-single-bit 122 // "Changing the nth bit to x" 123 number ^= (-static_cast<unsigned char>(x) ^ number) & (1U << n); 124 } 125 126 /** 127 * A namespace for geometric information on reference cells. 128 */ 129 namespace Info 130 { 131 /** 132 * Interface to be used in TriaAccessor/TriaCellAccessor to access 133 * sub-entities of dimension d' of geometric entities of dimension d, with 134 * 0<=d'<d<=3. 135 */ 136 struct Base 137 { 138 /** 139 * Destructor. 140 */ 141 virtual ~Base() = default; 142 143 /** 144 * Number of vertices. 145 */ 146 virtual unsigned int n_verticesBase147 n_vertices() const 148 { 149 Assert(false, ExcNotImplemented()); 150 return 0; 151 } 152 153 /** 154 * Number of lines. 155 */ 156 virtual unsigned int n_linesBase157 n_lines() const 158 { 159 Assert(false, ExcNotImplemented()); 160 return 0; 161 } 162 163 164 /** 165 * Number of faces. 166 */ 167 virtual unsigned int n_facesBase168 n_faces() const 169 { 170 Assert(false, ExcNotImplemented()); 171 return 0; 172 } 173 174 /** 175 * Return an object that can be thought of as an array containing all 176 * indices from zero to n_vertices(). 177 */ 178 inline std_cxx20::ranges::iota_view<unsigned int, unsigned int> vertex_indicesBase179 vertex_indices() const 180 { 181 return {0U, n_vertices()}; 182 } 183 184 /** 185 * Return an object that can be thought of as an array containing all 186 * indices from zero to n_lines(). 187 */ 188 inline std_cxx20::ranges::iota_view<unsigned int, unsigned int> line_indicesBase189 line_indices() const 190 { 191 return {0U, n_lines()}; 192 } 193 194 /** 195 * Return an object that can be thought of as an array containing all 196 * indices from zero to n_faces(). 197 */ 198 inline std_cxx20::ranges::iota_view<unsigned int, unsigned int> face_indicesBase199 face_indices() const 200 { 201 return {0U, n_faces()}; 202 } 203 204 /** 205 * Standard decomposition of vertex index into face and face-vertex 206 * index. 207 */ 208 virtual std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexBase209 standard_vertex_to_face_and_vertex_index( 210 const unsigned int vertex) const 211 { 212 Assert(false, ExcNotImplemented()); 213 214 (void)vertex; 215 216 return {{0u, 0u}}; 217 } 218 219 /** 220 * Standard decomposition of line index into face and face-line index. 221 */ 222 virtual std::array<unsigned int, 2> standard_line_to_face_and_line_indexBase223 standard_line_to_face_and_line_index(const unsigned int line) const 224 { 225 Assert(false, ExcNotImplemented()); 226 227 (void)line; 228 229 return {{0, 0}}; 230 } 231 232 /** 233 * Correct vertex index depending on face orientation. 234 */ 235 virtual unsigned int standard_to_real_face_vertexBase236 standard_to_real_face_vertex(const unsigned int vertex, 237 const unsigned int face, 238 const unsigned char face_orientation) const 239 { 240 Assert(false, ExcNotImplemented()); 241 242 (void)vertex; 243 (void)face; 244 (void)face_orientation; 245 246 return 0; 247 } 248 249 /** 250 * Correct line index depending on face orientation. 251 */ 252 virtual unsigned int standard_to_real_face_lineBase253 standard_to_real_face_line(const unsigned int line, 254 const unsigned int face, 255 const unsigned char face_orientation) const 256 { 257 Assert(false, ExcNotImplemented()); 258 259 (void)line; 260 (void)face; 261 (void)face_orientation; 262 263 return 0; 264 } 265 266 /** 267 * Combine face and line orientation. 268 */ 269 virtual bool combine_face_and_line_orientationBase270 combine_face_and_line_orientation( 271 const unsigned int line, 272 const unsigned char face_orientation, 273 const unsigned char line_orientation) const 274 { 275 Assert(false, ExcNotImplemented()); 276 277 (void)line; 278 (void)face_orientation; 279 (void)line_orientation; 280 281 return true; 282 } 283 284 /** 285 * Return reference-cell type of face @p face_no. 286 */ 287 virtual ReferenceCell::Type face_reference_cell_typeBase288 face_reference_cell_type(const unsigned int face_no) const 289 { 290 Assert(false, ExcNotImplemented()); 291 (void)face_no; 292 293 return ReferenceCell::Type::Invalid; 294 } 295 296 /** 297 * Map face line number to cell line number. 298 */ 299 virtual unsigned int face_to_cell_linesBase300 face_to_cell_lines(const unsigned int face, 301 const unsigned int line, 302 const unsigned char face_orientation) const 303 { 304 Assert(false, ExcNotImplemented()); 305 (void)face; 306 (void)line; 307 (void)face_orientation; 308 309 return 0; 310 } 311 312 /** 313 * Map face vertex number to cell vertex number. 314 */ 315 virtual unsigned int face_to_cell_verticesBase316 face_to_cell_vertices(const unsigned int face, 317 const unsigned int vertex, 318 const unsigned char face_orientation) const 319 { 320 Assert(false, ExcNotImplemented()); 321 (void)face; 322 (void)vertex; 323 (void)face_orientation; 324 325 return 0; 326 } 327 }; 328 329 330 /** 331 * Base class for tensor-product geometric entities. 332 */ 333 template <int dim> 334 struct TensorProductBase : Base 335 { 336 unsigned int n_verticesTensorProductBase337 n_vertices() const override 338 { 339 return GeometryInfo<dim>::vertices_per_cell; 340 } 341 342 unsigned int n_linesTensorProductBase343 n_lines() const override 344 { 345 return GeometryInfo<dim>::lines_per_cell; 346 } 347 348 unsigned int n_facesTensorProductBase349 n_faces() const override 350 { 351 return GeometryInfo<dim>::faces_per_cell; 352 } 353 354 unsigned int face_to_cell_linesTensorProductBase355 face_to_cell_lines(const unsigned int face, 356 const unsigned int line, 357 const unsigned char face_orientation) const override 358 { 359 return GeometryInfo<dim>::face_to_cell_lines( 360 face, 361 line, 362 get_bit(face_orientation, 0), 363 get_bit(face_orientation, 2), 364 get_bit(face_orientation, 1)); 365 } 366 367 unsigned int face_to_cell_verticesTensorProductBase368 face_to_cell_vertices( 369 const unsigned int face, 370 const unsigned int vertex, 371 const unsigned char face_orientation) const override 372 { 373 return GeometryInfo<dim>::face_to_cell_vertices( 374 face, 375 vertex, 376 get_bit(face_orientation, 0), 377 get_bit(face_orientation, 2), 378 get_bit(face_orientation, 1)); 379 } 380 }; 381 382 383 384 /* 385 * Vertex. 386 */ 387 struct Vertex : public TensorProductBase<0> 388 { 389 ReferenceCell::Type face_reference_cell_typeVertex390 face_reference_cell_type(const unsigned int face_no) const override 391 { 392 (void)face_no; 393 return ReferenceCell::Type::Invalid; 394 } 395 }; 396 397 398 399 /* 400 * Line. 401 */ 402 struct Line : public TensorProductBase<1> 403 { 404 ReferenceCell::Type face_reference_cell_typeLine405 face_reference_cell_type(const unsigned int face_no) const override 406 { 407 (void)face_no; 408 return ReferenceCell::Type::Vertex; 409 } 410 }; 411 412 413 414 /** 415 * Triangle. 416 */ 417 struct Tri : public Base 418 { 419 unsigned int n_verticesTri420 n_vertices() const override 421 { 422 return 3; 423 } 424 425 unsigned int n_linesTri426 n_lines() const override 427 { 428 return 3; 429 } 430 431 unsigned int n_facesTri432 n_faces() const override 433 { 434 return this->n_lines(); 435 } 436 437 std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexTri438 standard_vertex_to_face_and_vertex_index( 439 const unsigned int vertex) const override 440 { 441 AssertIndexRange(vertex, 3); 442 443 static const std::array<std::array<unsigned int, 2>, 3> table = { 444 {{{0, 0}}, {{0, 1}}, {{1, 1}}}}; 445 446 return table[vertex]; 447 } 448 449 unsigned int standard_to_real_face_vertexTri450 standard_to_real_face_vertex( 451 const unsigned int vertex, 452 const unsigned int face, 453 const unsigned char line_orientation) const override 454 { 455 (void)face; 456 457 static const std::array<std::array<unsigned int, 2>, 2> table = { 458 {{{1, 0}}, {{0, 1}}}}; 459 460 return table[line_orientation][vertex]; 461 } 462 463 ReferenceCell::Type face_reference_cell_typeTri464 face_reference_cell_type(const unsigned int face_no) const override 465 { 466 (void)face_no; 467 468 AssertIndexRange(face_no, n_faces()); 469 470 return ReferenceCell::Type::Line; 471 } 472 473 unsigned int face_to_cell_linesTri474 face_to_cell_lines(const unsigned int face, 475 const unsigned int line, 476 const unsigned char face_orientation) const override 477 { 478 AssertIndexRange(face, n_faces()); 479 AssertDimension(line, 0); 480 481 (void)line; 482 (void)face_orientation; 483 484 return face; 485 } 486 487 unsigned int face_to_cell_verticesTri488 face_to_cell_vertices( 489 const unsigned int face, 490 const unsigned int vertex, 491 const unsigned char face_orientation) const override 492 { 493 static const std::array<std::array<unsigned int, 2>, 3> table = { 494 {{{0, 1}}, {{1, 2}}, {{2, 0}}}}; 495 496 return table[face][face_orientation ? vertex : (1 - vertex)]; 497 } 498 }; 499 500 501 502 /** 503 * Quad. 504 */ 505 struct Quad : public TensorProductBase<2> 506 { 507 std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexQuad508 standard_vertex_to_face_and_vertex_index( 509 const unsigned int vertex) const override 510 { 511 return GeometryInfo<2>::standard_quad_vertex_to_line_vertex_index( 512 vertex); 513 } 514 515 unsigned int standard_to_real_face_vertexQuad516 standard_to_real_face_vertex( 517 const unsigned int vertex, 518 const unsigned int face, 519 const unsigned char line_orientation) const override 520 { 521 (void)face; 522 523 return GeometryInfo<2>::standard_to_real_line_vertex( 524 vertex, line_orientation); 525 } 526 527 ReferenceCell::Type face_reference_cell_typeQuad528 face_reference_cell_type(const unsigned int face_no) const override 529 { 530 (void)face_no; 531 return ReferenceCell::Type::Line; 532 } 533 }; 534 535 536 537 /** 538 * Tet. 539 */ 540 struct Tet : public TensorProductBase<3> 541 { 542 unsigned int n_verticesTet543 n_vertices() const override 544 { 545 return 4; 546 } 547 548 unsigned int n_linesTet549 n_lines() const override 550 { 551 return 6; 552 } 553 554 unsigned int n_facesTet555 n_faces() const override 556 { 557 return 4; 558 } 559 560 std::array<unsigned int, 2> standard_line_to_face_and_line_indexTet561 standard_line_to_face_and_line_index( 562 const unsigned int line) const override 563 { 564 static const std::array<unsigned int, 2> table[6] = { 565 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}}; 566 567 return table[line]; 568 } 569 570 unsigned int standard_to_real_face_lineTet571 standard_to_real_face_line( 572 const unsigned int line, 573 const unsigned int face, 574 const unsigned char face_orientation) const override 575 { 576 (void)face; 577 578 static const std::array<std::array<unsigned int, 3>, 6> table = { 579 {{{2, 1, 0}}, 580 {{0, 1, 2}}, 581 {{1, 2, 0}}, 582 {{0, 2, 1}}, 583 {{1, 0, 2}}, 584 {{2, 0, 1}}}}; 585 586 return table[face_orientation][line]; 587 } 588 589 bool combine_face_and_line_orientationTet590 combine_face_and_line_orientation( 591 const unsigned int line, 592 const unsigned char face_orientation_raw, 593 const unsigned char line_orientation) const override 594 { 595 (void)line; 596 (void)face_orientation_raw; 597 598 return line_orientation; 599 } 600 601 std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexTet602 standard_vertex_to_face_and_vertex_index( 603 const unsigned int vertex) const override 604 { 605 AssertIndexRange(vertex, 4); 606 607 static const std::array<unsigned int, 2> table[4] = {{{0, 0}}, 608 {{0, 1}}, 609 {{0, 2}}, 610 {{1, 2}}}; 611 612 return table[vertex]; 613 } 614 615 unsigned int standard_to_real_face_vertexTet616 standard_to_real_face_vertex( 617 const unsigned int vertex, 618 const unsigned int face, 619 const unsigned char face_orientation) const override 620 { 621 AssertIndexRange(face_orientation, 6); 622 (void)face; 623 624 static const std::array<std::array<unsigned int, 3>, 6> table = { 625 {{{0, 2, 1}}, 626 {{0, 1, 2}}, 627 {{1, 2, 0}}, 628 {{1, 0, 2}}, 629 {{2, 1, 0}}, 630 {{2, 0, 1}}}}; 631 632 return table[face_orientation][vertex]; 633 } 634 635 ReferenceCell::Type face_reference_cell_typeTet636 face_reference_cell_type(const unsigned int face_no) const override 637 { 638 (void)face_no; 639 640 AssertIndexRange(face_no, n_faces()); 641 642 return ReferenceCell::Type::Tri; 643 } 644 645 unsigned int face_to_cell_linesTet646 face_to_cell_lines(const unsigned int face, 647 const unsigned int line, 648 const unsigned char face_orientation) const override 649 { 650 AssertIndexRange(face, n_faces()); 651 652 const static std::array<std::array<unsigned int, 3>, 4> table = { 653 {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}}; 654 655 return table[face][standard_to_real_face_line( 656 line, face, face_orientation)]; 657 } 658 659 unsigned int face_to_cell_verticesTet660 face_to_cell_vertices( 661 const unsigned int face, 662 const unsigned int vertex, 663 const unsigned char face_orientation) const override 664 { 665 static const std::array<std::array<unsigned int, 3>, 4> table = { 666 {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}}; 667 668 return table[face][standard_to_real_face_vertex( 669 vertex, face, face_orientation)]; 670 } 671 }; 672 673 674 675 /** 676 * Pyramid. 677 */ 678 struct Pyramid : public TensorProductBase<3> 679 { 680 unsigned int n_verticesPyramid681 n_vertices() const override 682 { 683 return 5; 684 } 685 686 unsigned int n_linesPyramid687 n_lines() const override 688 { 689 return 8; 690 } 691 692 unsigned int n_facesPyramid693 n_faces() const override 694 { 695 return 5; 696 } 697 698 std::array<unsigned int, 2> standard_line_to_face_and_line_indexPyramid699 standard_line_to_face_and_line_index( 700 const unsigned int line) const override 701 { 702 Assert(false, ExcNotImplemented()); 703 704 static const std::array<unsigned int, 2> table[6] = { 705 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}}; 706 707 return table[line]; 708 } 709 710 unsigned int standard_to_real_face_linePyramid711 standard_to_real_face_line( 712 const unsigned int line, 713 const unsigned int face, 714 const unsigned char face_orientation) const override 715 { 716 Assert(false, ExcNotImplemented()); 717 718 (void)face; 719 720 static const std::array<std::array<unsigned int, 3>, 6> table = { 721 {{{2, 1, 0}}, 722 {{0, 1, 2}}, 723 {{1, 2, 0}}, 724 {{0, 2, 1}}, 725 {{1, 0, 2}}, 726 {{2, 0, 1}}}}; 727 728 return table[face_orientation][line]; 729 } 730 731 bool combine_face_and_line_orientationPyramid732 combine_face_and_line_orientation( 733 const unsigned int line, 734 const unsigned char face_orientation_raw, 735 const unsigned char line_orientation) const override 736 { 737 (void)line; 738 (void)face_orientation_raw; 739 740 return line_orientation; 741 } 742 743 std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexPyramid744 standard_vertex_to_face_and_vertex_index( 745 const unsigned int vertex) const override 746 { 747 static const std::array<unsigned int, 2> table[5] = { 748 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}}; 749 750 return table[vertex]; 751 } 752 753 unsigned int standard_to_real_face_vertexPyramid754 standard_to_real_face_vertex( 755 const unsigned int vertex, 756 const unsigned int face, 757 const unsigned char face_orientation) const override 758 { 759 if (face == 0) // Quad 760 { 761 return GeometryInfo<3>::standard_to_real_face_vertex( 762 vertex, 763 get_bit(face_orientation, 0), 764 get_bit(face_orientation, 2), 765 get_bit(face_orientation, 1)); 766 } 767 else // Tri 768 { 769 static const std::array<std::array<unsigned int, 3>, 6> table = { 770 {{{0, 2, 1}}, 771 {{0, 1, 2}}, 772 {{1, 2, 0}}, 773 {{1, 0, 2}}, 774 {{2, 1, 0}}, 775 {{2, 0, 1}}}}; 776 777 return table[face_orientation][vertex]; 778 } 779 } 780 781 ReferenceCell::Type face_reference_cell_typePyramid782 face_reference_cell_type(const unsigned int face_no) const override 783 { 784 AssertIndexRange(face_no, n_faces()); 785 786 if (face_no == 1) 787 return ReferenceCell::Type::Quad; 788 else 789 return ReferenceCell::Type::Tri; 790 } 791 }; 792 793 794 795 /** 796 * Wedge. 797 */ 798 struct Wedge : public TensorProductBase<3> 799 { 800 unsigned int n_verticesWedge801 n_vertices() const override 802 { 803 return 6; 804 } 805 806 unsigned int n_linesWedge807 n_lines() const override 808 { 809 return 9; 810 } 811 812 unsigned int n_facesWedge813 n_faces() const override 814 { 815 return 6; 816 } 817 818 std::array<unsigned int, 2> standard_line_to_face_and_line_indexWedge819 standard_line_to_face_and_line_index( 820 const unsigned int line) const override 821 { 822 Assert(false, ExcNotImplemented()); 823 824 static const std::array<unsigned int, 2> table[6] = { 825 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}}; 826 827 return table[line]; 828 } 829 830 unsigned int standard_to_real_face_lineWedge831 standard_to_real_face_line( 832 const unsigned int line, 833 const unsigned int face, 834 const unsigned char face_orientation) const override 835 { 836 Assert(false, ExcNotImplemented()); 837 838 (void)face; 839 840 static const std::array<std::array<unsigned int, 3>, 6> table = { 841 {{{2, 1, 0}}, 842 {{0, 1, 2}}, 843 {{1, 2, 0}}, 844 {{0, 2, 1}}, 845 {{1, 0, 2}}, 846 {{2, 0, 1}}}}; 847 848 return table[face_orientation][line]; 849 } 850 851 bool combine_face_and_line_orientationWedge852 combine_face_and_line_orientation( 853 const unsigned int line, 854 const unsigned char face_orientation_raw, 855 const unsigned char line_orientation) const override 856 { 857 (void)line; 858 (void)face_orientation_raw; 859 860 return line_orientation; 861 } 862 863 std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexWedge864 standard_vertex_to_face_and_vertex_index( 865 const unsigned int vertex) const override 866 { 867 static const std::array<std::array<unsigned int, 2>, 6> table = { 868 {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}}; 869 870 return table[vertex]; 871 } 872 873 unsigned int standard_to_real_face_vertexWedge874 standard_to_real_face_vertex( 875 const unsigned int vertex, 876 const unsigned int face, 877 const unsigned char face_orientation) const override 878 { 879 if (face > 1) // QUAD 880 { 881 return GeometryInfo<3>::standard_to_real_face_vertex( 882 vertex, 883 get_bit(face_orientation, 0), 884 get_bit(face_orientation, 2), 885 get_bit(face_orientation, 1)); 886 } 887 else // TRI 888 { 889 static const std::array<std::array<unsigned int, 3>, 6> table = { 890 {{{0, 2, 1}}, 891 {{0, 1, 2}}, 892 {{1, 2, 0}}, 893 {{1, 0, 2}}, 894 {{2, 1, 0}}, 895 {{2, 0, 1}}}}; 896 897 return table[face_orientation][vertex]; 898 } 899 } 900 901 ReferenceCell::Type face_reference_cell_typeWedge902 face_reference_cell_type(const unsigned int face_no) const override 903 { 904 AssertIndexRange(face_no, n_faces()); 905 906 if (face_no > 1) 907 return ReferenceCell::Type::Quad; 908 else 909 return ReferenceCell::Type::Tri; 910 } 911 }; 912 913 914 915 /** 916 * Hex. 917 */ 918 struct Hex : public TensorProductBase<3> 919 { 920 std::array<unsigned int, 2> standard_line_to_face_and_line_indexHex921 standard_line_to_face_and_line_index( 922 const unsigned int line) const override 923 { 924 return GeometryInfo<3>::standard_hex_line_to_quad_line_index(line); 925 } 926 927 unsigned int standard_to_real_face_lineHex928 standard_to_real_face_line( 929 const unsigned int line, 930 const unsigned int face, 931 const unsigned char face_orientation) const override 932 { 933 (void)face; 934 935 return GeometryInfo<3>::standard_to_real_face_line( 936 line, 937 get_bit(face_orientation, 0), 938 get_bit(face_orientation, 2), 939 get_bit(face_orientation, 1)); 940 } 941 942 bool combine_face_and_line_orientationHex943 combine_face_and_line_orientation( 944 const unsigned int line, 945 const unsigned char face_orientation_raw, 946 const unsigned char line_orientation) const override 947 { 948 static const bool bool_table[2][2][2][2] = { 949 {{{true, false}, // lines 0/1, face_orientation=false, 950 // face_flip=false, face_rotation=false and true 951 {false, true}}, // lines 0/1, face_orientation=false, 952 // face_flip=true, face_rotation=false and true 953 {{true, true}, // lines 0/1, face_orientation=true, 954 // face_flip=false, face_rotation=false and true 955 {false, false}}}, // lines 0/1, face_orientation=true, 956 // face_flip=true, face_rotation=false and true 957 958 {{{true, true}, // lines 2/3 ... 959 {false, false}}, 960 {{true, false}, {false, true}}}}; 961 962 const bool face_orientation = get_bit(face_orientation_raw, 0); 963 const bool face_flip = get_bit(face_orientation_raw, 2); 964 const bool face_rotation = get_bit(face_orientation_raw, 1); 965 966 return ( 967 static_cast<bool>(line_orientation) == 968 bool_table[line / 2][face_orientation][face_flip][face_rotation]); 969 } 970 971 std::array<unsigned int, 2> standard_vertex_to_face_and_vertex_indexHex972 standard_vertex_to_face_and_vertex_index( 973 const unsigned int vertex) const override 974 { 975 return GeometryInfo<3>::standard_hex_vertex_to_quad_vertex_index( 976 vertex); 977 } 978 979 unsigned int standard_to_real_face_vertexHex980 standard_to_real_face_vertex( 981 const unsigned int vertex, 982 const unsigned int face, 983 const unsigned char face_orientation) const override 984 { 985 (void)face; 986 987 return GeometryInfo<3>::standard_to_real_face_vertex( 988 vertex, 989 get_bit(face_orientation, 0), 990 get_bit(face_orientation, 2), 991 get_bit(face_orientation, 1)); 992 } 993 994 ReferenceCell::Type face_reference_cell_typeHex995 face_reference_cell_type(const unsigned int face_no) const override 996 { 997 (void)face_no; 998 return ReferenceCell::Type::Quad; 999 } 1000 }; 1001 1002 /** 1003 * Return for a given reference-cell type the right Info. 1004 */ 1005 inline const ReferenceCell::internal::Info::Base & get_cell(const ReferenceCell::Type & type)1006 get_cell(const ReferenceCell::Type &type) 1007 { 1008 static const std:: 1009 array<std::unique_ptr<ReferenceCell::internal::Info::Base>, 8> 1010 gei{{std::make_unique<ReferenceCell::internal::Info::Vertex>(), 1011 std::make_unique<ReferenceCell::internal::Info::Line>(), 1012 std::make_unique<ReferenceCell::internal::Info::Tri>(), 1013 std::make_unique<ReferenceCell::internal::Info::Quad>(), 1014 std::make_unique<ReferenceCell::internal::Info::Tet>(), 1015 std::make_unique<ReferenceCell::internal::Info::Pyramid>(), 1016 std::make_unique<ReferenceCell::internal::Info::Wedge>(), 1017 std::make_unique<ReferenceCell::internal::Info::Hex>()}}; 1018 AssertIndexRange(static_cast<std::uint8_t>(type), 8); 1019 return *gei[static_cast<std::uint8_t>(type)]; 1020 } 1021 1022 /** 1023 * Return for a given reference-cell type @p and face number @p face_no the 1024 * right Info of the @p face_no-th face. 1025 */ 1026 inline const ReferenceCell::internal::Info::Base & get_face(const ReferenceCell::Type & type,const unsigned int face_no)1027 get_face(const ReferenceCell::Type &type, const unsigned int face_no) 1028 { 1029 return get_cell(get_cell(type).face_reference_cell_type(face_no)); 1030 } 1031 1032 } // namespace Info 1033 } // namespace internal 1034 } // namespace ReferenceCell 1035 1036 1037 DEAL_II_NAMESPACE_CLOSE 1038 1039 #endif 1040