]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
changed the old point_difference and point_value to wrapper functions
authorrschulz <rschulz@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 May 2006 11:31:08 +0000 (11:31 +0000)
committerrschulz <rschulz@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 May 2006 11:31:08 +0000 (11:31 +0000)
git-svn-id: https://svn.dealii.org/trunk@13112 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_tools.h
deal.II/deal.II/include/numerics/vectors.h
deal.II/deal.II/include/numerics/vectors.templates.h
deal.II/deal.II/source/grid/grid_tools.cc
deal.II/doc/news/changes.html

index 768106a8b2088ad99ee08cc7c9cdd3b60260ed91..9b210e423e25f5e53901df66db601cf6c8b98a50 100644 (file)
@@ -16,7 +16,6 @@
 
 #include <base/config.h>
 #include <fe/mapping.h>
-#include <fe/mapping_q1.h>
 #include <grid/tria.h>
 #include <grid/tria_accessor.h>
 #include <grid/tria_iterator.h>
@@ -586,15 +585,6 @@ void GridTools::transform (const Predicate    &predicate,
 }
 
 
-template <int dim, typename Container>
-inline
-typename Container::active_cell_iterator
-GridTools::find_active_cell_around_point (const Container  &container,
-                                          const Point<dim> &p)
-{
-   return find_active_cell_around_point(StaticMappingQ1<dim>::mapping, container, p).first;
-}
-
 
 /*----------------------------   grid_tools.h     ---------------------------*/
 /* end of #ifndef __deal2__grid_tools_H */
index 267aa2f7eab6260a952fd0abf0c5952f20a8867b..bbea9f130db7777ee4ce98cc6ff35357fe56ad2d 100644 (file)
@@ -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 <int dim, class InVector>
     static void point_difference (const DoFHandler<dim>& 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 <int dim, class InVector>
     static void point_difference (const Mapping<dim>    &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 <int dim, class InVector>
     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 <int dim, class InVector>
     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 <int dim, class InVector>
     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 <int dim, class InVector>
     static
index f23a7c8f4d65bf31cf86274b581a6200cb216997..fab9a727dd2d04f64b97ef40513658447f77221f 100644 (file)
@@ -1772,39 +1772,12 @@ VectorTools::point_difference (const DoFHandler<dim> &dof,
                               Vector<double>        &difference,
                               const Point<dim>      &point)
 {
-  static const MappingQ1<dim> mapping;
-  const FiniteElement<dim>& 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<dim>::active_cell_iterator
-    cell = GridTools::find_active_cell_around_point (dof, point);
-
-  const Point<dim> unit_point
-    = mapping.transform_real_to_unit_cell(cell, point);
-  Assert (GeometryInfo<dim>::is_inside_unit_cell (unit_point),
-          ExcInternalError());
-
-  const Quadrature<dim> quadrature (unit_point);
-  FEValues<dim> 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<Vector<double> > u_value(1, Vector<double> (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<difference.size(); ++i)
-    difference(i) -= u_value[0](i);
+   point_difference(StaticMappingQ1<dim>::mapping,
+                    dof,
+                    fe_function,
+                    exact_function,
+                    difference,
+                    point);
 }
 
 
@@ -1857,33 +1830,11 @@ VectorTools::point_value (const DoFHandler<dim> &dof,
                          const Point<dim>      &point,
                          Vector<double>        &value)
 {
-  static const MappingQ1<dim> mapping;
-  const FiniteElement<dim>& 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<dim>::active_cell_iterator
-    cell = GridTools::find_active_cell_around_point (dof, point);
-
-  const Point<dim> unit_point
-    = mapping.transform_real_to_unit_cell(cell, point);
-  Assert (GeometryInfo<dim>::is_inside_unit_cell (unit_point),
-          ExcInternalError());
-
-  const Quadrature<dim> quadrature (unit_point);
-  FEValues<dim> 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<Vector<double> > u_value(1, Vector<double> (fe.n_components()));
-  fe_values.get_function_values(fe_function, u_value);
-
-  value = u_value[0];
+   point_value(StaticMappingQ1<dim>::mapping,
+               dof,
+               fe_function,
+               point,
+               value);
 }
 
 
@@ -1894,33 +1845,10 @@ VectorTools::point_value (const DoFHandler<dim> &dof,
                          const InVector        &fe_function,
                          const Point<dim>      &point)
 {
-  static const MappingQ1<dim> mapping;
-  const FiniteElement<dim>& 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<dim>::active_cell_iterator
-    cell = GridTools::find_active_cell_around_point (dof, point);
-
-  const Point<dim> unit_point
-    = mapping.transform_real_to_unit_cell(cell, point);
-  Assert (GeometryInfo<dim>::is_inside_unit_cell (unit_point),
-          ExcInternalError());
-
-  const Quadrature<dim> quadrature (unit_point);
-  FEValues<dim> 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<double> u_value(1);
-  fe_values.get_function_values(fe_function, u_value);
-
-  return u_value[0];
+   return point_value(StaticMappingQ1<dim>::mapping,
+                      dof,
+                      fe_function,
+                      point);
 }
 
 template <int dim, class InVector>
index d375a5036ca6402f56150fb80ea6112f584f8599..8b44eb1df61bd6c84b284ec9dea9e0cd1af792c4 100644 (file)
@@ -23,6 +23,7 @@
 #include <dofs/dof_accessor.h>
 #include <dofs/dof_tools.h>
 #include <fe/fe_dgq.h>
+#include <fe/mapping_q1.h>
 #include <multigrid/mg_dof_handler.h>
 #include <multigrid/mg_dof_accessor.h>
 
@@ -554,6 +555,15 @@ GridTools::find_cells_adjacent_to_vertex(const Container<dim> &container,
 }  
 
 
+template <int dim, typename Container>
+typename Container::active_cell_iterator
+GridTools::find_active_cell_around_point (const Container  &container,
+                                          const Point<dim> &p)
+{
+   return find_active_cell_around_point(StaticMappingQ1<dim>::mapping, container, p).first;
+}
+
+
 template <int dim, template <int> class Container>
 std::pair<typename Container<dim>::active_cell_iterator, Point<dim> >
 GridTools::find_active_cell_around_point (const Mapping<dim>   &mapping,
index 45b36b6298caeef7a2cd6a8979afead56a965bd5..d9fc3b0109098a9055999626bfa9a55b2a55799d 100644 (file)
@@ -646,6 +646,15 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       Changed: Functions <code class="class">VectorTools</code>::<code
+       class="member">point_value</code> and <code class="class">VectorTools</code>::<code
+       class="member">point_difference</code> using the old interface without
+       boundary mapping were replaced by wrapper functions calling the new
+       versions.
+       <br> 
+       (Ralf B. Schulz, 2006/05/15)
+       </p>
   <li> <p>
        Changed: The old version of <code class="class">GridTools</code>::<code
        class="member">find_active_cell_around_point</code> has been replaced

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.