From e2181a4e1031624d496feccb62baa536eef86055 Mon Sep 17 00:00:00 2001 From: Siarhei Uzunbajakau Date: Mon, 7 Oct 2024 13:08:13 +0200 Subject: [PATCH] Nedelec orientation fix nr.1. Update fe_nedelec.h Indent errors fix A comment correction Moving the swap tables to the source file Resorting to c-tyle arrays. Nedelec face edge orientation test Nedelec face edge orientation test Nedelec face edge orientation test Nedelec face edge orientation test Nedelec face edge orientation test Changed misleading variable name Changes after review --- include/deal.II/fe/fe_nedelec.h | 5 +- source/fe/fe_nedelec.cc | 356 +++++++++++++- source/fe/fe_poly_tensor.cc | 135 +++++- tests/fe/fe_nedelec_face_edge_orientation.h | 451 ++++++++++++++++++ .../fe/fe_nedelec_face_edge_orientation_0.cc | 32 ++ .../fe_nedelec_face_edge_orientation_0.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_1.cc | 32 ++ .../fe_nedelec_face_edge_orientation_1.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_2.cc | 32 ++ .../fe_nedelec_face_edge_orientation_2.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_3.cc | 32 ++ .../fe_nedelec_face_edge_orientation_3.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_4.cc | 32 ++ .../fe_nedelec_face_edge_orientation_4.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_5.cc | 32 ++ .../fe_nedelec_face_edge_orientation_5.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_6.cc | 32 ++ .../fe_nedelec_face_edge_orientation_6.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_7.cc | 32 ++ .../fe_nedelec_face_edge_orientation_7.output | 101 ++++ .../fe/fe_nedelec_face_edge_orientation_8.cc | 41 ++ .../fe_nedelec_face_edge_orientation_8.output | 41 ++ 22 files changed, 2069 insertions(+), 24 deletions(-) create mode 100644 tests/fe/fe_nedelec_face_edge_orientation.h create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_0.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_0.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_1.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_1.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_2.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_2.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_3.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_3.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_4.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_4.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_5.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_5.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_6.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_6.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_7.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_7.output create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_8.cc create mode 100644 tests/fe/fe_nedelec_face_edge_orientation_8.output diff --git a/include/deal.II/fe/fe_nedelec.h b/include/deal.II/fe/fe_nedelec.h index 67fbcb7d54..078cb14ad9 100644 --- a/include/deal.II/fe/fe_nedelec.h +++ b/include/deal.II/fe/fe_nedelec.h @@ -644,9 +644,8 @@ private: /** * Initialize the permutation pattern and the pattern of sign change. * - * @note This function is not fully filled with the correct implementation - * yet. It needs to be consistently implemented in a future release to work - * on meshes that contain cells with flipped faces. + * @note Currently this function is implemented for the finite elements of + * the order k < 5. */ void initialize_quad_dof_index_permutation_and_sign_change(); diff --git a/source/fe/fe_nedelec.cc b/source/fe/fe_nedelec.cc index cc3e1ac90e..ef7f533b19 100644 --- a/source/fe/fe_nedelec.cc +++ b/source/fe/fe_nedelec.cc @@ -215,11 +215,361 @@ template void FE_Nedelec::initialize_quad_dof_index_permutation_and_sign_change() { - // for 1d and 2d, do nothing - if (dim < 3) + // The order of the Nedelec elements equals the tensor degree minus one, + // k = n - 1. In the three-dimensional space the Nedelec elements of the + // lowermost order, k = 0, have only 12 line (edge) dofs. The Nedelec + // elements of the higher orders, k > 0, have 3*(k+1)*(k+2)^2 dofs in + // total if dim=3. The dofs in a cell are distributed between lines + // (edges), quads (faces), and the hex (the interior of the cell) as the + // following: + // + // 12*(k+1) line dofs; (k+1) dofs per line. + // 2*6*k*(k+1) quad dofs; 2*k*(k+1) dofs per quad. + // 3*(k+2)^2*(k+1) hex dofs; + // + // The dofs are indexed in the following order: first all line dofs, + // then all quad dofs, and then all hex dofs. + // + // The line dofs need only sign adjustments. No permutation of line + // dofs is needed. The line dofs are treated by + // internal::FE_PolyTensor::get_dof_sign_change_nedelec(...) + // in fe_poly_tensor.cc. + // + // The hex dofs need no adjustments: they are not shared between + // neighbouring mesh cells. + // + // The two-dimensional Nedelec finite elements share no quad dofs between + // neighbouring mesh cells. The zero-order three-dimensional Nedelec + // finite elements have no quad dofs. Consequently, here we treat only + // quad dofs of the three-dimensional Nedelec finite elements of the + // higher orders, k>0. The questions how the curl looks like in the + // higher-dimensional spaces and what does it mean to be curl-conforming + // if dim>3 we leave unanswered. + // + // In this function we need to change some entries in the following two + // vectors of tables: + // adjust_quad_dof_index_for_face_orientation_table + // and + // adjust_quad_dof_sign_for_face_orientation_table. + // These tables specify the permutations and sign adjustments of the quad + // dofs only. The tables are already filled with zeros meaning no + // permutations or sign change are required. We need to change some + // entries of the tables such that the shape functions that correspond to + // the quad dofs and are shared between neighbouring cells have consistent + // orientations. + // + // The swap tables below describe the dof permutations and sign changes + // that need to be done. The function + // FE_Nedelec::initialize_quad_dof_index_permutation_and_sign_change() + // simply reads the information in the swap tables below and puts it into + // tables + // adjust_quad_dof_index_for_face_orientation_table + // and + // adjust_quad_dof_sign_for_face_orientation_table. + // A good question is: why don't we put the information into the tables of + // deal.II right away? The answer is the following. The information on the + // necessary dof permutations and sign changes is derived by plotting the + // shape functions and observing them on faces of different orientations. + // It is convenient to put the observations first in the format of the + // swap tables below and then convert the swap tables into the format used + // by deal.II. + // + // The dofs on a quad are indexed as the following: + // + // | x0, x1, x2, x3, ..., xk | y0, y1, y2, y3 ..., yk | + // | | | + // |<------ k*(k+1) --------->|<------ k*(k+1) -------->| + // | | + // |<------------------- 2*k*(k+1) -------------------->| + // + // Only one type of dof permutation is needed: swap between two dofs; one + // dof being xi, another yj. That is, if x4 is replaced with y7, + // then y7 must be replaced with x4. Such swaps can be ordered as + // illustrated by the following example: + // + // * + // y0, y9, y1, y2, ..., yk + // --------------------------- (swap) + // x0, x1, x2, x3, ..., xk + // * + // + // An x-dof below the line is swapped with the corresponding y-dof above + // the line. A dof marked by the asterisk must change its sign before the + // swap. + // + // The x-dofs are assumed to have the normal order. There is no need to + // encode it. Therefore, the swap tables need to encode the following + // information: indices of the y-dofs, the sign change of the x-dofs, and + // sign change of the y-dofs. The swap above is encoded as the following: + // + // swap = { 0, 9, 1, 2, ...., yk, // indices of the y-dofs + // 1, 0, 0, 0, ...., 0, // sign change of the x-dofs, + // 0, 1, 0, 0, ...., 0}; // sign change of the y-dofs. + // + // If no swap is needed, -1 is placed instead of the y-dof index. + // + // Such swaps are assembled into the swap table: + // + // swap_table = {swap_0, swap_1, ... swap_7}; + // + // Each swap table contains eight swaps - one swap for each possible quad + // orientation. The deal.II encodes the orientation of a quad using + // three boolean parameters: + // face_orientation - true if face is in standard orientation + // and false otherwise; + // face_rotation - rotation by 90 deg counterclockwise if true; + // face_flip - rotation by 180 deg counterclockwise if true. + // See the documentation of GeometryInfo. + // + // The combined face orientation is computes as + // orientation_no = face_flip*4 + face_rotation*2 + face_orientation*1; + // See tria_orientation.h. + // + // The parameter orientation_no (0...7) indexes the swaps in a swap table. + // + // Nedelec elements of order k have their own swap table, swap_table_k. + // Recall, the swap_table_0 is empty as the Nedelec finite elements of the + // lowermost order have no quad dofs. + + static const int c_swap_table_0 = 0; + + static const int c_swap_table_1[8][3][2] = { // 0 1 + {{0, 1}, // 0 + {0, 0}, + {0, 0}}, + {{-1, -1}, // 1 + {0, 0}, + {0, 0}}, + {{-1, -1}, // 2 + {0, 0}, + {1, 0}}, + {{0, 1}, // 3 + {1, 0}, + {0, 0}}, + {{0, 1}, // 4 + {1, 0}, + {1, 0}}, + {{-1, -1}, // 5 + {1, 0}, + {1, 0}}, + {{-1, -1}, // 6 + {1, 0}, + {0, 0}}, + {{0, 1}, // 7 + {0, 0}, + {1, 0}}}; + + static const int c_swap_table_2[8][3][6] = {// 0 1 2 3 4 5 + {{0, 3, 1, 4, 2, 5}, // 0 + {0, 0, 0, 0, 0, 0}, + {0, 0, 0, 0, 0, 0}}, + {{-1, -1, -1, -1, -1, -1}, // 1 + {0, 0, 0, 0, 0, 0}, + {0, 0, 0, 0, 0, 0}}, + {{-1, -1, -1, -1, -1, -1}, // 2 + {0, 1, 0, 1, 0, 1}, + {1, 0, 1, 1, 0, 1}}, + {{0, 3, 1, 4, 2, 5}, // 3 + {1, 1, 0, 0, 1, 1}, + {0, 0, 0, 1, 1, 1}}, + {{0, 3, 1, 4, 2, 5}, // 4 + {1, 0, 0, 1, 1, 0}, + {1, 0, 1, 0, 1, 0}}, + {{-1, -1, -1, -1, -1, -1}, // 5 + {1, 0, 0, 1, 1, 0}, + {1, 0, 1, 0, 1, 0}}, + {{-1, -1, -1, -1, -1, -1}, // 6 + {1, 1, 0, 0, 1, 1}, + {0, 0, 0, 1, 1, 1}}, + {{0, 3, 1, 4, 2, 5}, // 7 + {0, 1, 0, 1, 0, 1}, + {1, 0, 1, 1, 0, 1}}}; + + static const int c_swap_table_3[8][3][12] = { + // 0 1 2 3 4 5 6 7 8 9 10 11 + {{0, 4, 8, 1, 5, 9, 2, 6, 10, 3, 7, 11}, // 0 + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 1 + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 2 + {0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0}, + {1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0}}, + {{0, 4, 8, 1, 5, 9, 2, 6, 10, 3, 7, 11}, // 3 + {1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0}, + {0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0}}, + {{0, 4, 8, 1, 5, 9, 2, 6, 10, 3, 7, 11}, // 4 + {1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0}, + {1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 5 + {1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0}, + {1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 6 + {1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0}, + {0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0}}, + {{0, 4, 8, 1, 5, 9, 2, 6, 10, 3, 7, 11}, // 7 + {0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0}, + {1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0}}}; + + static const int c_swap_table_4[8][3][20] = { + // Swap sign_X and sign_Y rows if k=4. Why?... + // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 + // 19 + {{0, 5, 10, 15, 1, 6, 11, 16, 2, 7, + 12, 17, 3, 8, 13, 18, 4, 9, 14, 19}, // 0 + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, + -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 1 + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, + {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, + -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 2 + {1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1}, + {0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1}}, + {{0, 5, 10, 15, 1, 6, 11, 16, 2, 7, + 12, 17, 3, 8, 13, 18, 4, 9, 14, 19}, // 3 + {0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1}, + {1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1}}, + {{0, 5, 10, 15, 1, 6, 11, 16, 2, 7, + 12, 17, 3, 8, 13, 18, 4, 9, 14, 19}, // 4 + {1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0}, + {1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, + -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 5 + {1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0}, + {1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 1, 0}}, + {{-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, + -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}, // 6 + {0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1}, + {1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1}}, + {{0, 5, 10, 15, 1, 6, 11, 16, 2, 7, + 12, 17, 3, 8, 13, 18, 4, 9, 14, 19}, // 7 + {1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1}, + {0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1}}}; + + static const int *swap_table_array[5] = {&c_swap_table_0, + &c_swap_table_1[0][0][0], + &c_swap_table_2[0][0][0], + &c_swap_table_3[0][0][0], + &c_swap_table_4[0][0][0]}; + + static const int row_length[5] = {0, 2, 6, 12, 20}; + static const int table_size[5] = { + 0, 8 * 3 * 2, 8 * 3 * 6, 8 * 3 * 12, 8 * 3 * 20}; + + // Only three-dimensional Nedelec finite elements are treated. The + // two-dimensional Nedelec finite elements only need sign adjustments of the + // line dofs. These adjustments are done by + // internal::FE_PolyTensor::get_dof_sign_change_nedelec(...) + // in fe_poly_tensor.cc. The notions of curl and curl-conforming finite + // elements in higher-dimensional spaces, dim >3, are somewhat unclear as + // curl, strictly peaking, exists only in the three-dimensional space. + if (dim != 3) return; - // TODO: Implement this for this class + const unsigned int k = this->tensor_degree() - 1; + + // The Nedelec finite elements of the lowermost order have no quad dofs. + if (k == 0) + return; + + // The finite element orders > 4 are not implemented. + AssertThrow(k < 5, ExcNotImplemented()); + + // TODO: the implementation makes the assumption that all quads have the + // same number of dofs + AssertDimension(this->n_unique_faces(), 1); + const unsigned int face_no = 0; + + Assert( + this->adjust_quad_dof_index_for_face_orientation_table[0].n_elements() == + this->reference_cell().n_face_orientations(face_no) * + this->n_dofs_per_quad(face_no), + ExcInternalError()); + + Assert( + this->adjust_quad_dof_sign_for_face_orientation_table[0].n_elements() == + this->reference_cell().n_face_orientations(face_no) * + this->n_dofs_per_quad(face_no), + ExcInternalError()); + + // The 3D Nedelec finite elements have 2*k*(k+1) dofs per each quad. + Assert(2 * k * (k + 1) == this->n_dofs_per_quad(face_no), ExcInternalError()); + + const int *swap_table = swap_table_array[k]; + + const unsigned int half_dofs = k * (k + 1); // see below; + + const int rl = row_length[k]; + + int offset = table_size[0]; + // The assignment above is to prevent the compiler from complaining about the + // unused table_size. + + int value = 0; + + for (const bool face_orientation : {false, true}) + for (const bool face_rotation : {false, true}) + for (const bool face_flip : {false, true}) + { + const auto case_no = + internal::combined_face_orientation(face_orientation, + face_rotation, + face_flip); + + // The dofs on a quad are indexed as the following: + // + // | x0, x1, x2, x3, ..., xk | y0, y1, y2, y3 ..., yk | + // | | | + // |-- half_ dofs = k*(k+1) --|-- half_dofs = k*(k+1) --| + // | | + // |-------------------- 2*k*(k+1) ---------------------| + + for (unsigned int indx_x = 0; indx_x < half_dofs; indx_x++) + { + offset = 3 * rl * case_no + 0 * rl + indx_x; + Assert(offset < table_size[k], ExcInternalError()); + value = *(swap_table + offset); + Assert(value < table_size[k], ExcInternalError()); + Assert(value > -2, ExcInternalError()); + + if (value != -1) + { + const unsigned int indx_y = + half_dofs + static_cast(value); + + // dofs swap + this + ->adjust_quad_dof_index_for_face_orientation_table[face_no]( + indx_x, case_no) = indx_y - indx_x; + + this + ->adjust_quad_dof_index_for_face_orientation_table[face_no]( + indx_y, case_no) = indx_x - indx_y; + } + + // dof sign change + offset = 3 * rl * case_no + 1 * rl + indx_x; + Assert(offset < table_size[k], ExcInternalError()); + value = *(swap_table + offset); + Assert((value == 0) || (value == 1), ExcInternalError()); + + this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + indx_x, case_no) = static_cast(value); + + + offset = 3 * rl * case_no + 2 * rl + indx_x; + Assert(offset < table_size[k], ExcInternalError()); + value = *(swap_table + offset); + Assert((value == 0) || (value == 1), ExcInternalError()); + + this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + indx_x + half_dofs, case_no) = static_cast(value); + } + } + return; } diff --git a/source/fe/fe_poly_tensor.cc b/source/fe/fe_poly_tensor.cc index 5f27d831a0..9fda6d6a9a 100644 --- a/source/fe/fe_poly_tensor.cc +++ b/source/fe/fe_poly_tensor.cc @@ -134,42 +134,141 @@ namespace internal // nothing to do in 1d } - - template void get_dof_sign_change_nedelec( const typename dealii::Triangulation<2, spacedim>::cell_iterator &cell, - const FiniteElement<2, spacedim> & /*fe*/, + const FiniteElement<2, spacedim> &fe, const std::vector &mapping_kind, std::vector &line_dof_sign) { - const unsigned int dim = 2; - // TODO: This fixes only lowest order - for (unsigned int l = 0; l < GeometryInfo::lines_per_cell; ++l) - if ((cell->line_orientation(l) == - numbers::reverse_line_orientation) && + // The Nedelec finite elements in two spacial dimensions have two types + // of dofs: the line dofs and the quad dofs. The line dofs are + // associated with the edges. They are shared between the neighbouring + // cells and, consequently, need to be adjusted to compensate for a + // possible mismatch in the edge orientation. The quad dofs, they are + // associated with the interiors of the cells, are not shared between + // the cells in two spacial dimensions and need no adjustments. + // + // The Nedelec finite elements in two spacial dimensions have + // 2*(k+1)*(k+2) dofs per cell. All dofs are distributed between the + // line and quad dofs as the following: + // + // 4*(k+1) line dofs; (k+1) dofs per line. + // 2*k*(k+1) quad dofs. + // + // The dofs are indexed in the following order: first all line dofs, + // then all quad dofs. + // + // Here we adjust the line dofs. The sign of a line dof needs to be + // changed if the edge on which the dof resides points in the opposite + // direction. + const unsigned int k = fe.tensor_degree() - 1; + + for (unsigned int l = 0; l < GeometryInfo<2>::lines_per_cell; ++l) + if (!(cell->line_orientation(l)) && mapping_kind[0] == mapping_nedelec) - line_dof_sign[l] = -1.0; + { + if (k == 0) + { + // The lowest order element (k=0) is straightforward, because + // there is a single dof per edge, which needs to be flipped: + line_dof_sign[l] = -1.0; + } + else + { + // The case k > 0 is a bit more complicated. As we adjust + // only edge dofs in this function, we need to concern + // ourselves with the first 4*(k+1) entries in line_dof_sign + // vector ignoring the rest. There are (k+1) dofs per edge. + // Let us consider the local dof indices on one edge, + // local_line_dof = 0...k. The shape functions with even + // indices are asymmetric. The corresponding dofs need sign + // adjustment if the edge points in the opposite direction. + // The shape functions with odd indices are symmetric. The + // corresponding dofs need no sign adjustment even if the edge + // points in the opposite direction. In the current context + // the notion of symmetry of a shape function means that a + // shape function looks exactly the same if it is looked upon + // from the centers of the two neighbouring cells that share + // it. + for (unsigned int local_line_dof = 0; + local_line_dof < (k + 1); + local_line_dof++) + if (local_line_dof % 2 == 0) + line_dof_sign[local_line_dof + l * (k + 1)] = -1.0; + } + } } - template void get_dof_sign_change_nedelec( const typename dealii::Triangulation<3, spacedim>::cell_iterator &cell, - const FiniteElement<3, spacedim> & /*fe*/, + const FiniteElement<3, spacedim> &fe, const std::vector &mapping_kind, std::vector &line_dof_sign) { - const unsigned int dim = 3; - // TODO: This is probably only going to work for those elements for - // which all dofs are face dofs - for (unsigned int l = 0; l < GeometryInfo::lines_per_cell; ++l) - if ((cell->line_orientation(l) == - numbers::reverse_line_orientation) && + // This function does half of the job - it adjusts the sign of the + // line (edge) dofs. In the three-dimensional space the quad (face) dofs + // need to be adjusted as well. The quad dofs are treated by + // FE_Nedelec::initialize_quad_dof_index_permutation_and_sign_change() + // in fe_nedelec.cc. The dofs associated with the interior of the cells, + // the hex dofs, need no adjustments as they are not shared between the + // neighboring cells. + + const unsigned int k = fe.tensor_degree() - 1; + // The order of the Nedelec elements equals the tensor degree minus one, + // k = n - 1. In the three-dimensional space the Nedelec elements of the + // lowermost order, k = 0, have only 12 line (edge) dofs. The Nedelec + // elements of the higher orders, k > 0, have 3*(k+1)*(k+2)^2 dofs in + // total if dim=3. The dofs in a cell are distributed between lines + // (edges), quads (faces), and the hex (the interior of the cell) as the + // following: + // + // 12*(k+1) line dofs; (k+1) dofs per line. + // 2*6*k*(k+1) quad dofs; 2*k*(k+1) dofs per quad. + // 3*(k+2)^2*(k+1) hex dofs. + // + // The dofs are indexed in the following order: first all line dofs, + // then all quad dofs, and then all hex dofs. + // + // Here we adjust only the line (edge) dofs. The line dofs need only + // sign adjustment. That is, no permutation of the line dofs is needed. + for (unsigned int l = 0; l < GeometryInfo<3>::lines_per_cell; ++l) + if (!(cell->line_orientation(l)) && mapping_kind[0] == mapping_nedelec) - line_dof_sign[l] = -1.0; + { + if (k == 0) + { + // The lowest order element (k=0) is straightforward, because + // there is a single dof per edge, which needs to be flipped: + line_dof_sign[l] = -1.0; + } + else + { + // The case k > 0 is a bit more complicated. As we adjust + // only edge dofs in this function, we need to concern + // ourselves with the first 12*(k+1) entries in line_dof_sign + // vector ignoring the rest. There are (k+1) dofs per edge. + // Let us consider the local dof indices on one edge, + // local_line_dof = 0...k. The shape functions with even + // indices are asymmetric. The corresponding dofs need sign + // adjustment if the edge points in the opposite direction. + // The shape functions with odd indices are symmetric. The + // corresponding dofs need no sign adjustment even if the edge + // points in the opposite direction. In the current context + // the notion of symmetry of a shape function means that a + // shape function looks exactly the same if it is looked upon + // from the centers of the two neighbouring cells that share + // it. + for (unsigned int local_line_dof = 0; + local_line_dof < (k + 1); + local_line_dof++) + if (local_line_dof % 2 == 0) + line_dof_sign[local_line_dof + l * (k + 1)] = -1.0; + } + } } } // namespace } // namespace FE_PolyTensor diff --git a/tests/fe/fe_nedelec_face_edge_orientation.h b/tests/fe/fe_nedelec_face_edge_orientation.h new file mode 100644 index 0000000000..402650c738 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation.h @@ -0,0 +1,451 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// In general, vertices of a face can be incorporated into two neighbouring +// cells such that the face appears to have different orientations in the two +// cells. The code of deal.II compensates for such mismatch in face orientation. +// This code tests if deal.II does the compensation well for the Nedelec finite +// elements in two- and three-dimensional spaces on globally refined meshes. +// +// In the two-dimensional space the Nedelec finite elements have dofs associated +// with lines (edges) and the quad (the interior of the cell). Only line +// dofs are shared between the cells in 2D. In the three-dimensional space the +// Nedelec finite elements have dofs that are associated with lines (edges), +// quads (faces), and the hex (the interior of the cell). The dofs associated +// with lines and quads are shared between the cells in 3D. Therefore, the dofs +// affected by the mismatch in face orientation are: +// (i) the line dofs in 2D, +// (ii) the line dofs in 3D, +// (iii) the quad dofs in 3D. +// +// The shape functions that correspond to these three kinds of dofs share the +// same feature: they exhibit reflection symmetry in the shared face if the +// two cells that share the face are cubes (squares) of equal volume (area). +// The code below exploits this feature. +// +// First, a mesh with a particular face mismatch is generated. Second, all shape +// functions are sampled at sampling points. Indices of the sampling points in +// one cell is the mirrored image of the indices in the other cell. Third, the +// sampled shape functions are sorted in two categories: (i) shared between two +// cells and (ii) locked inside a cell. The locked shaped functions are ignored. +// The reflection symmetry of the shared shape functions is tested by computing +// a score coefficient which is a cumulative l2_norm of a difference between a +// sampled vector and its reflection. The score coefficient should be zero. The +// program logs the score coefficients and two versions of the total amount of +// shared dofs - one is theoretically expected, another is actually counted +// during the execution of the code. +// +// In the three-dimensional space deal.II encodes the mismatch in face +// orientation by three boolean parameters: +// +// (i) face_orientation - true if face is in standard orientation. +// (ii) face_rotation - rotation by 90 deg counterclockwise if true; +// (iii) face_flip - rotation by 180 deg counterclockwise if true. +// +// See the documentation of GeometryInfo. +// +// The combined face orientation is computes as +// +// orientation_no = face_flip*4 + face_rotation*2 + face_orientation*1; +// +// See tria_orientation.h. +// +// In the two-dimensional space the face orientation is encoded by +// one boolean parameter line_orientation. The parameter equals true if the +// faces have the same orientation. +// +// In this program the one parameter is used to encode the face orientation in +// two- and three- dimensions. This parameter is +// +// unsigned char combined_face_orientation; +// +// Interpretation of this parameter in three dimensions is straightforward: +// it is the combined face orientation as it is known in deal.II. In the +// two-dimensional space this parameter is interpreted as the following. +// The odd parameter (LSB = 1) implies that the shared faces (lines) are +// aligned. The even parameter (LSB = 0) implies that the shared faces point in +// opposite directions. That is, the least significant bit (LSB) of the +// combined_face_orientation parameter is interpreted as the +// line_orientation parameter which is used in deal.II on two-dimensional +// meshes. +// +// Test meshes in two and three dimensions are generated by functions +// GridGenerator::non_standard_orientation_mesh(). +// The three-dimensional mesh is used as it is. Only the leftmost shared face +// is used in two dimensions. That is, shape functions are sampled only on +// the leftmost and the middle cells of the two-dimensional mesh. + +#include + +#include + +#include +#include + +#include +#include +#include +#include + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + +template +void +create_triangulation(Triangulation &triangulation, + const unsigned char combined_face_orientation); + +template +void +create_sampling_points(std::vector> &sampling_points_left, + std::vector> &sampling_points_right, + unsigned int N); + +template +void +save_sampling_points(std::vector> &sampling_points_left, + std::vector> &sampling_points_right); + +template +void +print_shared_face(Triangulation &triangulation); + +template <> +void +create_triangulation<2>(Triangulation<2> &triangulation, + const unsigned char combined_face_orientation) + +{ + // Only the leftmost shared face is considered. Even values of + // combined_face_orientation produce misaligned faces (the middle square is + // rotated by 2*[pi/2]). Odd values of combined_face_orientation produce + // aligned faces (the middle square is not rotated). + const unsigned int n_rotate_middle_square = + (combined_face_orientation % 2 == 0) ? 2 : 0; + + GridGenerator::non_standard_orientation_mesh(triangulation, + n_rotate_middle_square); +} + +template <> +void +create_triangulation<3>(Triangulation<3> &triangulation, + const unsigned char combined_face_orientation) + +{ + bool face_orientation; + bool face_rotation; + bool face_flip; + + std::tie(face_orientation, face_rotation, face_flip) = + dealii::internal::split_face_orientation(combined_face_orientation); + + GridGenerator::non_standard_orientation_mesh( + triangulation, face_orientation, face_flip, face_rotation, false); +} + +template <> +void +create_sampling_points<2>(std::vector> &sampling_points_left, + std::vector> &sampling_points_right, + unsigned int N) +{ + const double offset = 0.1; + const double step = (1.0 - 2 * offset) / (N - 1); + + unsigned int indx = 0; + + for (unsigned int j = 0; j < N; j++) + for (unsigned int i = 0; i < N; i++) + { + indx = i + N * j; + + sampling_points_left.at(indx)[0] = offset + step * i; + sampling_points_left.at(indx)[1] = offset + step * j; + + sampling_points_right.at(indx)[0] = + 2.0 - sampling_points_left.at(indx)[0]; + sampling_points_right.at(indx)[1] = sampling_points_left.at(indx)[1]; + } +} + +template <> +void +create_sampling_points<3>(std::vector> &sampling_points_left, + std::vector> &sampling_points_right, + unsigned int N) +{ + const double offset = 0.1; + const double step = (1.0 - 2 * offset) / (N - 1); + + unsigned int indx = 0; + + for (unsigned int k = 0; k < N; k++) + for (unsigned int j = 0; j < N; j++) + for (unsigned int i = 0; i < N; i++) + { + indx = i + N * j + N * N * k; + + sampling_points_left.at(indx)[0] = offset + step * i; + sampling_points_left.at(indx)[1] = offset + step * j; + sampling_points_left.at(indx)[2] = offset + step * k; + + sampling_points_right.at(indx)[0] = + 2.0 - sampling_points_left.at(indx)[0]; + sampling_points_right.at(indx)[1] = sampling_points_left.at(indx)[1]; + sampling_points_right.at(indx)[2] = sampling_points_left.at(indx)[2]; + } +} + +template +void +save_mesh(Triangulation &triangulation) +{ + GridOut gridout; + GridOutFlags::Msh msh_flags(true, true); + gridout.set_flags(msh_flags); + + std::ofstream ofs("mesh_" + std::to_string(dim) + "D.msh"); + gridout.write_msh(triangulation, ofs); +} + +template <> +void +save_sampling_points<2>(std::vector> &sampling_points_left, + std::vector> &sampling_points_right) +{ + std::ofstream ofs_left("sampling_points_left_2D.csv"); + std::ofstream ofs_right("sampling_points_right_2D.csv"); + + for (unsigned int i = 0; i < sampling_points_left.size(); i++) + { + ofs_left << sampling_points_left.at(i)[0] << " , " + << sampling_points_left.at(i)[1] << std::endl; + + ofs_right << sampling_points_right.at(i)[0] << " , " + << sampling_points_right.at(i)[1] << std::endl; + } + + ofs_left.close(); + ofs_right.close(); +} + +template <> +void +save_sampling_points<3>(std::vector> &sampling_points_left, + std::vector> &sampling_points_right) +{ + std::ofstream ofs_left("sampling_points_left_3D.csv"); + std::ofstream ofs_right("sampling_points_right_3D.csv"); + + for (unsigned int i = 0; i < sampling_points_left.size(); i++) + { + ofs_left << sampling_points_left.at(i)[0] << " , " + << sampling_points_left.at(i)[1] << " , " + << sampling_points_left.at(i)[2] << std::endl; + + ofs_right << sampling_points_right.at(i)[0] << " , " + << sampling_points_right.at(i)[1] << " , " + << sampling_points_right.at(i)[2] << std::endl; + } + + ofs_left.close(); + ofs_right.close(); +} + +template <> +void +print_shared_face(Triangulation<2> &triangulation) +{ + const double eps = 1e-12; + + for (const auto &cell : triangulation.active_cell_iterators()) + for (unsigned int f = 0; f < 4; f++) + if (cell->neighbor_index(f) != -1) + if (((cell->face(f)->center()(0) - 1.0) < eps) && + ((cell->face(f)->center()(1) - 0.5) < eps)) + { + std::cout << "Cell center = (" << cell->center() << ")" + << std::endl; + std::cout << "The local index of the shared face: " << f + << std::endl; + std::cout << "Shared face - orientation: " + << cell->line_orientation(f) << std::endl; + } +} + +template <> +void +print_shared_face<3>(Triangulation<3> &triangulation) +{ + for (const auto &cell : triangulation.active_cell_iterators()) + { + unsigned int face = -1; + + for (unsigned int f = 0; f < 6; f++) + if (cell->neighbor_index(f) != -1) + face = f; + + std::cout << "Cell id = " << cell->id() << ": " << std::endl; + + std::cout << "The local index of the shared face: " << face << std::endl; + + std::cout << "Shared face - flip: " << cell->face_flip(face) << std::endl; + std::cout << "Shared face - rotation: " << cell->face_rotation(face) + << std::endl; + std::cout << "Shared face - orientation: " << cell->face_orientation(face) + << std::endl; + std::cout + << "Shared face - combined orientation: " + << static_cast( + internal::combined_face_orientation(cell->face_orientation(face), + cell->face_rotation(face), + cell->face_flip(face))) + << std::endl; + } +} + +struct SFData +{ + bool shared; + bool error; + double cum_norm_left; + double cum_norm_right; + double cum_norm_diff; +}; + +void +classify(const std::vector> &shape_function_left, + const std::vector> &shape_function_right, + SFData &result) +{ + result.shared = false; + result.error = false; + result.cum_norm_left = 0.0; + result.cum_norm_right = 0.0; + result.cum_norm_diff = 0.0; + + double th = 1e-3; + double eps = 1e-12; + + unsigned int v_size = shape_function_left.size(); + for (unsigned int i = 0; i < v_size; i++) + { + result.cum_norm_left += shape_function_left.at(i).l2_norm(); + result.cum_norm_right += shape_function_right.at(i).l2_norm(); + + Vector diff(shape_function_left.at(i)); + + diff.add(-1.0, shape_function_right.at(i)); + + result.cum_norm_diff += diff.l2_norm(); + } + + double s_v_size = std::sqrt(v_size); + + result.cum_norm_left /= s_v_size; + result.cum_norm_right /= s_v_size; + result.cum_norm_diff /= s_v_size; + + if ((result.cum_norm_left > th) && (result.cum_norm_right > th)) + result.shared = true; + + if (result.shared) + if (result.cum_norm_diff > eps) + result.error = true; +} + +template +void +run(unsigned char combined_face_orientation) +{ + Triangulation triangulation; + + create_triangulation(triangulation, combined_face_orientation); + + // save_mesh(triangulation); + + // print_shared_face(triangulation); + + FE_Nedelec fe(order); + + DoFHandler dof_handler(triangulation); + + dof_handler.distribute_dofs(fe); + + unsigned int N = order + 2; + unsigned int M = static_cast(std::pow(N, dim)); + + std::vector> sampling_points_left(M, Point()); + std::vector> sampling_points_right(M, Point()); + + std::vector> shape_function_left(M, Vector(dim)); + std::vector> shape_function_right(M, Vector(dim)); + + create_sampling_points(sampling_points_left, sampling_points_right, N); + + // save_sampling_points(sampling_points_left, sampling_points_right); + + Vector solution(dof_handler.n_dofs()); + + dealii::Functions::FEFieldFunction field_function(dof_handler, solution); + + unsigned int n_shared_dofs_theory; + + if (dim == 2) + n_shared_dofs_theory = order + 1; + + if (dim == 3) + n_shared_dofs_theory = (2 * order + 1) * (order + 1); + + unsigned int Q = solution.size(); + unsigned int n_shared_dofs = 0; + for (unsigned int i = 0; i < Q; i++) + { + solution.reinit(Q); + solution(i) = 1.0; + + field_function.vector_value_list(sampling_points_left, + shape_function_left); + field_function.vector_value_list(sampling_points_right, + shape_function_right); + + SFData result; + classify(shape_function_left, shape_function_right, result); + + if (result.error) + deallog << "error: " << result.cum_norm_left << " " + << result.cum_norm_right << " " << result.cum_norm_diff + << std::endl; + + if ((result.shared) && (!result.error)) + { + n_shared_dofs++; + deallog << result.cum_norm_diff << std::endl; + } + + + if (n_shared_dofs == n_shared_dofs_theory) + break; + } + + deallog << n_shared_dofs << " : " << n_shared_dofs_theory << std::endl; +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_0.cc b/tests/fe/fe_nedelec_face_edge_orientation_0.cc new file mode 100644 index 0000000000..acbc941177 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_0.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 0. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 0; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_0.output b/tests/fe/fe_nedelec_face_edge_orientation_0.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_0.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_1.cc b/tests/fe/fe_nedelec_face_edge_orientation_1.cc new file mode 100644 index 0000000000..64b6be3b84 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_1.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 1. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 1; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_1.output b/tests/fe/fe_nedelec_face_edge_orientation_1.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_1.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_2.cc b/tests/fe/fe_nedelec_face_edge_orientation_2.cc new file mode 100644 index 0000000000..445ec81a02 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_2.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 2. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 2; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_2.output b/tests/fe/fe_nedelec_face_edge_orientation_2.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_2.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_3.cc b/tests/fe/fe_nedelec_face_edge_orientation_3.cc new file mode 100644 index 0000000000..6097be34a0 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_3.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 3. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 3; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_3.output b/tests/fe/fe_nedelec_face_edge_orientation_3.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_3.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_4.cc b/tests/fe/fe_nedelec_face_edge_orientation_4.cc new file mode 100644 index 0000000000..c2a037e820 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_4.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 4. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 4; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_4.output b/tests/fe/fe_nedelec_face_edge_orientation_4.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_4.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_5.cc b/tests/fe/fe_nedelec_face_edge_orientation_5.cc new file mode 100644 index 0000000000..464a8f2fef --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_5.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 5. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 5; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_5.output b/tests/fe/fe_nedelec_face_edge_orientation_5.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_5.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_6.cc b/tests/fe/fe_nedelec_face_edge_orientation_6.cc new file mode 100644 index 0000000000..1bd08dbfea --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_6.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 6. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 6; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_6.output b/tests/fe/fe_nedelec_face_edge_orientation_6.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_6.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_7.cc b/tests/fe/fe_nedelec_face_edge_orientation_7.cc new file mode 100644 index 0000000000..ff859f1423 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_7.cc @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a two-cells 3D mesh +// with a shared face orientation 7. See the notes in the header file for +// more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 7; + + initlog(); + run<3, 0>(combined_face_orientation); + run<3, 1>(combined_face_orientation); + run<3, 2>(combined_face_orientation); + run<3, 3>(combined_face_orientation); + run<3, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_7.output b/tests/fe/fe_nedelec_face_edge_orientation_7.output new file mode 100644 index 0000000000..a4a8ecb81b --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_7.output @@ -0,0 +1,101 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::6 : 6 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::15 : 15 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::28 : 28 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::45 : 45 diff --git a/tests/fe/fe_nedelec_face_edge_orientation_8.cc b/tests/fe/fe_nedelec_face_edge_orientation_8.cc new file mode 100644 index 0000000000..647c351adf --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_8.cc @@ -0,0 +1,41 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2013 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Tests FE_Nedelec finite elements of orders 0...4 on a simple 2D mesh +// with a shared face orientations 0 and 1 (the are only two possible face +// orientations in 2D). See the notes in the header file for more detail. + +#include "fe_nedelec_face_edge_orientation.h" + +int +main() +{ + unsigned char combined_face_orientation = 0; + + initlog(); + + run<2, 0>(combined_face_orientation); + run<2, 1>(combined_face_orientation); + run<2, 2>(combined_face_orientation); + run<2, 3>(combined_face_orientation); + run<2, 4>(combined_face_orientation); + + combined_face_orientation = 1; + + run<2, 0>(combined_face_orientation); + run<2, 1>(combined_face_orientation); + run<2, 2>(combined_face_orientation); + run<2, 3>(combined_face_orientation); + run<2, 4>(combined_face_orientation); +} diff --git a/tests/fe/fe_nedelec_face_edge_orientation_8.output b/tests/fe/fe_nedelec_face_edge_orientation_8.output new file mode 100644 index 0000000000..5018a3b997 --- /dev/null +++ b/tests/fe/fe_nedelec_face_edge_orientation_8.output @@ -0,0 +1,41 @@ + +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::2 : 2 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::3 : 3 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::4 : 4 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::5 : 5 +DEAL::0 +DEAL::1 : 1 +DEAL::0 +DEAL::0 +DEAL::2 : 2 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::3 : 3 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::4 : 4 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::0 +DEAL::5 : 5 -- 2.39.5