]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
point evaluation implemented, not tested
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 14 Nov 2002 07:48:38 +0000 (07:48 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 14 Nov 2002 07:48:38 +0000 (07:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@6762 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/vectors.h
deal.II/deal.II/source/numerics/vectors.cc

index 160619dc763e30279cf5f3584a53d321d29fcae5..be6ec9107fe0a436e8730eb6b4f8bfcf9265135b 100644 (file)
@@ -21,6 +21,7 @@
 #include <vector>
 #include <set>
 
+template <int dim> class Point;
 template <int dim> class Function;
 template <int dim> class FunctionMap;
 template <int dim> class Quadrature;
@@ -820,6 +821,27 @@ class VectorTools
                                      const Function<dim>   *weight=0,
                                      const double exponent = 2.);
 
+                                    /**
+                                     * Point error evaluation. Find
+                                     * the first cell containing the
+                                     * given point and compute the
+                                     * difference of a finite element
+                                     * function and a continuous
+                                     * function 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 yet.
+                                     */
+    template <int dim, class InVector>
+    static void point_difference (const DoFHandler<dim>& dof,
+                                 const InVector&        fe_function,
+                                 const Function<dim>&   exact_solution,
+                                 Vector<double>&        difference,
+                                 const Point<dim>&      point);
+    
                                     /**
                                      * Mean-value filter for Stokes.
                                      * The pressure in Stokes'
index b0f9fa7f5cc93314005d728ee0d74990207bca7f..d806af77686d1e3394771d7692be268fa6d29d24 100644 (file)
@@ -1358,6 +1358,55 @@ VectorTools::integrate_difference (const DoFHandler<dim>    &dof,
 
 
 
+template <int dim, class InVector>
+void
+VectorTools::point_difference (const DoFHandler<dim>& dof,
+                              const InVector&        fe_function,
+                              const Function<dim>&   exact_function,
+                              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()));
+
+  DoFHandler<dim>::active_cell_iterator cell = dof_handler->begin_active();
+  const DoFHandler<dim>::active_cell_iterator endc = dof_handler->end();
+
+  std::vector<Vector<double> > u_value(1, Vector<double> (fe.n_components()));
+  
+  for (;cell != endc; ++cell)
+    {
+      if (cell->point_inside(point))
+       {
+         Point<dim> unit_point = mapping->transform_real_to_unit_cell(
+            cell, evaluation_point);
+          for (unsigned int i=0; i<dim; ++i)
+            Assert(0<=unit_point(i) && unit_point(i)<=1, ExcInternalError());
+          vector<Point<dim> > quadrature_point(1, unit_point);
+          vector<double> weight(1, 1.);
+          Quadrature<dim> quadrature(quadrature_point, weight);
+          FEValues<dim> fe_values(
+            mapping, fe, quadrature, UpdateFlags(update_values));
+          fe_values.reinit(cell);
+          fe_values.get_function_values(u, u_value);
+         break;
+       }
+    }
+  if (fe.n_components() == 1)
+    {
+      const double value = exact_function.value(point);
+      difference(0) = value;
+    } else {
+      exact_function.vector_value(point, difference);
+    }
+  for (unsigned int i=0;i<difference.size();++i)
+    difference(i) -= u_value[0](i);
+}
+
+
 
 template <int dim, class InVector>
 double

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.