]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
more cleaning of FE_RT
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Jan 2009 21:21:16 +0000 (21:21 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Jan 2009 21:21:16 +0000 (21:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@18166 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_raviart_thomas.h
deal.II/deal.II/source/fe/fe_raviart_thomas.cc
deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc

index fde8131c31288db2c87d2b85b1c9a7c276ae2f89..28f609f87ed883a4c4edd12265cb0ad9315fc321 100644 (file)
@@ -387,6 +387,17 @@ class FE_RaviartThomasNodal
                                      */
     static std::vector<bool>
     get_ria_vector (const unsigned int degree);
+                                    /**
+                                     * Check whether a shape function
+                                     * may be non-zero on a face.
+                                     *
+                                     * Right now, this is only
+                                     * implemented for RT0 in
+                                     * 1D. Otherwise, returns always
+                                     * @p true.
+                                     */
+    virtual bool has_support_on_face (const unsigned int shape_index,
+                                     const unsigned int face_index) const;
                                     /**
                                      * Initialize the
                                      * FiniteElement<dim>::generalized_support_points
index cf158ff4b2460b178ea57ed75ea9c445e7470c20..4ab6e9f1b0408a311a1a068fd0df49be493d7940 100644 (file)
@@ -416,10 +416,12 @@ FE_RaviartThomas<dim>::initialize_restriction()
 
 template <>
 std::vector<unsigned int>
-FE_RaviartThomas<1>::get_dpo_vector (const unsigned int)
+FE_RaviartThomas<1>::get_dpo_vector (const unsigned int deg)
 {
-  Assert (false, ExcImpossibleInDim(1));
-  return std::vector<unsigned int>();
+  std::vector<unsigned int> dpo(2);
+  dpo[0] = 1;
+  dpo[1] = deg;
+  return dpo;
 }
 
 #endif
@@ -427,20 +429,17 @@ FE_RaviartThomas<1>::get_dpo_vector (const unsigned int)
 
 template <int dim>
 std::vector<unsigned int>
-FE_RaviartThomas<dim>::get_dpo_vector (const unsigned int rt_order)
+FE_RaviartThomas<dim>::get_dpo_vector (const unsigned int deg)
 {
-                                   // the element is face-based (not
-                                   // to be confused with George
-                                   // W. Bush's Faith Based
-                                   // Initiative...), and we have
-                                   // (rt_order+1)^(dim-1) DoFs per face
+                                   // the element is face-based and we have
+                                   // (deg+1)^(dim-1) DoFs per face
   unsigned int dofs_per_face = 1;
-  for (unsigned int d=0; d<dim-1; ++d)
-    dofs_per_face *= rt_order+1;
+  for (unsigned int d=1; d<dim; ++d)
+    dofs_per_face *= deg+1;
 
                                    // and then there are interior dofs
   const unsigned int
-    interior_dofs = dim*rt_order*dofs_per_face;
+    interior_dofs = dim*deg*dofs_per_face;
   
   std::vector<unsigned int> dpo(dim+1);
   dpo[dim-1] = dofs_per_face;
@@ -498,7 +497,6 @@ FE_RaviartThomas<dim>::has_support_on_face (const unsigned int shape_index,
 }
 
 
-
 template <int dim>
 void
 FE_RaviartThomas<dim>::interpolate(
@@ -509,7 +507,6 @@ FE_RaviartThomas<dim>::interpolate(
 }
 
 
-
 template <int dim>
 void
 FE_RaviartThomas<dim>::interpolate(
index 8eb0fa196c73eb3504fc8144bf5133a69271d74d..738fad80ff1282ab86560195bcffe13d737a4352 100644 (file)
@@ -272,6 +272,35 @@ FE_RaviartThomasNodal<dim>::get_ria_vector (const unsigned int deg)
 }
 
 
+template <int dim>
+bool
+FE_RaviartThomasNodal<dim>::has_support_on_face (
+  const unsigned int shape_index,
+  const unsigned int face_index) const
+{
+  Assert (shape_index < this->dofs_per_cell,
+         ExcIndexRange (shape_index, 0, this->dofs_per_cell));
+  Assert (face_index < GeometryInfo<dim>::faces_per_cell,
+         ExcIndexRange (face_index, 0, GeometryInfo<dim>::faces_per_cell));
+
+                                  // The first degrees of freedom are
+                                  // on the faces and each face has
+                                  // degree degrees.
+  const unsigned int support_face = shape_index / this->degree;
+  
+                                  // The only thing we know for sure
+                                  // is that shape functions with
+                                  // support on one face are zero on
+                                  // the opposite face.
+  if (support_face < GeometryInfo<dim>::faces_per_cell)
+    return (face_index != GeometryInfo<dim>::opposite_face[support_face]);
+  
+                                  // In all other cases, return true,
+                                  // which is safe
+  return true;
+}
+
+
 template <int dim>
 void
 FE_RaviartThomasNodal<dim>::interpolate(
@@ -330,7 +359,6 @@ FE_RaviartThomasNodal<dim>::interpolate(
 }
 
 
-
 template <int dim>
 void
 FE_RaviartThomasNodal<dim>::interpolate(

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.