From: rschulz Date: Mon, 15 May 2006 11:31:08 +0000 (+0000) Subject: changed the old point_difference and point_value to wrapper functions X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=152111a64c5f2a55d01be685be8b429d85a78301;p=dealii-svn.git changed the old point_difference and point_value to wrapper functions git-svn-id: https://svn.dealii.org/trunk@13112 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/grid/grid_tools.h b/deal.II/deal.II/include/grid/grid_tools.h index 768106a8b2..9b210e423e 100644 --- a/deal.II/deal.II/include/grid/grid_tools.h +++ b/deal.II/deal.II/include/grid/grid_tools.h @@ -16,7 +16,6 @@ #include #include -#include #include #include #include @@ -586,15 +585,6 @@ void GridTools::transform (const Predicate &predicate, } -template -inline -typename Container::active_cell_iterator -GridTools::find_active_cell_around_point (const Container &container, - const Point &p) -{ - return find_active_cell_around_point(StaticMappingQ1::mapping, container, p).first; -} - /*---------------------------- grid_tools.h ---------------------------*/ /* end of #ifndef __deal2__grid_tools_H */ diff --git a/deal.II/deal.II/include/numerics/vectors.h b/deal.II/deal.II/include/numerics/vectors.h index 267aa2f7ea..bbea9f130d 100644 --- a/deal.II/deal.II/include/numerics/vectors.h +++ b/deal.II/deal.II/include/numerics/vectors.h @@ -964,11 +964,9 @@ class VectorTools * components as the finite * element) at this point. * - * Since the function uses a - * simple test for checking - * whether a point is in a cell, - * it is only implemented for - * Q1-mapping at this time. + * This is a wrapper function + * using a Q1-mapping for cell + * boundaries. */ template static void point_difference (const DoFHandler& dof, @@ -988,11 +986,8 @@ class VectorTools * components as the finite * element) at this point. * - * This function additionally - * expects a mapping and uses - * a different algorithm to - * determine the cell surrounding - * the given point. + * This function uses an arbitrary + * mapping to evaluate the difference. */ template static void point_difference (const Mapping &mapping, @@ -1011,6 +1006,9 @@ class VectorTools * the (vector) value of this * function through the last * argument. + * + * This is a wrapper function using + * a Q1-mapping for cells. */ template static @@ -1027,6 +1025,9 @@ class VectorTools * vector at the given point, and * return the value of this * function. + * + * This is a wrapper function using + * a Q1-mapping for cells. */ template static @@ -1043,11 +1044,10 @@ class VectorTools * the given point, and return * the (vector) value of this * function through the last - * argument. This function expects - * additionally a mapping, and - * uses a different algorithm to - * determine the location of the - * point (see GridTools::find_active_cell_around_point). + * argument. + * + * This function uses an arbitrary + * mapping to evaluate the difference. */ template static @@ -1064,11 +1064,10 @@ class VectorTools * the given DoFHandler and nodal * vector at the given point, and * return the value of this - * function. This function expects - * additionally a mapping, and - * uses a different algorithm to - * determine the location of the - * point (see GridTools::find_active_cell_around_point). + * function. + * + * This function uses an arbitrary + * mapping to evaluate the difference. */ template static diff --git a/deal.II/deal.II/include/numerics/vectors.templates.h b/deal.II/deal.II/include/numerics/vectors.templates.h index f23a7c8f4d..fab9a727dd 100644 --- a/deal.II/deal.II/include/numerics/vectors.templates.h +++ b/deal.II/deal.II/include/numerics/vectors.templates.h @@ -1772,39 +1772,12 @@ VectorTools::point_difference (const DoFHandler &dof, Vector &difference, const Point &point) { - static const MappingQ1 mapping; - const FiniteElement& fe = dof.get_fe(); - - Assert(difference.size() == fe.n_components(), - ExcDimensionMismatch(difference.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 typename DoFHandler::active_cell_iterator - cell = GridTools::find_active_cell_around_point (dof, point); - - const Point unit_point - = mapping.transform_real_to_unit_cell(cell, point); - Assert (GeometryInfo::is_inside_unit_cell (unit_point), - ExcInternalError()); - - const Quadrature quadrature (unit_point); - FEValues fe_values(mapping, fe, quadrature, update_values); - fe_values.reinit(cell); - - // 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); - - if (fe.n_components() == 1) - difference(0) = exact_function.value(point); - else - exact_function.vector_value(point, difference); - - for (unsigned int i=0; i::mapping, + dof, + fe_function, + exact_function, + difference, + point); } @@ -1857,33 +1830,11 @@ VectorTools::point_value (const DoFHandler &dof, const Point &point, Vector &value) { - static const MappingQ1 mapping; - const FiniteElement& 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 typename DoFHandler::active_cell_iterator - cell = GridTools::find_active_cell_around_point (dof, point); - - const Point unit_point - = mapping.transform_real_to_unit_cell(cell, point); - Assert (GeometryInfo::is_inside_unit_cell (unit_point), - ExcInternalError()); - - const Quadrature quadrature (unit_point); - FEValues fe_values(mapping, fe, quadrature, update_values); - fe_values.reinit(cell); - - // 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]; + point_value(StaticMappingQ1::mapping, + dof, + fe_function, + point, + value); } @@ -1894,33 +1845,10 @@ VectorTools::point_value (const DoFHandler &dof, const InVector &fe_function, const Point &point) { - static const MappingQ1 mapping; - const FiniteElement& 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 typename DoFHandler::active_cell_iterator - cell = GridTools::find_active_cell_around_point (dof, point); - - const Point unit_point - = mapping.transform_real_to_unit_cell(cell, point); - Assert (GeometryInfo::is_inside_unit_cell (unit_point), - ExcInternalError()); - - const Quadrature quadrature (unit_point); - FEValues fe_values(mapping, fe, quadrature, update_values); - fe_values.reinit(cell); - - // 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 point_value(StaticMappingQ1::mapping, + dof, + fe_function, + point); } template diff --git a/deal.II/deal.II/source/grid/grid_tools.cc b/deal.II/deal.II/source/grid/grid_tools.cc index d375a5036c..8b44eb1df6 100644 --- a/deal.II/deal.II/source/grid/grid_tools.cc +++ b/deal.II/deal.II/source/grid/grid_tools.cc @@ -23,6 +23,7 @@ #include #include #include +#include #include #include @@ -554,6 +555,15 @@ GridTools::find_cells_adjacent_to_vertex(const Container &container, } +template +typename Container::active_cell_iterator +GridTools::find_active_cell_around_point (const Container &container, + const Point &p) +{ + return find_active_cell_around_point(StaticMappingQ1::mapping, container, p).first; +} + + template class Container> std::pair::active_cell_iterator, Point > GridTools::find_active_cell_around_point (const Mapping &mapping, diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 45b36b6298..d9fc3b0109 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -646,6 +646,15 @@ inconvenience this causes.

deal.II

    +
  1. + Changed: Functions VectorTools::point_value and VectorTools::point_difference using the old interface without + boundary mapping were replaced by wrapper functions calling the new + versions. +
    + (Ralf B. Schulz, 2006/05/15) +

  2. Changed: The old version of GridTools::find_active_cell_around_point has been replaced