From 8827cdc40b8e61cce54a75038c67b3ba9a259c22 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Wed, 20 Jan 2016 16:28:24 -0500 Subject: [PATCH] rewrite normal vector computation --- .../numerics/error_estimator.templates.h | 33 ++++++++++++++----- 1 file changed, 25 insertions(+), 8 deletions(-) diff --git a/include/deal.II/numerics/error_estimator.templates.h b/include/deal.II/numerics/error_estimator.templates.h index 483cd2c5f0..82771dec1f 100644 --- a/include/deal.II/numerics/error_estimator.templates.h +++ b/include/deal.II/numerics/error_estimator.templates.h @@ -151,6 +151,11 @@ namespace internal */ std::vector > normal_vectors; + /** + * Normal vectors of the opposing face. + */ + std::vector > neighbor_normal_vectors; + /** * Two arrays needed for the values of coefficients in the jumps, if * they are given. @@ -254,6 +259,7 @@ namespace internal (face_quadratures.max_n_quadrature_points(), std::vector > (fe.n_components()))), normal_vectors (face_quadratures.max_n_quadrature_points()), + neighbor_normal_vectors (face_quadratures.max_n_quadrature_points()), coefficient_values1 (face_quadratures.max_n_quadrature_points()), coefficient_values (face_quadratures.max_n_quadrature_points(), dealii::Vector (fe.n_components())), @@ -275,6 +281,7 @@ namespace internal const unsigned int n_components = finite_element.n_components(); normal_vectors.resize(n_q_points); + neighbor_normal_vectors.resize(n_q_points); coefficient_values1.resize(n_q_points); coefficient_values.resize(n_q_points); JxW_values.resize(n_q_points); @@ -342,8 +349,7 @@ namespace internal integrate_over_face (ParallelData ¶llel_data, const typename DoFHandlerType::face_iterator &face, - dealii::hp::FEFaceValues &fe_face_values_cell, - const std::vector > & other_normals) + dealii::hp::FEFaceValues &fe_face_values_cell) { const unsigned int n_q_points = parallel_data.psi[0].size(), n_components = parallel_data.finite_element.n_components(), @@ -378,16 +384,19 @@ namespace internal parallel_data.phi[n][point][component] = (parallel_data.psi[n][point][component] * parallel_data.normal_vectors[point]); + if (face->at_boundary() == false) { // compute the jump in the gradients + for (unsigned int n=0; nface(neighbor_neighbor)] = - integrate_over_face (parallel_data, face, fe_face_values, - fe_subface_values.get_present_fe_values().get_all_normal_vectors()); + integrate_over_face (parallel_data, face, fe_face_values); for (unsigned int i = 0; i < local_face_integrals[neighbor_child->face(neighbor_neighbor)].size(); i++) local_face_integrals[neighbor_child->face(neighbor_neighbor)][i] *= factor; } -- 2.39.5