From 267500f57fa698c8eca90a44dccfe3fc13d847b3 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Mon, 4 May 2015 00:42:30 +0200 Subject: [PATCH] Add an entry to changes.h and reduce the number of actual implementations for point_value and point_gradient --- doc/news/changes.h | 6 ++ include/deal.II/numerics/vector_tools.h | 2 +- .../deal.II/numerics/vector_tools.templates.h | 98 +++---------------- 3 files changed, 19 insertions(+), 87 deletions(-) diff --git a/doc/news/changes.h b/doc/news/changes.h index 934a491e02..822adabc74 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -472,6 +472,12 @@ inconvenience this causes.

Specific improvements

    +
  1. New: The function VectorTools::point_gradient has been added to compute + the gradient of a given FE function. +
    + (Daniel Arndt, 2015/05/03) +
  2. +
  3. New: dealii:Vector, dealii::BlockVector, TrilinosWrappers::MPI::Vector, TrilinosWrappers::MPI::BlockVector, PETScWrappers::MPI::Vector and PETScWrappers::MPI::BlockVector now have diff --git a/include/deal.II/numerics/vector_tools.h b/include/deal.II/numerics/vector_tools.h index 7ce3c469e6..2da0dec20c 100644 --- a/include/deal.II/numerics/vector_tools.h +++ b/include/deal.II/numerics/vector_tools.h @@ -2112,7 +2112,7 @@ namespace VectorTools * (vector) gradient of this function through the last argument. * * This is a wrapper function using a Q1-mapping for cell boundaries to call - * the other point_difference() function. + * the other point_gradient() function. * * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index 08ecab5eff..3ec9ff28b3 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -6729,28 +6729,10 @@ namespace VectorTools Assert(fe.n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); - // first find the cell in which this point - // is, initialize a quadrature rule with - // it, and then a FEValues object - const std::pair::active_cell_iterator, Point > - cell_point = GridTools::find_active_cell_around_point (mapping, dof, point); - - AssertThrow(cell_point.first->is_locally_owned(), - ExcPointNotAvailableHere()); - Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < 1e-10, - ExcInternalError()); + Vector value(1); + point_value(mapping, dof, fe_function, point, value); - const Quadrature - quadrature (GeometryInfo::project_to_unit_cell(cell_point.second)); - FEValues fe_values(mapping, fe, quadrature, update_values); - fe_values.reinit(cell_point.first); - - // then use this to get at the values of - // the given fe_function at this point - std::vector u_value(1); - fe_values.get_function_values(fe_function, u_value); - - return u_value[0]; + return value(0); } @@ -6766,29 +6748,10 @@ namespace VectorTools Assert(fe.n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); - // first find the cell in which this point - // is, initialize a quadrature rule with - // it, and then a FEValues object - const std::pair::active_cell_iterator, Point > - cell_point = GridTools::find_active_cell_around_point (mapping, dof, point); - - AssertThrow(cell_point.first->is_locally_owned(), - ExcPointNotAvailableHere()); - Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < 1e-10, - ExcInternalError()); - - const Quadrature - quadrature (GeometryInfo::project_to_unit_cell(cell_point.second)); - hp::FEValues hp_fe_values(mapping, fe, hp::QCollection(quadrature), update_values); - hp_fe_values.reinit(cell_point.first); - const FEValues &fe_values = hp_fe_values.get_present_fe_values(); - - // then use this to get at the values of - // the given fe_function at this point - std::vector u_value(1); - fe_values.get_function_values(fe_function, u_value); + Vector value(1); + point_value(mapping, dof, fe_function, point, value); - return u_value[0]; + return value(0); } @@ -6941,28 +6904,10 @@ namespace VectorTools Assert(fe.n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); - // first find the cell in which this point - // is, initialize a quadrature rule with - // it, and then a FEValues object - const std::pair::active_cell_iterator, Point > - cell_point = GridTools::find_active_cell_around_point (mapping, dof, point); - - AssertThrow(cell_point.first->is_locally_owned(), - ExcPointNotAvailableHere()); - Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < 1e-10, - ExcInternalError()); + std::vector > gradient(1); + point_gradient (mapping, dof, fe_function, point, gradient); - const Quadrature - quadrature (GeometryInfo::project_to_unit_cell(cell_point.second)); - FEValues fe_values(mapping, fe, quadrature, update_gradients); - fe_values.reinit(cell_point.first); - - // then use this to get the gradients of - // the given fe_function at this point - std::vector > u_gradient(1); - fe_values.get_function_gradients(fe_function, u_gradient); - - return u_gradient[0]; + return gradient[0]; } @@ -6978,29 +6923,10 @@ namespace VectorTools Assert(fe.n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); - // first find the cell in which this point - // is, initialize a quadrature rule with - // it, and then a FEValues object - const std::pair::active_cell_iterator, Point > - cell_point = GridTools::find_active_cell_around_point (mapping, dof, point); - - AssertThrow(cell_point.first->is_locally_owned(), - ExcPointNotAvailableHere()); - Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < 1e-10, - ExcInternalError()); - - const Quadrature - quadrature (GeometryInfo::project_to_unit_cell(cell_point.second)); - hp::FEValues hp_fe_values(mapping, fe, hp::QCollection(quadrature), update_gradients); - hp_fe_values.reinit(cell_point.first); - const FEValues &fe_values = hp_fe_values.get_present_fe_values(); - - // then use this to get the gradients of - // the given fe_function at this point - std::vector > u_gradient(1); - fe_values.get_function_gradients(fe_function, u_gradient); + std::vector > gradient(1); + point_gradient (mapping, dof, fe_function, point, gradient); - return u_gradient[0]; + return gradient[0]; } -- 2.39.5