From c0f4e3832fb1fc6658c7044cec960745cfecb728 Mon Sep 17 00:00:00 2001 From: Sebastian Kinnewig Date: Tue, 21 Sep 2021 13:12:32 +0200 Subject: [PATCH] Add hessians to FE_NedelecSZ --- doc/news/changes/minor/20220128Kinnewig | 3 + include/deal.II/fe/fe_nedelec_sz.h | 10 + source/fe/fe_nedelec_sz.cc | 875 +++++++++++++++++- ...e_nedelec_sz_hessian_divergence_theorem.cc | 172 ++++ ...delec_sz_hessian_divergence_theorem.output | 760 +++++++++++++++ 5 files changed, 1792 insertions(+), 28 deletions(-) create mode 100644 doc/news/changes/minor/20220128Kinnewig create mode 100644 tests/fe/fe_nedelec_sz_hessian_divergence_theorem.cc create mode 100644 tests/fe/fe_nedelec_sz_hessian_divergence_theorem.output diff --git a/doc/news/changes/minor/20220128Kinnewig b/doc/news/changes/minor/20220128Kinnewig new file mode 100644 index 0000000000..941bc53776 --- /dev/null +++ b/doc/news/changes/minor/20220128Kinnewig @@ -0,0 +1,3 @@ +New: Added support for hessians to FE_NedelecSZ +
+(Sebastian Kinnewig, 2022/01/28) diff --git a/include/deal.II/fe/fe_nedelec_sz.h b/include/deal.II/fe/fe_nedelec_sz.h index 50801f6366..11d93ed495 100644 --- a/include/deal.II/fe/fe_nedelec_sz.h +++ b/include/deal.II/fe/fe_nedelec_sz.h @@ -255,6 +255,16 @@ protected: */ mutable std::vector>> shape_grads; + /** + * Storage for shape function hessians on the reference element. We only + * pre-compute cell-based DoFs, as the edge- and face-based DoFs depend on + * the cell. + * + * Due to the cell-dependent DoFs, this variable is declared mutable. + */ + mutable std::vector>> + shape_hessians; + /** * Storage for all possible edge parameterization between vertices. These * are required in the computation of edge- and face-based DoFs, which are diff --git a/source/fe/fe_nedelec_sz.cc b/source/fe/fe_nedelec_sz.cc index 3e80a41035..2163c9c095 100644 --- a/source/fe/fe_nedelec_sz.cc +++ b/source/fe/fe_nedelec_sz.cc @@ -177,10 +177,12 @@ FE_NedelecSZ::get_data( std::vector>( n_q_points)); } - // Not implementing second derivatives yet: + if ((flags & update_hessians) != 0u) { - Assert(false, ExcNotImplemented()); + data.shape_hessians.resize(this->n_dofs_per_cell(), + std::vector>( + n_q_points)); } std::vector> p_list(n_q_points); @@ -401,7 +403,8 @@ FE_NedelecSZ::get_data( cell_type2_offset + degree * degree; const unsigned int cell_type3_offset2 = cell_type3_offset1 + degree; - if ((flags & (update_values | update_gradients)) != 0u) + if ((flags & (update_values | update_gradients | update_hessians)) != + 0u) { // compute all points we must evaluate the 1d polynomials at: std::vector> cell_points(n_q_points); @@ -424,12 +427,12 @@ FE_NedelecSZ::get_data( // We only need poly values and 1st derivative for // update_values, but need the 2nd derivative too for - // update_gradients. - // - // Note that this will need to be updated if we're supporting - // update_hessians. + // update_gradients. For update_hessians we also need the 3rd + // derivatives. const unsigned int poly_length( - (flags & update_gradients) != 0u ? 3 : 2); + (flags & update_hessians) != 0u ? + 4 : + ((flags & update_gradients) != 0u ? 3 : 2)); std::vector> polyx( degree, std::vector(poly_length)); @@ -534,6 +537,83 @@ FE_NedelecSZ::get_data( data.shape_grads[dof_index3_2][q][1][1] = 0.0; } } + if (flags & update_hessians) + { + for (unsigned int j = 0; j < degree; ++j) + { + const unsigned int shift_j(j * degree); + for (unsigned int i = 0; i < degree; ++i) + { + const unsigned int shift_ij(i + shift_j); + + // Type 1: + const unsigned int dof_index1(cell_type1_offset + + shift_ij); + data.shape_hessians[dof_index1][q][0][0][0] = + 8.0 * polyx[i][3] * polyy[j][0]; + data.shape_hessians[dof_index1][q][1][0][0] = + 8.0 * polyx[i][2] * polyy[j][1]; + + data.shape_hessians[dof_index1][q][0][1][0] = + data.shape_hessians[dof_index1][q][1][0][0]; + data.shape_hessians[dof_index1][q][1][1][0] = + 8.0 * polyx[i][1] * polyy[j][2]; + + data.shape_hessians[dof_index1][q][0][0][1] = + data.shape_hessians[dof_index1][q][1][0][0]; + data.shape_hessians[dof_index1][q][1][0][1] = + data.shape_hessians[dof_index1][q][1][1][0]; + + data.shape_hessians[dof_index1][q][0][1][1] = + data.shape_hessians[dof_index1][q][1][1][0]; + data.shape_hessians[dof_index1][q][1][1][1] = + 8.0 * polyx[i][0] * polyy[j][3]; + + + + // Type 2: + const unsigned int dof_index2(cell_type2_offset + + shift_ij); + for (unsigned int d = 0; d < dim; ++d) + { + data.shape_hessians[dof_index2][q][0][0][d] = + data.shape_hessians[dof_index1][q][0][0][d]; + data.shape_hessians[dof_index2][q][0][1][d] = + data.shape_hessians[dof_index1][q][0][1][d]; + data.shape_hessians[dof_index2][q][1][0][d] = + -1.0 * + data.shape_hessians[dof_index1][q][1][0][d]; + data.shape_hessians[dof_index2][q][1][1][d] = + -1.0 * + data.shape_hessians[dof_index1][q][1][1][d]; + } + } + // Type 3: + const unsigned int dof_index3_1(cell_type3_offset1 + + j); + data.shape_hessians[dof_index3_1][q][0][0][0] = 0.0; + data.shape_hessians[dof_index3_1][q][0][0][1] = 0.0; + data.shape_hessians[dof_index3_1][q][0][1][0] = 0.0; + data.shape_hessians[dof_index3_1][q][0][1][1] = + 4.0 * polyy[j][2]; + data.shape_hessians[dof_index3_1][q][1][0][0] = 0.0; + data.shape_hessians[dof_index3_1][q][1][0][1] = 0.0; + data.shape_hessians[dof_index3_1][q][1][1][0] = 0.0; + data.shape_hessians[dof_index3_1][q][1][1][1] = 0.0; + + const unsigned int dof_index3_2(cell_type3_offset2 + + j); + data.shape_hessians[dof_index3_2][q][0][0][0] = 0.0; + data.shape_hessians[dof_index3_2][q][0][0][1] = 0.0; + data.shape_hessians[dof_index3_2][q][0][1][0] = 0.0; + data.shape_hessians[dof_index3_2][q][0][1][1] = 0.0; + data.shape_hessians[dof_index3_2][q][1][0][0] = + 4.0 * polyx[j][2]; + data.shape_hessians[dof_index3_2][q][1][0][1] = 0.0; + data.shape_hessians[dof_index3_2][q][1][1][0] = 0.0; + data.shape_hessians[dof_index3_2][q][1][1][1] = 0.0; + } + } } } break; @@ -799,7 +879,8 @@ FE_NedelecSZ::get_data( // for cell-based shape functions: // these don't depend on the cell, so can precompute all here: - if ((flags & (update_values | update_gradients)) != 0u) + if ((flags & + (update_values | update_gradients | update_hessians)) != 0u) { // Cell-based shape functions: // @@ -851,11 +932,15 @@ FE_NedelecSZ::get_data( } } - - // only need poly values and 1st derivative for update_values, - // but need 2nd derivative too for update_gradients. + // We only need poly values and 1st derivative for + // update_values, but need the 2nd derivative too for + // update_gradients. For update_hessians we also need 3rd + // derivative. const unsigned int poly_length( - (flags & update_gradients) != 0u ? 3 : 2); + (flags & update_hessians) != 0u ? + 4 : + ((flags & update_gradients) != 0u ? 3 : 2)); + // Loop through quad points: for (unsigned int q = 0; q < n_q_points; ++q) { @@ -1091,6 +1176,284 @@ FE_NedelecSZ::get_data( } } } + if (flags & update_hessians) + { + for (unsigned int k = 0; k < degree; ++k) + { + const unsigned int shift_k(k * degree * degree); + const unsigned int shift_j( + k * degree); // Used below when subbing k for j + // type 3 + + for (unsigned int j = 0; j < degree; ++j) + { + const unsigned int shift_jk(j * degree + + shift_k); + for (unsigned int i = 0; i < degree; ++i) + { + const unsigned int shift_ijk(shift_jk + + i); + + // Type 1: + const unsigned int dof_index1( + cell_type1_offset + shift_ijk); + + data.shape_hessians[dof_index1][q][0][0] + [0] = + 8.0 * polyx[i][3] * polyy[j][0] * + polyz[k][0]; + data.shape_hessians[dof_index1][q][1][0] + [0] = + 8.0 * polyx[i][2] * polyy[j][1] * + polyz[k][0]; + data.shape_hessians[dof_index1][q][2][0] + [0] = + 8.0 * polyx[i][2] * polyy[j][0] * + polyz[k][1]; + + data.shape_hessians[dof_index1][q][0][1] + [0] = + data.shape_hessians[dof_index1][q][1][0] + [0]; + data.shape_hessians[dof_index1][q][1][1] + [0] = + 8.0 * polyx[i][1] * polyy[j][2] * + polyz[k][0]; + data.shape_hessians[dof_index1][q][2][1] + [0] = + 8.0 * polyx[i][1] * polyy[j][1] * + polyz[k][1]; + + data.shape_hessians[dof_index1][q][0][2] + [0] = + data.shape_hessians[dof_index1][q][2][0] + [0]; + data.shape_hessians[dof_index1][q][1][2] + [0] = + data.shape_hessians[dof_index1][q][2][1] + [0]; + data.shape_hessians[dof_index1][q][2][2] + [0] = + 8.0 * polyx[i][1] * polyy[j][0] * + polyz[k][2]; + + + data.shape_hessians[dof_index1][q][0][0] + [1] = + data.shape_hessians[dof_index1][q][1][0] + [0]; + data.shape_hessians[dof_index1][q][1][0] + [1] = + data.shape_hessians[dof_index1][q][1][1] + [0]; + data.shape_hessians[dof_index1][q][2][0] + [1] = + data.shape_hessians[dof_index1][q][2][1] + [0]; + + data.shape_hessians[dof_index1][q][0][1] + [1] = + data.shape_hessians[dof_index1][q][1][1] + [0]; + data.shape_hessians[dof_index1][q][1][1] + [1] = + 8.0 * polyx[i][0] * polyy[j][3] * + polyz[k][0]; + data.shape_hessians[dof_index1][q][2][1] + [1] = + 8.0 * polyx[i][0] * polyy[j][2] * + polyz[k][1]; + + data.shape_hessians[dof_index1][q][0][2] + [1] = + data.shape_hessians[dof_index1][q][2][1] + [0]; + data.shape_hessians[dof_index1][q][1][2] + [1] = + data.shape_hessians[dof_index1][q][2][1] + [1]; + data.shape_hessians[dof_index1][q][2][2] + [1] = + 8.0 * polyx[i][0] * polyy[j][1] * + polyz[k][2]; + + + data.shape_hessians[dof_index1][q][0][0] + [2] = + data.shape_hessians[dof_index1][q][2][0] + [0]; + data.shape_hessians[dof_index1][q][1][0] + [2] = + data.shape_hessians[dof_index1][q][2][1] + [0]; + data.shape_hessians[dof_index1][q][2][0] + [2] = + data.shape_hessians[dof_index1][q][2][2] + [0]; + + data.shape_hessians[dof_index1][q][0][1] + [2] = + data.shape_hessians[dof_index1][q][2][1] + [0]; + data.shape_hessians[dof_index1][q][1][1] + [2] = + data.shape_hessians[dof_index1][q][2][1] + [1]; + data.shape_hessians[dof_index1][q][2][1] + [2] = + data.shape_hessians[dof_index1][q][2][2] + [1]; + + data.shape_hessians[dof_index1][q][0][2] + [2] = + data.shape_hessians[dof_index1][q][2][2] + [0]; + data.shape_hessians[dof_index1][q][1][2] + [2] = + data.shape_hessians[dof_index1][q][2][2] + [1]; + data.shape_hessians[dof_index1][q][2][2] + [2] = + 8.0 * polyx[i][0] * polyy[j][0] * + polyz[k][3]; + + + // Type 2: + const unsigned int dof_index2_1( + cell_type2_offset1 + shift_ijk); + const unsigned int dof_index2_2( + cell_type2_offset2 + shift_ijk); + + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; + ++d2) + { + data + .shape_hessians[dof_index2_1][q] + [0][d1][d2] = + data + .shape_hessians[dof_index1][q] + [0][d1][d2]; + data + .shape_hessians[dof_index2_1][q] + [1][d1][d2] = + -1.0 * + data + .shape_hessians[dof_index1][q] + [1][d1][d2]; + data + .shape_hessians[dof_index2_1][q] + [2][d1][d2] = + data + .shape_hessians[dof_index1][q] + [2][d1][d2]; + + data + .shape_hessians[dof_index2_2][q] + [0][d1][d2] = + data + .shape_hessians[dof_index1][q] + [0][d1][d2]; + data + .shape_hessians[dof_index2_2][q] + [1][d1][d2] = + -1.0 * + data + .shape_hessians[dof_index1][q] + [1][d1][d2]; + data + .shape_hessians[dof_index2_2][q] + [2][d1][d2] = + -1.0 * + data + .shape_hessians[dof_index1][q] + [2][d1][d2]; + } + } + } + // Type 3: (note we re-use k and j for + // convenience): + const unsigned int shift_ij( + j + shift_j); // here we've subbed j for i, + // k for j. + const unsigned int dof_index3_1( + cell_type3_offset1 + shift_ij); + const unsigned int dof_index3_2( + cell_type3_offset2 + shift_ij); + const unsigned int dof_index3_3( + cell_type3_offset3 + shift_ij); + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; + ++d3) + { + data + .shape_hessians[dof_index3_1][q] + [d1][d2][d3] = + 0.0; + data + .shape_hessians[dof_index3_2][q] + [d1][d2][d3] = + 0.0; + data + .shape_hessians[dof_index3_3][q] + [d1][d2][d3] = + 0.0; + } + } + } + data + .shape_hessians[dof_index3_1][q][0][1][1] = + 4.0 * polyy[j][2] * polyz[k][0]; + data + .shape_hessians[dof_index3_1][q][0][1][2] = + 4.0 * polyy[j][1] * polyz[k][1]; + + data + .shape_hessians[dof_index3_1][q][0][2][1] = + data + .shape_hessians[dof_index3_1][q][0][1][2]; + data + .shape_hessians[dof_index3_1][q][0][2][2] = + 4.0 * polyy[j][0] * polyz[k][2]; + + + data + .shape_hessians[dof_index3_2][q][1][0][0] = + 4.0 * polyx[j][2] * polyz[k][0]; + data + .shape_hessians[dof_index3_2][q][1][0][2] = + 4.0 * polyx[j][1] * polyz[k][1]; + + data + .shape_hessians[dof_index3_2][q][1][2][0] = + data + .shape_hessians[dof_index3_2][q][1][0][2]; + data + .shape_hessians[dof_index3_2][q][1][2][2] = + 4.0 * polyx[j][0] * polyz[k][2]; + + + data + .shape_hessians[dof_index3_3][q][2][0][0] = + 4.0 * polyx[j][2] * polyy[k][0]; + data + .shape_hessians[dof_index3_3][q][2][0][1] = + 4.0 * polyx[j][1] * polyy[k][1]; + + data + .shape_hessians[dof_index3_3][q][2][1][0] = + data + .shape_hessians[dof_index3_3][q][2][0][1]; + data + .shape_hessians[dof_index3_3][q][2][1][1] = + 4.0 * polyx[j][0] * polyy[k][2]; + } + } + } } } } @@ -1155,7 +1518,8 @@ FE_NedelecSZ::fill_edge_values( { case 2: { - if ((flags & (update_values | update_gradients)) != 0u) + if ((flags & (update_values | update_gradients | update_hessians)) != + 0u) { // Define an edge numbering so that each edge, E_{m} = [e^{m}_{1}, // e^{m}_{2}] e1 = higher global numbering of the two local @@ -1182,7 +1546,7 @@ FE_NedelecSZ::fill_edge_values( } } - // Define \sigma_{m} = sigma_{e^{m}_{2}} - sigma_{e^{m}_{2}} + // Define \sigma_{m} = sigma_{e^{m}_{2}} - sigma_{e^{m}_{1}} // \lambda_{m} = \lambda_{e^{m}_{1}} + \lambda_{e^{m}_{2}} // // To help things, in fe_data, we have precomputed (sigma_{i} - @@ -1249,7 +1613,9 @@ FE_NedelecSZ::fill_edge_values( // derivatives of the 1d polynomials, but only first derivatives // for the shape values. const unsigned int poly_length( - (flags & update_gradients) != 0u ? 3 : 2); + (flags & update_hessians) != 0u ? + 4 : + ((flags & update_gradients) != 0u ? 3 : 2)); for (unsigned int m = 0; m < lines_per_cell; ++m) { @@ -1263,7 +1629,7 @@ FE_NedelecSZ::fill_edge_values( { // Compute all required 1d polynomials and their // derivatives, starting at degree 2. e.g. to access - // L'_{3}(2x-1) use polyx[1][1]. + // L'_{i+2}(edge_sigma) use polyx[i][1]. IntegratedLegendrePolynomials[i + 1].value( edge_sigma_values[m][q], poly[i - 1]); } @@ -1333,6 +1699,54 @@ FE_NedelecSZ::fill_edge_values( edge_lambda_values[m][q] * poly[poly_index][2]; } } + if (flags & update_hessians) + { + // Lowest order edge shape function + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + fe_data.shape_hessians[shift_m][q][d1][d2] + [d3] = 0; + } + } + } + + // Higher order edge shape function + for (unsigned int i = 0; i < degree; ++i) + { + const unsigned int dof_index(i + 1 + shift_m); + + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + fe_data.shape_hessians[dof_index][q] + [d1][d2][d3] = + edge_sigma_grads[m][d1] * + edge_sigma_grads[m][d2] * + edge_sigma_grads[m][d3] * + poly[i][3] * + edge_lambda_values[m][q] + + poly[i][2] * + (edge_sigma_grads[m][d1] * + edge_sigma_grads[m][d2] * + edge_lambda_grads[m][d3] + + edge_sigma_grads[m][d3] * + edge_sigma_grads[m][d1] * + edge_lambda_grads[m][d2] + + edge_sigma_grads[m][d3] * + edge_sigma_grads[m][d2] * + edge_lambda_grads[m][d1]); + } + } + } + } + } } } } @@ -1340,7 +1754,8 @@ FE_NedelecSZ::fill_edge_values( } case 3: { - if ((flags & (update_values | update_gradients)) != 0u) + if ((flags & (update_values | update_gradients | update_hessians)) != + 0u) { // Define an edge numbering so that each edge, E_{m} = [e^{m}_{1}, // e^{m}_{2}] e1 = higher global numbering of the two local @@ -1367,7 +1782,7 @@ FE_NedelecSZ::fill_edge_values( } } - // Define \sigma_{m} = sigma_{e^{m}_{2}} - sigma_{e^{m}_{2}} + // Define \sigma_{m} = sigma_{e^{m}_{1}} - sigma_{e^{m}_{2}} // \lambda_{m} = \lambda_{e^{m}_{1}} + \lambda_{e^{m}_{2}} // // To help things, in fe_data, we have precomputed (sigma_{i} - @@ -1437,7 +1852,10 @@ FE_NedelecSZ::fill_edge_values( // derivatives of the 1d polynomials, but only first derivatives // for the shape values. const unsigned int poly_length( - (flags & update_gradients) != 0u ? 3 : 2); + (flags & update_hessians) != 0u ? + 4 : + ((flags & update_gradients) != 0u ? 3 : 2)); + std::vector> poly( degree, std::vector(poly_length)); for (unsigned int m = 0; m < lines_per_cell; ++m) @@ -1446,6 +1864,7 @@ FE_NedelecSZ::fill_edge_values( for (unsigned int q = 0; q < n_q_points; ++q) { // precompute values of all 1d polynomials required: + // for example poly[i][1] = L'_{i+2}(edge_sigma_values) if (degree > 0) { for (unsigned int i = 0; i < degree; ++i) @@ -1456,13 +1875,14 @@ FE_NedelecSZ::fill_edge_values( } if ((flags & update_values) != 0u) { - // Lowest order shape functions: + // Lowest order edge shape functions: for (unsigned int d = 0; d < dim; ++d) { fe_data.shape_values[shift_m][q][d] = 0.5 * edge_sigma_grads[m][d] * edge_lambda_values[m][q]; } + // Higher order edge shape functions if (degree > 0) { for (unsigned int i = 0; i < degree; ++i) @@ -1480,7 +1900,7 @@ FE_NedelecSZ::fill_edge_values( } if ((flags & update_gradients) != 0u) { - // Lowest order shape functions: + // Lowest order edge shape functions: for (unsigned int d1 = 0; d1 < dim; ++d1) { for (unsigned int d2 = 0; d2 < dim; ++d2) @@ -1490,6 +1910,7 @@ FE_NedelecSZ::fill_edge_values( edge_lambda_grads[m][q][d2]; } } + // Higher order edge shape functions if (degree > 0) { for (unsigned int i = 0; i < degree; ++i) @@ -1519,6 +1940,74 @@ FE_NedelecSZ::fill_edge_values( } } } + if (flags & update_hessians) + { + // Lowest order edge shape functions: + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + fe_data.shape_hessians[shift_m][q][d1][d2] + [d3] = + 0.5 * edge_sigma_grads[m][d1] * + edge_lambda_gradgrads_3d[m][d3][d2]; + } + } + } + + // Higher order edge shape functions + if (degree > 0) + { + for (unsigned int i = 0; i < degree; ++i) + { + const unsigned int dof_index(i + 1 + shift_m); + + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; + ++d3) + { + fe_data + .shape_hessians[dof_index][q] + [d1][d2][d3] = + edge_sigma_grads[m][d1] * + edge_sigma_grads[m][d2] * + edge_sigma_grads[m][d3] * + poly[i][3] * + edge_lambda_values[m][q] + + poly[i][2] * + (edge_sigma_grads[m][d1] * + edge_sigma_grads[m][d2] * + edge_lambda_grads[m][q] + [d3] + + edge_sigma_grads[m][d3] * + edge_sigma_grads[m][d1] * + edge_lambda_grads[m][q] + [d2] + + edge_sigma_grads[m][d3] * + edge_sigma_grads[m][d2] * + edge_lambda_grads[m][q] + [d1]) + + poly[i][1] * + (edge_sigma_grads[m][d1] * + edge_lambda_gradgrads_3d + [m][d3][d2] + + edge_sigma_grads[m][d2] * + edge_lambda_gradgrads_3d + [m][d3][d1] + + edge_sigma_grads[m][d3] * + edge_lambda_gradgrads_3d + [m][d2][d1]); + } + } + } + } + } + } } } } @@ -1562,7 +2051,7 @@ FE_NedelecSZ::fill_face_values( { const UpdateFlags flags(fe_data.update_each); - if ((flags & (update_values | update_gradients)) != 0u) + if ((flags & (update_values | update_gradients | update_hessians)) != 0u) { const unsigned int n_q_points = quadrature.size(); @@ -1688,8 +2177,10 @@ FE_NedelecSZ::fill_face_values( } } // Now can generate the basis - const unsigned int poly_length((flags & update_gradients) != 0u ? 3 : + const unsigned int poly_length((flags & update_hessians) != 0u ? 4 : + (flags & update_gradients) != 0u ? 3 : 2); + std::vector> polyxi( degree, std::vector(poly_length)); std::vector> polyeta( @@ -1899,13 +2390,210 @@ FE_NedelecSZ::fill_face_values( } } } + if ((flags & update_hessians) != 0u) + { + for (unsigned int j = 0; j < degree; ++j) + { + const unsigned int shift_j(j * degree); + for (unsigned int i = 0; i < degree; ++i) + { + const unsigned int shift_ij(shift_j + i); + + // Type 1: + const unsigned int dof_index1(face_type1_offset + + shift_ij); + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + fe_data.shape_hessians[dof_index1][q] + [d1][d2][d3] = + polyxi[i][1] * + face_xi_grads[m][d3] * + (face_eta_grads[m][d1] * + (polyeta[j][2] * + face_eta_grads[m][d2] * + face_lambda_values[m][q] + + polyeta[j][1] * + face_lambda_grads[m][d2]) + + polyeta[j][1] * + face_eta_grads[m][d2] * + face_lambda_grads[m][d1]) + + polyxi[i][0] * + (polyeta[j][3] * + face_eta_grads[m][d1] * + face_eta_grads[m][d2] * + face_eta_grads[m][d3] * + face_lambda_values[m][q] + + polyeta[j][2] * + (face_eta_grads[m][d1] * + face_eta_grads[m][d2] * + face_lambda_grads[m][d3] + + face_eta_grads[m][d3] * + (face_eta_grads[m][d1] * + face_lambda_grads[m] + [d2] + + face_eta_grads[m][d2] * + face_lambda_grads + [m][d1]))) + + (polyxi[i][1] * polyeta[j][1] * + face_eta_grads[m][d3] + + polyxi[i][2] * polyeta[j][0] * + face_xi_grads[m][d3]) * + (face_xi_grads[m][d1] * + face_lambda_grads[m][d2] + + face_xi_grads[m][d2] * + face_lambda_grads[m][d1]) + + face_lambda_grads[m][d3] * + (polyxi[i][2] * polyeta[j][0] * + face_xi_grads[m][d1] * + face_xi_grads[m][d2] + + polyxi[i][1] * polyeta[j][1] * + (face_xi_grads[m][d1] * + face_eta_grads[m][d2] + + face_xi_grads[m][d2] * + face_eta_grads[m][d1])) + + face_lambda_values[m][q] * + (polyxi[i][3] * polyeta[j][0] * + face_xi_grads[m][d1] * + face_xi_grads[m][d2] * + face_xi_grads[m][d3] + + polyxi[i][1] * polyeta[j][2] * + face_eta_grads[m][d3] * + (face_xi_grads[m][d1] * + face_eta_grads[m][d2] + + face_xi_grads[m][d2] * + face_eta_grads[m][d1]) + + polyxi[i][2] * polyeta[j][1] * + (face_xi_grads[m][d3] * + face_xi_grads[m][d2] * + face_eta_grads[m][d1] + + face_xi_grads[m][d1] * + (face_xi_grads[m][d2] * + face_eta_grads[m][d3] + + face_xi_grads[m][d3] * + face_eta_grads[m][d2]))); + } + } + } + + // Type 2: + const unsigned int dof_index2(face_type2_offset + + shift_ij); + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + fe_data.shape_hessians[dof_index2][q] + [d1][d2][d3] = + face_xi_grads[m][d1] * + (polyxi[i][1] * polyeta[j][1] * + (face_eta_grads[m][d2] * + face_lambda_grads[m][d3] + + face_eta_grads[m][d3] * + face_lambda_grads[m][d2]) + + polyxi[i][2] * polyeta[j][0] * + (face_xi_grads[m][d2] * + face_lambda_grads[m][d3] + + face_xi_grads[m][d3] * + face_lambda_grads[m][d2]) + + face_lambda_values[m][q] * + (face_eta_grads[m][d2] * + (polyxi[i][1] * + polyeta[j][2] * + face_eta_grads[m][d3] + + polyxi[i][2] * + polyeta[j][1] * + face_xi_grads[m][d3]) + + face_xi_grads[m][d2] * + (polyxi[i][2] * + polyeta[j][1] * + face_eta_grads[m][d3] + + polyxi[i][3] * + polyeta[j][0] * + face_xi_grads[m][d3]))) - + polyxi[i][0] * + face_eta_grads[m][d1] * + (face_eta_grads[m][d2] * + (polyeta[j][3] * + face_eta_grads[m][d3] * + face_lambda_values[m][q] + + polyeta[j][2] * + face_lambda_grads[m][d3]) + + polyeta[j][2] * + face_eta_grads[m][d3] * + face_lambda_grads[m][d2]) - + face_eta_grads[m][d1] * + (polyxi[i][1] * + face_xi_grads[m][d3] * + (polyeta[j][2] * + face_eta_grads[m][d2] * + face_lambda_values[m][q] + + polyeta[j][1] * + face_lambda_grads[m][d2]) + + face_xi_grads[m][d2] * + (polyxi[i][1] * + (polyeta[j][2] * + face_eta_grads[m][d3] * + face_lambda_values[m] + [q] + + polyeta[j][1] * + face_lambda_grads[m] + [d3]) + + polyxi[i][2] * polyeta[j][1] * + face_xi_grads[m][d3] * + face_lambda_values[m][q])); + } + } + } + } + // Type 3: + const unsigned int dof_index3_1(face_type3_offset1 + + j); + const unsigned int dof_index3_2(face_type3_offset2 + + j); + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + fe_data.shape_hessians[dof_index3_1][q] + [d1][d2][d3] = + face_xi_grads[m][d1] * + (face_eta_grads[m][d2] * + (polyeta[j][2] * + face_eta_grads[m][d3] * + face_lambda_values[m][q] + + polyeta[j][1] * + face_lambda_grads[m][d3]) + + face_lambda_grads[m][d2] * + polyeta[j][1] * + face_eta_grads[m][d3]); + + fe_data.shape_hessians[dof_index3_2][q] + [d1][d2][d3] = + face_eta_grads[m][d1] * + (face_xi_grads[m][d2] * + (polyxi[j][2] * + face_xi_grads[m][d3] * + face_lambda_values[m][q] + + polyxi[j][1] * + face_lambda_grads[m][d3]) + + face_lambda_grads[m][d2] * + polyxi[j][1] * face_xi_grads[m][d3]); + } + } + } + } + } } } } - if ((flags & update_hessians) != 0u) - { - Assert(false, ExcNotImplemented()); - } } } @@ -1976,6 +2664,7 @@ FE_NedelecSZ::fill_fe_values( } } } + if ((flags & update_gradients) != 0u) { // Now have all shape_grads stored on the reference cell. @@ -2021,6 +2710,72 @@ FE_NedelecSZ::fill_fe_values( } } } + + if (flags & update_hessians) + { + // Now have all shape_grads stored on the reference cell. + // Must now transform to the physical cell. + std::vector> input(n_q_points); + std::vector> transformed_shape_hessians(n_q_points); + for (unsigned int dof = 0; dof < this->n_dofs_per_cell(); ++dof) + { + for (unsigned int q = 0; q < n_q_points; ++q) + { + input[q] = fe_data.shape_hessians[dof][q]; + } + mapping.transform(make_array_view(input), + mapping_covariant_hessian, + mapping_internal, + make_array_view(transformed_shape_hessians)); + + const unsigned int first = + data.shape_function_to_row_table[dof * this->n_components() + + this->get_nonzero_components(dof) + .first_selected_component()]; + + for (unsigned int q = 0; q < n_q_points; ++q) + { + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + for (unsigned int d4 = 0; d4 < dim; ++d4) + { + transformed_shape_hessians[q][d1][d3][d4] -= + (data.shape_values(first + d2, q) * + mapping_data + .jacobian_pushed_forward_2nd_derivatives + [q][d2][d1][d3][d4]) + + (data.shape_gradients[first + d1][q][d2] * + mapping_data + .jacobian_pushed_forward_grads[q][d2][d3] + [d4]) + + (data.shape_gradients[first + d2][q][d3] * + mapping_data + .jacobian_pushed_forward_grads[q][d2][d1] + [d4]) + + (data.shape_gradients[first + d2][q][d4] * + mapping_data + .jacobian_pushed_forward_grads[q][d2][d3] + [d1]); + } + } + } + } + } + + for (unsigned int q = 0; q < n_q_points; ++q) + { + for (unsigned int d = 0; d < dim; ++d) + { + data.shape_hessians[first + d][q] = + transformed_shape_hessians[q][d]; + } + } + } + } } @@ -2157,6 +2912,70 @@ FE_NedelecSZ::fill_fe_face_values( } } } + if (flags & update_hessians) + { + // Now have all shape_grads stored on the reference cell. + // Must now transform to the physical cell. + std::vector> input(n_q_points); + std::vector> transformed_shape_hessians(n_q_points); + for (unsigned int dof = 0; dof < this->n_dofs_per_cell(); ++dof) + { + for (unsigned int q = 0; q < n_q_points; ++q) + input[q] = fe_data.shape_hessians[dof][offset + q]; + + mapping.transform(input, + mapping_covariant_hessian, + mapping_internal, + make_array_view(transformed_shape_hessians)); + + const unsigned int first = + data.shape_function_to_row_table[dof * this->n_components() + + this->get_nonzero_components(dof) + .first_selected_component()]; + + for (unsigned int q = 0; q < n_q_points; ++q) + { + for (unsigned int d1 = 0; d1 < dim; ++d1) + { + for (unsigned int d2 = 0; d2 < dim; ++d2) + { + for (unsigned int d3 = 0; d3 < dim; ++d3) + { + for (unsigned int d4 = 0; d4 < dim; ++d4) + { + transformed_shape_hessians[q][d1][d3][d4] -= + (data.shape_values(first + d2, q) * + mapping_data + .jacobian_pushed_forward_2nd_derivatives + [q][d2][d1][d3][d4]) + + (data.shape_gradients[first + d1][q][d2] * + mapping_data + .jacobian_pushed_forward_grads[q][d2][d3] + [d4]) + + (data.shape_gradients[first + d2][q][d3] * + mapping_data + .jacobian_pushed_forward_grads[q][d2][d1] + [d4]) + + (data.shape_gradients[first + d2][q][d4] * + mapping_data + .jacobian_pushed_forward_grads[q][d2][d3] + [d1]); + } + } + } + } + } + + for (unsigned int q = 0; q < n_q_points; ++q) + { + for (unsigned int d = 0; d < dim; ++d) + { + data.shape_hessians[first + d][q] = + transformed_shape_hessians[q][d]; + } + } + } + } } diff --git a/tests/fe/fe_nedelec_sz_hessian_divergence_theorem.cc b/tests/fe/fe_nedelec_sz_hessian_divergence_theorem.cc new file mode 100644 index 0000000000..acc4afac5e --- /dev/null +++ b/tests/fe/fe_nedelec_sz_hessian_divergence_theorem.cc @@ -0,0 +1,172 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2003 - 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Observe how the values of the shape functions change as we make a +// cell smaller and smaller. + +#include +#include + +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +template +Tensor<1, dim> +ones() +{ + Tensor<1, dim> result; + for (unsigned int i = 0; i < dim; ++i) + result[i] = 1.0; + return result; +} + +template +void +test(const Triangulation &tr, + const FiniteElement &fe, + const double tolerance) +{ + DoFHandler dof(tr); + dof.distribute_dofs(fe); + + std::stringstream ss; + + deallog << "FE=" << fe.get_name() << std::endl; + + const QGauss quadrature(4); + FEValues fe_values(fe, + quadrature, + update_hessians | update_quadrature_points | + update_JxW_values); + + const QGauss face_quadrature(4); + FEFaceValues fe_face_values(fe, + face_quadrature, + update_gradients | update_quadrature_points | + update_normal_vectors | update_JxW_values); + + for (const auto &cell : dof.active_cell_iterators()) + { + fe_values.reinit(cell); + + deallog << "Cell nodes:" << std::endl; + for (const unsigned int i : GeometryInfo::vertex_indices()) + deallog << i << ": ( " << cell->vertex(i) << " )" << std::endl; + + bool cell_ok = true; + + for (unsigned int c = 0; c < fe.n_components(); ++c) + { + FEValuesExtractors::Scalar single_component(c); + + for (unsigned int i = 0; i < fe_values.dofs_per_cell; ++i) + { + ss << "component=" << c << ", dof=" << i << std::endl; + + Tensor<2, dim> bulk_integral; + for (const auto q : fe_values.quadrature_point_indices()) + { + bulk_integral += fe_values[single_component].hessian(i, q) * + fe_values.JxW(q); + } + + Tensor<2, dim> boundary_integral; + for (const unsigned int face : GeometryInfo::face_indices()) + { + fe_face_values.reinit(cell, face); + for (const auto q : fe_face_values.quadrature_point_indices()) + { + Tensor<1, dim> gradient = + fe_face_values[single_component].gradient(i, q); + Tensor<2, dim> gradient_normal_outer_prod = + outer_product(gradient, + fe_face_values.normal_vector(q)); + boundary_integral += + gradient_normal_outer_prod * fe_face_values.JxW(q); + } + } + + if ((bulk_integral - boundary_integral).norm_square() > + tolerance * (bulk_integral.norm() + boundary_integral.norm())) + { + deallog << "Failed:" << std::endl; + deallog << ss.str() << std::endl; + deallog << " bulk integral=" << bulk_integral << std::endl; + deallog << "boundary integral=" << boundary_integral + << std::endl; + deallog + << "Error! difference between bulk and surface integrals is greater than " + << tolerance << "!\n\n" + << std::endl; + cell_ok = false; + } + + ss.str(""); + } + } + + deallog << (cell_ok ? "OK: cell bulk and boundary integrals match...\n" : + "Failed divergence test...\n") + << std::endl; + } +} + + + +template +void +test_hyper_ball(const double tolerance) +{ + Triangulation tr; + GridGenerator::hyper_ball(tr); + + static const SphericalManifold boundary; + tr.set_manifold(0, boundary); + + tr.refine_global(1); + + FE_NedelecSZ fe(1); + test(tr, fe, tolerance); +} + + +int +main() +{ + initlog(); + deallog << std::setprecision(8); + + test_hyper_ball<2>(1e-6); + test_hyper_ball<3>(1e-6); + + deallog << "done..." << std::endl; +} diff --git a/tests/fe/fe_nedelec_sz_hessian_divergence_theorem.output b/tests/fe/fe_nedelec_sz_hessian_divergence_theorem.output new file mode 100644 index 0000000000..08332ff31e --- /dev/null +++ b/tests/fe/fe_nedelec_sz_hessian_divergence_theorem.output @@ -0,0 +1,760 @@ + +DEAL::FE=FE_NedelecSZ<2>(1) +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 -0.70710678 ) +DEAL::1: ( -1.1102230e-16 -1.0000000 ) +DEAL::2: ( -0.50000000 -0.50000000 ) +DEAL::3: ( -5.5511151e-17 -0.64644661 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -1.1102230e-16 -1.0000000 ) +DEAL::1: ( 0.70710678 -0.70710678 ) +DEAL::2: ( -5.5511151e-17 -0.64644661 ) +DEAL::3: ( 0.50000000 -0.50000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.50000000 -0.50000000 ) +DEAL::1: ( -5.5511151e-17 -0.64644661 ) +DEAL::2: ( -0.29289322 -0.29289322 ) +DEAL::3: ( 0.0000000 -0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -5.5511151e-17 -0.64644661 ) +DEAL::1: ( 0.50000000 -0.50000000 ) +DEAL::2: ( 0.0000000 -0.29289322 ) +DEAL::3: ( 0.29289322 -0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 -0.70710678 ) +DEAL::1: ( -0.50000000 -0.50000000 ) +DEAL::2: ( -1.0000000 -1.1102230e-16 ) +DEAL::3: ( -0.64644661 -5.5511151e-17 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.50000000 -0.50000000 ) +DEAL::1: ( -0.29289322 -0.29289322 ) +DEAL::2: ( -0.64644661 -5.5511151e-17 ) +DEAL::3: ( -0.29289322 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -1.0000000 -1.1102230e-16 ) +DEAL::1: ( -0.64644661 -5.5511151e-17 ) +DEAL::2: ( -0.70710678 0.70710678 ) +DEAL::3: ( -0.50000000 0.50000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.64644661 -5.5511151e-17 ) +DEAL::1: ( -0.29289322 0.0000000 ) +DEAL::2: ( -0.50000000 0.50000000 ) +DEAL::3: ( -0.29289322 0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.29289322 -0.29289322 ) +DEAL::1: ( 0.0000000 -0.29289322 ) +DEAL::2: ( -0.29289322 0.0000000 ) +DEAL::3: ( 0.0000000 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.29289322 ) +DEAL::1: ( 0.29289322 -0.29289322 ) +DEAL::2: ( 0.0000000 0.0000000 ) +DEAL::3: ( 0.29289322 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.29289322 0.0000000 ) +DEAL::1: ( 0.0000000 0.0000000 ) +DEAL::2: ( -0.29289322 0.29289322 ) +DEAL::3: ( 0.0000000 0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.0000000 ) +DEAL::1: ( 0.29289322 0.0000000 ) +DEAL::2: ( 0.0000000 0.29289322 ) +DEAL::3: ( 0.29289322 0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.70710678 -0.70710678 ) +DEAL::1: ( 1.0000000 -1.1102230e-16 ) +DEAL::2: ( 0.50000000 -0.50000000 ) +DEAL::3: ( 0.64644661 -5.5511151e-17 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 1.0000000 -1.1102230e-16 ) +DEAL::1: ( 0.70710678 0.70710678 ) +DEAL::2: ( 0.64644661 -5.5511151e-17 ) +DEAL::3: ( 0.50000000 0.50000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.50000000 -0.50000000 ) +DEAL::1: ( 0.64644661 -5.5511151e-17 ) +DEAL::2: ( 0.29289322 -0.29289322 ) +DEAL::3: ( 0.29289322 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.64644661 -5.5511151e-17 ) +DEAL::1: ( 0.50000000 0.50000000 ) +DEAL::2: ( 0.29289322 0.0000000 ) +DEAL::3: ( 0.29289322 0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 0.70710678 ) +DEAL::1: ( -0.50000000 0.50000000 ) +DEAL::2: ( -1.1102230e-16 1.0000000 ) +DEAL::3: ( -5.5511151e-17 0.64644661 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.50000000 0.50000000 ) +DEAL::1: ( -0.29289322 0.29289322 ) +DEAL::2: ( -5.5511151e-17 0.64644661 ) +DEAL::3: ( 0.0000000 0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -1.1102230e-16 1.0000000 ) +DEAL::1: ( -5.5511151e-17 0.64644661 ) +DEAL::2: ( 0.70710678 0.70710678 ) +DEAL::3: ( 0.50000000 0.50000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -5.5511151e-17 0.64644661 ) +DEAL::1: ( 0.0000000 0.29289322 ) +DEAL::2: ( 0.50000000 0.50000000 ) +DEAL::3: ( 0.29289322 0.29289322 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::FE=FE_NedelecSZ<3>(1) +DEAL::Cell nodes: +DEAL::0: ( -0.21132487 -0.21132487 -0.21132487 ) +DEAL::1: ( 0.0000000 -0.21132487 -0.21132487 ) +DEAL::2: ( -0.21132487 0.0000000 -0.21132487 ) +DEAL::3: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::4: ( -0.21132487 -0.21132487 0.0000000 ) +DEAL::5: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::6: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::7: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.21132487 -0.21132487 ) +DEAL::1: ( 0.21132487 -0.21132487 -0.21132487 ) +DEAL::2: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::3: ( 0.21132487 0.0000000 -0.21132487 ) +DEAL::4: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::5: ( 0.21132487 -0.21132487 0.0000000 ) +DEAL::6: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::7: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.21132487 0.0000000 -0.21132487 ) +DEAL::1: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::2: ( -0.21132487 0.21132487 -0.21132487 ) +DEAL::3: ( 0.0000000 0.21132487 -0.21132487 ) +DEAL::4: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::5: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::6: ( -0.21132487 0.21132487 0.0000000 ) +DEAL::7: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::1: ( 0.21132487 0.0000000 -0.21132487 ) +DEAL::2: ( 0.0000000 0.21132487 -0.21132487 ) +DEAL::3: ( 0.21132487 0.21132487 -0.21132487 ) +DEAL::4: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::5: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::6: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::7: ( 0.21132487 0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.21132487 -0.21132487 0.0000000 ) +DEAL::1: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::2: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::3: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::4: ( -0.21132487 -0.21132487 0.21132487 ) +DEAL::5: ( 0.0000000 -0.21132487 0.21132487 ) +DEAL::6: ( -0.21132487 0.0000000 0.21132487 ) +DEAL::7: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::1: ( 0.21132487 -0.21132487 0.0000000 ) +DEAL::2: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::3: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::4: ( 0.0000000 -0.21132487 0.21132487 ) +DEAL::5: ( 0.21132487 -0.21132487 0.21132487 ) +DEAL::6: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::7: ( 0.21132487 0.0000000 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::1: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::2: ( -0.21132487 0.21132487 0.0000000 ) +DEAL::3: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::4: ( -0.21132487 0.0000000 0.21132487 ) +DEAL::5: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::6: ( -0.21132487 0.21132487 0.21132487 ) +DEAL::7: ( 0.0000000 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.0000000 -1.3877788e-17 ) +DEAL::1: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::2: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::3: ( 0.21132487 0.21132487 0.0000000 ) +DEAL::4: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::5: ( 0.21132487 0.0000000 0.21132487 ) +DEAL::6: ( 0.0000000 0.21132487 0.21132487 ) +DEAL::7: ( 0.21132487 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.57735027 -0.57735027 -0.57735027 ) +DEAL::1: ( 0.0000000 -0.70710678 -0.70710678 ) +DEAL::2: ( -0.70710678 0.0000000 -0.70710678 ) +DEAL::3: ( 0.0000000 0.0000000 -1.0000000 ) +DEAL::4: ( -0.39433757 -0.39433757 -0.39433757 ) +DEAL::5: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::6: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::7: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.70710678 -0.70710678 ) +DEAL::1: ( 0.57735027 -0.57735027 -0.57735027 ) +DEAL::2: ( 0.0000000 0.0000000 -1.0000000 ) +DEAL::3: ( 0.70710678 0.0000000 -0.70710678 ) +DEAL::4: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::5: ( 0.39433757 -0.39433757 -0.39433757 ) +DEAL::6: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::7: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 0.0000000 -0.70710678 ) +DEAL::1: ( 0.0000000 0.0000000 -1.0000000 ) +DEAL::2: ( -0.57735027 0.57735027 -0.57735027 ) +DEAL::3: ( 0.0000000 0.70710678 -0.70710678 ) +DEAL::4: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::5: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::6: ( -0.39433757 0.39433757 -0.39433757 ) +DEAL::7: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.0000000 -1.0000000 ) +DEAL::1: ( 0.70710678 0.0000000 -0.70710678 ) +DEAL::2: ( 0.0000000 0.70710678 -0.70710678 ) +DEAL::3: ( 0.57735027 0.57735027 -0.57735027 ) +DEAL::4: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::5: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::6: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::7: ( 0.39433757 0.39433757 -0.39433757 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.39433757 -0.39433757 -0.39433757 ) +DEAL::1: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::2: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::3: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::4: ( -0.21132487 -0.21132487 -0.21132487 ) +DEAL::5: ( 0.0000000 -0.21132487 -0.21132487 ) +DEAL::6: ( -0.21132487 0.0000000 -0.21132487 ) +DEAL::7: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::1: ( 0.39433757 -0.39433757 -0.39433757 ) +DEAL::2: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::3: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::4: ( 0.0000000 -0.21132487 -0.21132487 ) +DEAL::5: ( 0.21132487 -0.21132487 -0.21132487 ) +DEAL::6: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::7: ( 0.21132487 0.0000000 -0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::1: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::2: ( -0.39433757 0.39433757 -0.39433757 ) +DEAL::3: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::4: ( -0.21132487 0.0000000 -0.21132487 ) +DEAL::5: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::6: ( -0.21132487 0.21132487 -0.21132487 ) +DEAL::7: ( 0.0000000 0.21132487 -0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.0000000 -0.60566243 ) +DEAL::1: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::2: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::3: ( 0.39433757 0.39433757 -0.39433757 ) +DEAL::4: ( 0.0000000 0.0000000 -0.21132487 ) +DEAL::5: ( 0.21132487 0.0000000 -0.21132487 ) +DEAL::6: ( 0.0000000 0.21132487 -0.21132487 ) +DEAL::7: ( 0.21132487 0.21132487 -0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.57735027 -0.57735027 -0.57735027 ) +DEAL::1: ( 0.70710678 0.0000000 -0.70710678 ) +DEAL::2: ( 0.39433757 -0.39433757 -0.39433757 ) +DEAL::3: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::4: ( 0.70710678 -0.70710678 0.0000000 ) +DEAL::5: ( 1.0000000 0.0000000 0.0000000 ) +DEAL::6: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::7: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.70710678 0.0000000 -0.70710678 ) +DEAL::1: ( 0.57735027 0.57735027 -0.57735027 ) +DEAL::2: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::3: ( 0.39433757 0.39433757 -0.39433757 ) +DEAL::4: ( 1.0000000 0.0000000 0.0000000 ) +DEAL::5: ( 0.70710678 0.70710678 0.0000000 ) +DEAL::6: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::7: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.39433757 -0.39433757 -0.39433757 ) +DEAL::1: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::2: ( 0.21132487 -0.21132487 -0.21132487 ) +DEAL::3: ( 0.21132487 0.0000000 -0.21132487 ) +DEAL::4: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::5: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::6: ( 0.21132487 -0.21132487 0.0000000 ) +DEAL::7: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.45921582 0.0000000 -0.45921582 ) +DEAL::1: ( 0.39433757 0.39433757 -0.39433757 ) +DEAL::2: ( 0.21132487 0.0000000 -0.21132487 ) +DEAL::3: ( 0.21132487 0.21132487 -0.21132487 ) +DEAL::4: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::5: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::6: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::7: ( 0.21132487 0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.70710678 -0.70710678 0.0000000 ) +DEAL::1: ( 1.0000000 0.0000000 0.0000000 ) +DEAL::2: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::3: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::4: ( 0.57735027 -0.57735027 0.57735027 ) +DEAL::5: ( 0.70710678 0.0000000 0.70710678 ) +DEAL::6: ( 0.39433757 -0.39433757 0.39433757 ) +DEAL::7: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 1.0000000 0.0000000 0.0000000 ) +DEAL::1: ( 0.70710678 0.70710678 0.0000000 ) +DEAL::2: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::3: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::4: ( 0.70710678 0.0000000 0.70710678 ) +DEAL::5: ( 0.57735027 0.57735027 0.57735027 ) +DEAL::6: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::7: ( 0.39433757 0.39433757 0.39433757 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::1: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::2: ( 0.21132487 -0.21132487 0.0000000 ) +DEAL::3: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::4: ( 0.39433757 -0.39433757 0.39433757 ) +DEAL::5: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::6: ( 0.21132487 -0.21132487 0.21132487 ) +DEAL::7: ( 0.21132487 0.0000000 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.60566243 0.0000000 -2.7755576e-17 ) +DEAL::1: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::2: ( 0.21132487 0.0000000 0.0000000 ) +DEAL::3: ( 0.21132487 0.21132487 0.0000000 ) +DEAL::4: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::5: ( 0.39433757 0.39433757 0.39433757 ) +DEAL::6: ( 0.21132487 0.0000000 0.21132487 ) +DEAL::7: ( 0.21132487 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.57735027 -0.57735027 0.57735027 ) +DEAL::1: ( 0.0000000 -0.70710678 0.70710678 ) +DEAL::2: ( -0.39433757 -0.39433757 0.39433757 ) +DEAL::3: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::4: ( -0.70710678 0.0000000 0.70710678 ) +DEAL::5: ( 0.0000000 0.0000000 1.0000000 ) +DEAL::6: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::7: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.70710678 0.70710678 ) +DEAL::1: ( 0.57735027 -0.57735027 0.57735027 ) +DEAL::2: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::3: ( 0.39433757 -0.39433757 0.39433757 ) +DEAL::4: ( 0.0000000 0.0000000 1.0000000 ) +DEAL::5: ( 0.70710678 0.0000000 0.70710678 ) +DEAL::6: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::7: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.39433757 -0.39433757 0.39433757 ) +DEAL::1: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::2: ( -0.21132487 -0.21132487 0.21132487 ) +DEAL::3: ( 0.0000000 -0.21132487 0.21132487 ) +DEAL::4: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::5: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::6: ( -0.21132487 0.0000000 0.21132487 ) +DEAL::7: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::1: ( 0.39433757 -0.39433757 0.39433757 ) +DEAL::2: ( 0.0000000 -0.21132487 0.21132487 ) +DEAL::3: ( 0.21132487 -0.21132487 0.21132487 ) +DEAL::4: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::5: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::6: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::7: ( 0.21132487 0.0000000 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 0.0000000 0.70710678 ) +DEAL::1: ( 0.0000000 0.0000000 1.0000000 ) +DEAL::2: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::3: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::4: ( -0.57735027 0.57735027 0.57735027 ) +DEAL::5: ( 0.0000000 0.70710678 0.70710678 ) +DEAL::6: ( -0.39433757 0.39433757 0.39433757 ) +DEAL::7: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.0000000 1.0000000 ) +DEAL::1: ( 0.70710678 0.0000000 0.70710678 ) +DEAL::2: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::3: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::4: ( 0.0000000 0.70710678 0.70710678 ) +DEAL::5: ( 0.57735027 0.57735027 0.57735027 ) +DEAL::6: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::7: ( 0.39433757 0.39433757 0.39433757 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::1: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::2: ( -0.21132487 0.0000000 0.21132487 ) +DEAL::3: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::4: ( -0.39433757 0.39433757 0.39433757 ) +DEAL::5: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::6: ( -0.21132487 0.21132487 0.21132487 ) +DEAL::7: ( 0.0000000 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -2.7755576e-17 -2.7755576e-17 0.60566243 ) +DEAL::1: ( 0.45921582 0.0000000 0.45921582 ) +DEAL::2: ( 0.0000000 0.0000000 0.21132487 ) +DEAL::3: ( 0.21132487 0.0000000 0.21132487 ) +DEAL::4: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::5: ( 0.39433757 0.39433757 0.39433757 ) +DEAL::6: ( 0.0000000 0.21132487 0.21132487 ) +DEAL::7: ( 0.21132487 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.57735027 -0.57735027 -0.57735027 ) +DEAL::1: ( -0.39433757 -0.39433757 -0.39433757 ) +DEAL::2: ( -0.70710678 0.0000000 -0.70710678 ) +DEAL::3: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::4: ( -0.70710678 -0.70710678 0.0000000 ) +DEAL::5: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::6: ( -1.0000000 0.0000000 0.0000000 ) +DEAL::7: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.39433757 -0.39433757 -0.39433757 ) +DEAL::1: ( -0.21132487 -0.21132487 -0.21132487 ) +DEAL::2: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::3: ( -0.21132487 0.0000000 -0.21132487 ) +DEAL::4: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::5: ( -0.21132487 -0.21132487 0.0000000 ) +DEAL::6: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::7: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 0.0000000 -0.70710678 ) +DEAL::1: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::2: ( -0.57735027 0.57735027 -0.57735027 ) +DEAL::3: ( -0.39433757 0.39433757 -0.39433757 ) +DEAL::4: ( -1.0000000 0.0000000 0.0000000 ) +DEAL::5: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::6: ( -0.70710678 0.70710678 0.0000000 ) +DEAL::7: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.45921582 0.0000000 -0.45921582 ) +DEAL::1: ( -0.21132487 0.0000000 -0.21132487 ) +DEAL::2: ( -0.39433757 0.39433757 -0.39433757 ) +DEAL::3: ( -0.21132487 0.21132487 -0.21132487 ) +DEAL::4: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::5: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::6: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::7: ( -0.21132487 0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 -0.70710678 0.0000000 ) +DEAL::1: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::2: ( -1.0000000 0.0000000 0.0000000 ) +DEAL::3: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::4: ( -0.57735027 -0.57735027 0.57735027 ) +DEAL::5: ( -0.39433757 -0.39433757 0.39433757 ) +DEAL::6: ( -0.70710678 0.0000000 0.70710678 ) +DEAL::7: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::1: ( -0.21132487 -0.21132487 0.0000000 ) +DEAL::2: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::3: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::4: ( -0.39433757 -0.39433757 0.39433757 ) +DEAL::5: ( -0.21132487 -0.21132487 0.21132487 ) +DEAL::6: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::7: ( -0.21132487 0.0000000 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -1.0000000 0.0000000 0.0000000 ) +DEAL::1: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::2: ( -0.70710678 0.70710678 0.0000000 ) +DEAL::3: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::4: ( -0.70710678 0.0000000 0.70710678 ) +DEAL::5: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::6: ( -0.57735027 0.57735027 0.57735027 ) +DEAL::7: ( -0.39433757 0.39433757 0.39433757 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.60566243 0.0000000 1.1102230e-16 ) +DEAL::1: ( -0.21132487 0.0000000 0.0000000 ) +DEAL::2: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::3: ( -0.21132487 0.21132487 0.0000000 ) +DEAL::4: ( -0.45921582 0.0000000 0.45921582 ) +DEAL::5: ( -0.21132487 0.0000000 0.21132487 ) +DEAL::6: ( -0.39433757 0.39433757 0.39433757 ) +DEAL::7: ( -0.21132487 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.57735027 -0.57735027 -0.57735027 ) +DEAL::1: ( 0.0000000 -0.70710678 -0.70710678 ) +DEAL::2: ( -0.39433757 -0.39433757 -0.39433757 ) +DEAL::3: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::4: ( -0.70710678 -0.70710678 0.0000000 ) +DEAL::5: ( 0.0000000 -1.0000000 0.0000000 ) +DEAL::6: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::7: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.70710678 -0.70710678 ) +DEAL::1: ( 0.57735027 -0.57735027 -0.57735027 ) +DEAL::2: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::3: ( 0.39433757 -0.39433757 -0.39433757 ) +DEAL::4: ( 0.0000000 -1.0000000 0.0000000 ) +DEAL::5: ( 0.70710678 -0.70710678 0.0000000 ) +DEAL::6: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::7: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.39433757 -0.39433757 -0.39433757 ) +DEAL::1: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::2: ( -0.21132487 -0.21132487 -0.21132487 ) +DEAL::3: ( 0.0000000 -0.21132487 -0.21132487 ) +DEAL::4: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::5: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::6: ( -0.21132487 -0.21132487 0.0000000 ) +DEAL::7: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -0.45921582 -0.45921582 ) +DEAL::1: ( 0.39433757 -0.39433757 -0.39433757 ) +DEAL::2: ( 0.0000000 -0.21132487 -0.21132487 ) +DEAL::3: ( 0.21132487 -0.21132487 -0.21132487 ) +DEAL::4: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::5: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::6: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::7: ( 0.21132487 -0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 -0.70710678 0.0000000 ) +DEAL::1: ( 0.0000000 -1.0000000 0.0000000 ) +DEAL::2: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::3: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::4: ( -0.57735027 -0.57735027 0.57735027 ) +DEAL::5: ( 0.0000000 -0.70710678 0.70710678 ) +DEAL::6: ( -0.39433757 -0.39433757 0.39433757 ) +DEAL::7: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 -1.0000000 0.0000000 ) +DEAL::1: ( 0.70710678 -0.70710678 0.0000000 ) +DEAL::2: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::3: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::4: ( 0.0000000 -0.70710678 0.70710678 ) +DEAL::5: ( 0.57735027 -0.57735027 0.57735027 ) +DEAL::6: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::7: ( 0.39433757 -0.39433757 0.39433757 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.45921582 -0.45921582 0.0000000 ) +DEAL::1: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::2: ( -0.21132487 -0.21132487 0.0000000 ) +DEAL::3: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::4: ( -0.39433757 -0.39433757 0.39433757 ) +DEAL::5: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::6: ( -0.21132487 -0.21132487 0.21132487 ) +DEAL::7: ( 0.0000000 -0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 2.7755576e-17 -0.60566243 -5.5511151e-17 ) +DEAL::1: ( 0.45921582 -0.45921582 0.0000000 ) +DEAL::2: ( 0.0000000 -0.21132487 0.0000000 ) +DEAL::3: ( 0.21132487 -0.21132487 0.0000000 ) +DEAL::4: ( 0.0000000 -0.45921582 0.45921582 ) +DEAL::5: ( 0.39433757 -0.39433757 0.39433757 ) +DEAL::6: ( 0.0000000 -0.21132487 0.21132487 ) +DEAL::7: ( 0.21132487 -0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.57735027 0.57735027 -0.57735027 ) +DEAL::1: ( -0.39433757 0.39433757 -0.39433757 ) +DEAL::2: ( 0.0000000 0.70710678 -0.70710678 ) +DEAL::3: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::4: ( -0.70710678 0.70710678 0.0000000 ) +DEAL::5: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::6: ( 0.0000000 1.0000000 0.0000000 ) +DEAL::7: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.39433757 0.39433757 -0.39433757 ) +DEAL::1: ( -0.21132487 0.21132487 -0.21132487 ) +DEAL::2: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::3: ( 0.0000000 0.21132487 -0.21132487 ) +DEAL::4: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::5: ( -0.21132487 0.21132487 0.0000000 ) +DEAL::6: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::7: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.70710678 -0.70710678 ) +DEAL::1: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::2: ( 0.57735027 0.57735027 -0.57735027 ) +DEAL::3: ( 0.39433757 0.39433757 -0.39433757 ) +DEAL::4: ( 0.0000000 1.0000000 0.0000000 ) +DEAL::5: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::6: ( 0.70710678 0.70710678 0.0000000 ) +DEAL::7: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 0.45921582 -0.45921582 ) +DEAL::1: ( 0.0000000 0.21132487 -0.21132487 ) +DEAL::2: ( 0.39433757 0.39433757 -0.39433757 ) +DEAL::3: ( 0.21132487 0.21132487 -0.21132487 ) +DEAL::4: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::5: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::6: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::7: ( 0.21132487 0.21132487 0.0000000 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.70710678 0.70710678 0.0000000 ) +DEAL::1: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::2: ( 0.0000000 1.0000000 0.0000000 ) +DEAL::3: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::4: ( -0.57735027 0.57735027 0.57735027 ) +DEAL::5: ( -0.39433757 0.39433757 0.39433757 ) +DEAL::6: ( 0.0000000 0.70710678 0.70710678 ) +DEAL::7: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( -0.45921582 0.45921582 0.0000000 ) +DEAL::1: ( -0.21132487 0.21132487 0.0000000 ) +DEAL::2: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::3: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::4: ( -0.39433757 0.39433757 0.39433757 ) +DEAL::5: ( -0.21132487 0.21132487 0.21132487 ) +DEAL::6: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::7: ( 0.0000000 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 0.0000000 1.0000000 0.0000000 ) +DEAL::1: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::2: ( 0.70710678 0.70710678 0.0000000 ) +DEAL::3: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::4: ( 0.0000000 0.70710678 0.70710678 ) +DEAL::5: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::6: ( 0.57735027 0.57735027 0.57735027 ) +DEAL::7: ( 0.39433757 0.39433757 0.39433757 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::Cell nodes: +DEAL::0: ( 2.7755576e-17 0.60566243 5.5511151e-17 ) +DEAL::1: ( 0.0000000 0.21132487 0.0000000 ) +DEAL::2: ( 0.45921582 0.45921582 0.0000000 ) +DEAL::3: ( 0.21132487 0.21132487 0.0000000 ) +DEAL::4: ( 0.0000000 0.45921582 0.45921582 ) +DEAL::5: ( 0.0000000 0.21132487 0.21132487 ) +DEAL::6: ( 0.39433757 0.39433757 0.39433757 ) +DEAL::7: ( 0.21132487 0.21132487 0.21132487 ) +DEAL::OK: cell bulk and boundary integrals match... + +DEAL::done... -- 2.39.5