From e3f62a2567b713a6de97fe72fd29437e5b810c41 Mon Sep 17 00:00:00 2001 From: goll Date: Thu, 26 Jan 2012 09:47:37 +0000 Subject: [PATCH] point_value now also works with hp::dofhandler. git-svn-id: https://svn.dealii.org/trunk@24932 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/numerics/vectors.h | 45 ++++++++ .../deal.II/numerics/vectors.templates.h | 100 ++++++++++++++++++ 2 files changed, 145 insertions(+) diff --git a/deal.II/include/deal.II/numerics/vectors.h b/deal.II/include/deal.II/numerics/vectors.h index 38b8a6cb89..0bf46f634b 100644 --- a/deal.II/include/deal.II/numerics/vectors.h +++ b/deal.II/include/deal.II/numerics/vectors.h @@ -1874,6 +1874,18 @@ namespace VectorTools const Point &point, Vector &value); + /** + * Same as above for hp. + */ + + template + void + point_value (const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point, + Vector &value); + + /** * Evaluate a scalar finite * element function defined by @@ -1898,6 +1910,16 @@ namespace VectorTools const InVector &fe_function, const Point &point); + /** + * Same as above for hp. + */ + + template + double + point_value (const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point); + /** * Evaluate a possibly * vector-valued finite element @@ -1921,6 +1943,18 @@ namespace VectorTools const Point &point, Vector &value); + /** + * Same as above for hp. + */ + + template + void + point_value (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point, + Vector &value); + /** * Evaluate a scalar finite * element function defined by @@ -1941,6 +1975,17 @@ namespace VectorTools const InVector &fe_function, const Point &point); + /** + * Same as above for hp. + */ + + template + double + point_value (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point); + //@} /** * Mean value operations diff --git a/deal.II/include/deal.II/numerics/vectors.templates.h b/deal.II/include/deal.II/numerics/vectors.templates.h index ff66a2d8c3..578c21b52b 100644 --- a/deal.II/include/deal.II/numerics/vectors.templates.h +++ b/deal.II/include/deal.II/numerics/vectors.templates.h @@ -5359,6 +5359,20 @@ namespace VectorTools } + template + void + point_value (const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point, + Vector &value) + { + point_value(hp::StaticMappingQ1::mapping_collection, + dof, + fe_function, + point, + value); + } + template double @@ -5372,6 +5386,20 @@ namespace VectorTools point); } + + template + double + point_value (const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point) + { + return point_value(hp::StaticMappingQ1::mapping_collection, + dof, + fe_function, + point); + } + + template void point_value (const Mapping &mapping, @@ -5410,6 +5438,42 @@ namespace VectorTools } + template + void + point_value (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point, + Vector &value) + { + const hp::FECollection& fe = dof.get_fe(); + + Assert(value.size() == fe.n_components(), + ExcDimensionMismatch(value.size(), fe.n_components())); + + // 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); + + 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, Vector (fe.n_components())); + fe_values.get_function_values(fe_function, u_value); + + value = u_value[0]; + } + template double @@ -5446,6 +5510,42 @@ namespace VectorTools } + template + double + point_value (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const InVector &fe_function, + const Point &point) + { + const hp::FECollection& fe = dof.get_fe(); + + 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); + + 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); + + return u_value[0]; + } + + template double compute_mean_value (const Mapping &mapping, -- 2.39.5