From: Timo Heister Date: Thu, 21 Jan 2016 16:01:19 +0000 (-0500) Subject: fix 1d Kelly estimator X-Git-Tag: v8.4.0-rc2~66^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=94ddae16cd6295abf6f1979fbbc705a2cf1519aa;p=dealii.git fix 1d Kelly estimator --- diff --git a/source/numerics/error_estimator_1d.cc b/source/numerics/error_estimator_1d.cc index 1868c6d78e..675797954a 100644 --- a/source/numerics/error_estimator_1d.cc +++ b/source/numerics/error_estimator_1d.cc @@ -319,6 +319,8 @@ estimate (const Mapping<1,spacedim> &mapping, const QTrapez<1> quadrature; const hp::QCollection<1> q_collection(quadrature); + const QGauss<0> face_quadrature(1); + const hp::QCollection<0> q_face_collection(face_quadrature); const hp::FECollection<1,spacedim> fe (dof_handler.get_fe()); @@ -327,6 +329,8 @@ estimate (const Mapping<1,spacedim> &mapping, hp::FEValues<1,spacedim> fe_values (mapping_collection, fe, q_collection, update_gradients); + hp::FEFaceValues<1,spacedim> fe_face_values (/*mapping_collection,*/ fe, q_face_collection, + update_normal_vectors); // loop over all cells and do something on the cells which we're told to // work on. note that the error indicator is only a sum over the two @@ -345,6 +349,11 @@ estimate (const Mapping<1,spacedim> &mapping, for (unsigned int n=0; nactive_cell_index()) = 0; + fe_values.reinit (cell); + for (unsigned int s=0; s &mapping, while (neighbor->has_children()) neighbor = neighbor->child(n==0 ? 1 : 0); - // now get the gradients on the both sides of the point - fe_values.reinit (cell); - - for (unsigned int s=0; s normal = + fe_face_values.get_present_fe_values().get_all_normal_vectors()[0]; if (neighbor.state() == IteratorState::valid) { @@ -371,12 +377,15 @@ estimate (const Mapping<1,spacedim> &mapping, .get_function_gradients (*solutions[s], gradients_neighbor[s]); - // extract the gradients of all the components. [0] means: - // x-derivative, which is the only one here + fe_face_values.reinit (neighbor, n==0 ? 1 : 0); + Tensor<1,spacedim> neighbor_normal = + fe_face_values.get_present_fe_values().get_all_normal_vectors()[0]; + + // extract the gradient in normal direction of all the components. for (unsigned int s=0; s &mapping, for (unsigned int component=0; component