]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make sure we don't actually call a function that isn't meant to be
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 14 Aug 2010 03:09:47 +0000 (03:09 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 14 Aug 2010 03:09:47 +0000 (03:09 +0000)
called for this particular space dimension. We use the idiom:

  void foo (Quadrature<3> &q)
  {...}

  template <int dim>
  void foo (Quadrature<dim> &q)
  { Assert (false, ExcInternalError()); }

for functions that are supposed to be called only in 3d.

git-svn-id: https://svn.dealii.org/trunk@21656 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/vectors.templates.h

index 033ff60afa5e110de8652762f6c8282adceab79a..e40974fb8f7de5f569a139c30556521495992626 100644 (file)
@@ -2930,17 +2930,19 @@ namespace internals
                                     // This function computes the
                                     // projection of the boundary
                                     // function on edges for 3D.
-    template<int dim, typename cell_iterator>
+    template<typename cell_iterator>
     void
     compute_edge_projection (const cell_iterator& cell,
                             const unsigned int face,
                             const unsigned int line,
-                            FEValues<dim>& fe_values,
-                            const Quadrature<dim>& quadrature,
-                            const Function<dim>& boundary_function,
+                            FEValues<3>& fe_values,
+                            const Quadrature<3>& quadrature,
+                            const Function<3>& boundary_function,
                             const unsigned int first_vector_component,
                             std::vector<double>& dof_values)
     {
+      const unsigned int dim = 3;
+
       fe_values.reinit (cell);
 
                                       // Initialize the required
@@ -3149,20 +3151,41 @@ namespace internals
     }
 
 
+                                    // dummy implementation of above
+                                    // function for all other
+                                    // dimensions
+    template<int dim, typename cell_iterator>
+    void
+    compute_edge_projection (const cell_iterator&,
+                            const unsigned int,
+                            const unsigned int,
+                            FEValues<dim>&,
+                            const Quadrature<dim>&,
+                            const Function<dim>&,
+                            const unsigned int,
+                            std::vector<double>&)
+    {
+      Assert (false, ExcInternalError());
+    }
+
+
+
 
                                     // This function computes the
                                     // projection of the boundary
                                     // function on the interior of
                                     // faces in 3D.
-    template<int dim, typename cell_iterator>
+    template <typename cell_iterator>
     void
     compute_face_projection (const cell_iterator& cell,
                             const unsigned int face,
-                            FEValues<dim>& fe_values,
-                            const Function<dim>& boundary_function,
+                            FEValues<3>& fe_values,
+                            const Function<3>& boundary_function,
                             const unsigned int first_vector_component,
                             std::vector<double>& dof_values)
     {
+      const unsigned dim = 3;
+
       fe_values.reinit (cell);
 
                                       // Initialize the required objects.
@@ -3397,20 +3420,36 @@ namespace internals
 
 
 
+                                    // dummy implementation of above
+                                    // function for dim != 3
+    template<int dim, typename cell_iterator>
+    void
+    compute_face_projection (const cell_iterator&,
+                            const unsigned int,
+                            FEValues<dim>&,
+                            const Function<dim>&,
+                            const unsigned int,
+                            std::vector<double>&)
+    {
+      Assert (false, ExcInternalError ());
+    }
+
 
                                     // This function computes the
                                     // projection of the boundary
                                     // function on the faces in 2D.
-    template<int dim, typename cell_iterator>
+    template<typename cell_iterator>
     void
     compute_face_projection (const cell_iterator& cell,
                             const unsigned int face,
-                            FEValues<dim>& fe_values,
-                            const Quadrature<dim>& quadrature,
-                            const Function<dim>& boundary_function,
+                            FEValues<2>& fe_values,
+                            const Quadrature<2>& quadrature,
+                            const Function<2>& boundary_function,
                             const unsigned int first_vector_component,
                             std::vector<double>& dof_values)
     {
+      const unsigned int dim = 2;
+
       fe_values.reinit (cell);
 
                                       // Initialize the required objects.
@@ -3570,6 +3609,22 @@ namespace internals
            dof_values[i + 1] = solution (i);
        }
     }
+
+
+                                    // dummy implementation of above
+                                    // function for dim != 2
+    template<int dim, typename cell_iterator>
+    void
+    compute_face_projection (const cell_iterator&,
+                            const unsigned int,
+                            FEValues<dim>&,
+                            const Quadrature<dim>& quadrature,
+                            const Function<dim>&,
+                            const unsigned int,
+                            std::vector<double>&)
+    {
+      Assert (false, ExcInternalError ());
+    }
   }
 }
 

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.