]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add the weight function argument also in the scalar #integrate_difference# function
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 5 Jul 1999 11:19:10 +0000 (11:19 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 5 Jul 1999 11:19:10 +0000 (11:19 +0000)
git-svn-id: https://svn.dealii.org/trunk@1533 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 95fda10a0ec70bf999a7288fcdc38df70d884dcd..2e4a8f1e375154b56d309c214b40c649c2bd5b10 100644 (file)
@@ -434,20 +434,6 @@ class VectorTools //<dim>
                                      * accuracy of the #double# data type is
                                      * used.
                                      *
-                                     * See the general documentation of this
-                                     * class for more information.
-                                     */
-    static void integrate_difference (const DoFHandler<dim>    &dof,
-                                     const Vector<double>     &fe_function,
-                                     const Function<dim>      &exact_solution,
-                                     Vector<float>            &difference,
-                                     const Quadrature<dim>    &q,
-                                     const NormType           &norm);
-
-                                    /**
-                                     * Compute the error for the solution of a system.
-                                     * See the other #integrate_difference#.
-                                     *
                                      * The additional argument #weight# allows
                                      * to evaluate weighted norms. This is useful
                                      * for weighting the error of different parts
@@ -460,6 +446,21 @@ class VectorTools //<dim>
                                      * solutions.
                                      * By default, no weighting function is given,
                                      * i.e. weight=1 in the whole domain.
+                                     *
+                                     * See the general documentation of this
+                                     * class for more information.
+                                     */
+    static void integrate_difference (const DoFHandler<dim>    &dof,
+                                     const Vector<double>     &fe_function,
+                                     const Function<dim>      &exact_solution,
+                                     Vector<float>            &difference,
+                                     const Quadrature<dim>    &q,
+                                     const NormType           &norm,
+                                     const Function<dim>      *weight=0);
+
+                                    /**
+                                     * Compute the error for the solution of a system.
+                                     * See the other #integrate_difference#.
                                      */
     static void integrate_difference (const DoFHandler<dim>   &dof,
                                      const Vector<double>     &fe_function,
index 2dc3a163e0e26167d99442da310295d6177a1fc6..5bcdaeb20e3172fd9cb9fb35b9b9cdbc6fe517c3 100644 (file)
@@ -693,7 +693,9 @@ void VectorTools<dim>::integrate_difference (const DoFHandler<dim>    &dof,
                                             const Function<dim>      &exact_solution,
                                             Vector<float>            &difference,
                                             const Quadrature<dim>    &q,
-                                            const NormType           &norm) {
+                                            const NormType           &norm,
+                                            const Function<dim>      *weight=0)
+{
   const FiniteElement<dim> &fe = dof.get_fe();
     
   difference.reinit (dof.get_tria().n_active_cells());
@@ -782,6 +784,17 @@ void VectorTools<dim>::integrate_difference (const DoFHandler<dim>    &dof,
                      Assert (false, ExcNotImplemented());
              };
 
+                                            // now weight the values
+                                            // at the quadrature points due
+                                            // to the weighting function
+           if (weight)
+             {
+               vector<double> w(n_q_points);
+               weight->value_list(fe_values.get_quadrature_points(),w);
+               for (unsigned int q=0; q<n_q_points; ++q)
+                 psi[q]*=w[q];
+             }
+
                                             // ok, now we have the integrand,
                                             // let's compute the integral,
                                             // which is
@@ -852,6 +865,17 @@ void VectorTools<dim>::integrate_difference (const DoFHandler<dim>    &dof,
            for (unsigned int i=0; i<n_q_points; ++i)
              psi_square[i] = sqr_point(psi[i]);
 
+                                            // now weight the values
+                                            // at the quadrature points due
+                                            // to the weighting function
+           if (weight)
+             {
+               vector<double> w(n_q_points);
+               weight->value_list(fe_values.get_quadrature_points(),w);
+               for (unsigned int q=0; q<n_q_points; ++q)
+                 psi_square[q]*=w[q];
+             }
+
                                             // add seminorm to L_2 norm or
                                             // to zero
            diff += inner_product (psi_square.begin(), psi_square.end(),

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.