]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Patch by Joshua White: Implementation of more functions of FE_Nothing.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 1 Oct 2009 20:08:07 +0000 (20:08 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 1 Oct 2009 20:08:07 +0000 (20:08 +0000)
git-svn-id: https://svn.dealii.org/trunk@19649 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_nothing.h
deal.II/deal.II/source/fe/fe_nothing.cc

index 02a17f4f7381cdb96e6c369b0a1b182934294456..a980d7cb400ea61e2b7aaa23e7c9d94308ad7953 100644 (file)
@@ -40,9 +40,11 @@ class FE_Nothing : public FiniteElement<dim>
   public:
 
                                     /**
-                                      * Constructor.  
+                                      * Constructor. Argument denotes the
+                                      * number of components to give this
+                                      * finite element (default = 1).  
                                       */
-    FE_Nothing ();
+    FE_Nothing (unsigned int n_components = 1);
     
                                      /**
                                       * A sort of virtual copy
@@ -277,6 +279,65 @@ class FE_Nothing : public FiniteElement<dim>
     FiniteElementDomination::Domination
     compare_for_face_domination (const FiniteElement<dim> & fe_other) const;
     
+    
+    
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
+    
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
+    
+    virtual
+    std::vector<std::pair<unsigned int, unsigned int> >
+    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+    
+    virtual
+    bool
+    hp_constraints_are_implemented () const;
+
+                                      /**
+                                      * Return the matrix
+                                      * interpolating from a face of
+                                      * of one element to the face of
+                                      * the neighboring element. 
+                                      * The size of the matrix is
+                                      * then <tt>source.#dofs_per_face</tt> times
+                                      * <tt>this->#dofs_per_face</tt>.
+                                      *
+                                      * Since the current finite element has no
+                                      * degrees of freedom, the interpolation
+                                      * matrix is necessarily empty.
+                                      */
+
+    virtual
+    void
+    get_face_interpolation_matrix (const FiniteElement<dim> &source_fe,
+                                   FullMatrix<double>       &interpolation_matrix) const;
+                                   
+                                   
+                                     /**
+                                      * Return the matrix
+                                      * interpolating from a face of
+                                      * of one element to the subface of
+                                      * the neighboring element. 
+                                      * The size of the matrix is
+                                      * then <tt>source.#dofs_per_face</tt> times
+                                      * <tt>this->#dofs_per_face</tt>.
+                                      *
+                                      * Since the current finite element has no
+                                      * degrees of freedom, the interpolation
+                                      * matrix is necessarily empty.
+                                      */
+                                      
+    virtual
+    void
+    get_subface_interpolation_matrix (const FiniteElement<dim> & source_fe,
+                                      const unsigned int index,
+                                      FullMatrix<double>  &interpolation_matrix) const;
+
+    
 };
 
 
index eb86ec5cf148b47a9950fc0a60421269c436c0de..78ea19564c0988218e69e4fa9dc34af0810de11e 100644 (file)
@@ -25,11 +25,11 @@ namespace
 
 
 template <int dim>
-FE_Nothing<dim>::FE_Nothing ()
+FE_Nothing<dim>::FE_Nothing (const unsigned n_components)
                 :
                 FiniteElement<dim>
                (FiniteElementData<dim>(std::vector<unsigned>(dim+1,0),
-                                       1, 0,
+                                       n_components, 0,
                                        FiniteElementData<dim>::unknown),
                 std::vector<bool>(),
                 std::vector<std::vector<bool> >() )
@@ -120,6 +120,9 @@ FE_Nothing<dim>::get_data (const UpdateFlags  /*flags*/,
                           const Mapping<dim> & /*mapping*/,
                           const Quadrature<dim> & /*quadrature*/) const
 {
+                // Create a default data object.  Normally we would then
+                // need to resize things to hold the appropriate numbers
+                // of dofs, but in this case all data fields are empty.
   typename Mapping<dim>::InternalDataBase* data = new typename Mapping<dim>::InternalDataBase();
   return data;
 }
@@ -137,7 +140,7 @@ fill_fe_values (const Mapping<dim> & /*mapping*/,
                FEValuesData<dim,dim> & /*data*/,
                CellSimilarity::Similarity & /*cell_similarity*/) const
 {
-  Assert(false,ExcMessage(zero_dof_message));
+                // leave data fields empty
 }
 
 
@@ -153,11 +156,9 @@ fill_fe_face_values (const Mapping<dim> & /*mapping*/,
                     typename Mapping<dim>::InternalDataBase & /*fedata*/,
                     FEValuesData<dim,dim> & /*data*/) const
 {
-  Assert(false,ExcMessage(zero_dof_message));
+                // leave data fields empty
 }
 
-
-
 template <int dim>
 void
 FE_Nothing<dim>::
@@ -170,7 +171,7 @@ fill_fe_subface_values (const Mapping<dim> & /*mapping*/,
                        typename Mapping<dim>::InternalDataBase & /*fedata*/,
                        FEValuesData<dim,dim> & /*data*/) const
 {
-  Assert(false,ExcMessage(zero_dof_message));
+                // leave data fields empty
 }
 
 
@@ -181,10 +182,103 @@ compare_for_face_domination (const FiniteElement<dim> & fe_other) const
 {
   if(dynamic_cast<const FE_Nothing<dim>*>(&fe_other) != 0)
     return FiniteElementDomination::either_element_can_dominate;
+    
+              // Question: Best to assume this element always dominates,
+              // since we will always no how to do deal with other elements
+              // regardless of type, or do we leave it up to the other 
+              // element to handle the FENothing appropriately?  Set
+              // to the second choice for now.
   else
     return FiniteElementDomination::other_element_dominates;
+    //return FiniteElementDomination::this_element_dominates;
+}
+   
+   
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Nothing<dim> ::
+hp_vertex_dof_identities (const FiniteElement<dim> &/*fe_other*/) const
+{
+               // the FE_Nothing has no
+                                      // degrees of freedom, so there
+                                      // are no equivalencies to be
+                                      // recorded
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
 }
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Nothing<dim> ::
+hp_line_dof_identities (const FiniteElement<dim> &/*fe_other*/) const
+{
+               // the FE_Nothing has no
+                                      // degrees of freedom, so there
+                                      // are no equivalencies to be
+                                      // recorded
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+   
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Nothing<dim> ::
+hp_quad_dof_identities (const FiniteElement<dim> &/*fe_other*/) const
+{
+               // the FE_Nothing has no
+                                      // degrees of freedom, so there
+                                      // are no equivalencies to be
+                                      // recorded
+      return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
     
+template <int dim>
+bool
+FE_Nothing<dim> ::
+hp_constraints_are_implemented () const
+{
+  return true;
+}
+  
+        
+template <int dim>
+void
+FE_Nothing<dim>::
+get_face_interpolation_matrix (const FiniteElement<dim> &/*source_fe*/,
+                               FullMatrix<double>       &interpolation_matrix) const
+{
+                                   // since this element has no face dofs, the
+                                   // interpolation matrix is necessarily empty
+                                   
+  Assert (interpolation_matrix.m() == 0,
+          ExcDimensionMismatch (interpolation_matrix.m(),
+                                0));
+  Assert (interpolation_matrix.n() == 0,
+          ExcDimensionMismatch (interpolation_matrix.m(),
+                                0));
+}
+
+
+template <int dim>
+void
+FE_Nothing<dim>::
+get_subface_interpolation_matrix (const FiniteElement<dim> & /*source_fe*/,
+                                  const unsigned int /*index*/,
+                                  FullMatrix<double>  &interpolation_matrix) const
+{
+                                   // since this element has no face dofs, the
+                                   // interpolation matrix is necessarily empty
+
+  Assert (interpolation_matrix.m() == 0,
+          ExcDimensionMismatch (interpolation_matrix.m(),
+                                0));
+  Assert (interpolation_matrix.n() == 0,
+          ExcDimensionMismatch (interpolation_matrix.m(),
+                                0));
+}
+
+       
 
 template class FE_Nothing<deal_II_dimension>;
 

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.