From 4e5d9f52d5536f9c8c3e72b59c4fdf83ca4d7ce9 Mon Sep 17 00:00:00 2001 From: Nicola Giuliani Date: Thu, 23 Jan 2020 08:57:12 +0100 Subject: [PATCH] added description of integration of distorted cell and modified test --- source/grid/tria_accessor.cc | 7 ++++++- .../measure_of_distorted_codim_one_face.cc | 14 ++++---------- 2 files changed, 10 insertions(+), 11 deletions(-) diff --git a/source/grid/tria_accessor.cc b/source/grid/tria_accessor.cc index 09015c006a..a26edc3eea 100644 --- a/source/grid/tria_accessor.cc +++ b/source/grid/tria_accessor.cc @@ -1312,7 +1312,12 @@ namespace ((v03 * v03) * (v01 * v01) * (v02 * v02))) >= 1e-24) { // If the vectors are non planar we integrate the norm of the normal - // vector using a numerical Gauss scheme of order 4. + // vector using a numerical Gauss scheme of order 4. In particular we + // consider a bilinear quad x(u,v) = (1-v)((1-u)v_0 + u v_1) + + // v((1-u)v_2 + u v_3), consequently we compute the normal vector as + // n(u,v) = t_u x t_v = w_1 + u w_2 + v w_3. The integrand function is + // || n(u,v) || = sqrt(a + b u^2 + c v^2 + d u + e v + f uv). + // We integrate it using a QGauss<2> (4) computed explicitly. const Tensor<1, 3> w_1 = cross_product_3d(accessor.vertex(1) - accessor.vertex(0), accessor.vertex(2) - accessor.vertex(0)); diff --git a/tests/codim_one/measure_of_distorted_codim_one_face.cc b/tests/codim_one/measure_of_distorted_codim_one_face.cc index 342a64f100..9ae5eedc3a 100644 --- a/tests/codim_one/measure_of_distorted_codim_one_face.cc +++ b/tests/codim_one/measure_of_distorted_codim_one_face.cc @@ -43,11 +43,7 @@ test() std::vector> cells; SubCellData subcelldata; - double tol = 1e-12; - // vertices.push_back(Point<3>{10, 56, 0}); - // vertices.push_back(Point<3>{22, 1, 0}); - // vertices.push_back(Point<3>{15, 44, 0}); - // vertices.push_back(Point<3>{1, 1, 1}); + const double tol = 1e-12; vertices.push_back(Point<3>{0, 0, 1}); vertices.push_back(Point<3>{1, 0, -10}); vertices.push_back(Point<3>{0, 1, -1}); @@ -70,12 +66,10 @@ test() MappingQ1<2, 3> mapping; FEValues<2, 3> fe_values(mapping, fe, quadrature_formula, update_JxW_values); - double sum_1 = 0; - double sum_2 = 0; - auto cell = dof_handler.begin_active(); - auto accessor = tria.begin_active(); + double sum_2 = 0; + auto cell = dof_handler.begin_active(); - sum_1 = cell->measure(); + const double sum_1 = cell->measure(); fe_values.reinit(cell); for (unsigned int q = 0; q < quadrature_formula.size(); ++q) { -- 2.39.5