]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Patch by Joshua White: add documentation, implement face domination.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 26 Sep 2009 00:47:21 +0000 (00:47 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 26 Sep 2009 00:47:21 +0000 (00:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@19569 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 512e7222355d746fce6c70d20271bc9c2c7e9d50..02a17f4f7381cdb96e6c369b0a1b182934294456 100644 (file)
@@ -39,40 +39,145 @@ class FE_Nothing : public FiniteElement<dim>
 {
   public:
 
+                                    /**
+                                      * Constructor.  
+                                      */
     FE_Nothing ();
-
+    
+                                     /**
+                                      * A sort of virtual copy
+                                      * constructor. Some places in
+                                      * the library, for example the
+                                      * constructors of FESystem as
+                                      * well as the hp::FECollection
+                                      * class, need to make copied of
+                                      * finite elements without
+                                      * knowing their exact type. They
+                                      * do so through this function.
+                                      */
     virtual
     FiniteElement<dim> *
     clone() const;
 
+                                     /**
+                                      * Return a string that uniquely
+                                      * identifies a finite
+                                      * element. In this case it is
+                                      * "FE_Nothing<dim>()
+                                      */
     virtual
     std::string
     get_name() const;
 
+                                     /**
+                                      * Number of base elements in a
+                                      * mixed discretization. In this case
+                                      * we only have one.
+                                      */
     virtual
     unsigned int
     n_base_elements () const;
-
+    
+                                     /**
+                                      * Access to base element
+                                      * objects. Since the element is
+                                      * scalar, then
+                                      * <code>base_element(0)</code> is
+                                      * @p this.
+                                      */
     virtual
     const FiniteElement<dim> &
     base_element (const unsigned int index) const;
 
+                                     /**
+                                      * This index denotes how often
+                                      * the base element @p index is
+                                      * used in a composed element. If
+                                      * the element is scalar, then
+                                      * the result is always equal to
+                                      * one. See the documentation for
+                                      * the n_base_elements()
+                                      * function for more details.
+                                      */
     virtual
     unsigned int
     element_multiplicity (const unsigned int index) const;
 
+                                     /**
+                                      * Determine the values a finite
+                                      * element should compute on
+                                      * initialization of data for
+                                      * FEValues.
+                                      *
+                                      * Given a set of flags
+                                      * indicating what quantities are
+                                      * requested from a FEValues
+                                      * object, update_once() and
+                                      * update_each() compute which
+                                      * values must really be
+                                      * computed. Then, the
+                                      * <tt>fill_*_values</tt> functions
+                                      * are called with the result of
+                                      * these.
+                                      *
+                                      * In this case, since the element 
+                                      * has zero degrees of freedom and
+                                      * no information can be computed on
+                                      * it, this function simply returns
+                                      * the default (empty) set of update
+                                      * flags.
+                                      */
+
     virtual
     UpdateFlags
     update_once (const UpdateFlags flags) const;
 
+                                     /**
+                                      * Complementary function for
+                                      * update_once().
+                                      *
+                                      * While update_once() returns
+                                      * the values to be computed on
+                                      * the unit cell for yielding the
+                                      * required data, this function
+                                      * determines the values that
+                                      * must be recomputed on each
+                                      * cell.
+                                      *
+                                      * Refer to update_once() for
+                                      * more details.
+                                      */
     virtual
     UpdateFlags
     update_each (const UpdateFlags flags) const;
 
+                                     /**
+                                      * Return the value of the
+                                      * @p ith shape function at the
+                                      * point @p p. @p p is a point
+                                      * on the reference element. Because the
+                                      * current element has no degrees of freedom,
+                                      * this function should obviously not be
+                                      * called in practice.  All this function
+                                      * really does, therefore, is trigger an
+                                      * exception.
+                                      */
     virtual
     double
     shape_value (const unsigned int i, const Point<dim> &p) const;
 
+                                     /**
+                                      * Fill the fields of
+                                      * FEValues. This function
+                                      * performs all the operations
+                                      * needed to compute the data of an
+                                      * FEValues object.
+                                      *
+                                      * In the current case, this function
+                                      * returns no meaningful information, 
+                                      * since the element has no degrees of
+                                      * freedom.
+                                      */
     virtual
     void
     fill_fe_values (const Mapping<dim> & mapping,
@@ -83,6 +188,18 @@ class FE_Nothing : public FiniteElement<dim>
                    FEValuesData<dim,dim> & data,
                    CellSimilarity::Similarity & cell_similarity) const;
 
+                                     /**
+                                      * Fill the fields of
+                                      * FEFaceValues. This function
+                                      * performs all the operations
+                                      * needed to compute the data of an
+                                      * FEFaceValues object.
+                                      *
+                                      * In the current case, this function
+                                      * returns no meaningful information, 
+                                      * since the element has no degrees of
+                                      * freedom.
+                                      */
     virtual
     void
     fill_fe_face_values (const Mapping<dim> & mapping,
@@ -93,6 +210,18 @@ class FE_Nothing : public FiniteElement<dim>
                         typename Mapping<dim> :: InternalDataBase & fedata,
                         FEValuesData<dim,dim> & data) const;
 
+                                     /**
+                                      * Fill the fields of
+                                      * FESubFaceValues. This function
+                                      * performs all the operations
+                                      * needed to compute the data of an
+                                      * FESubFaceValues object.
+                                      *
+                                      * In the current case, this function
+                                      * returns no meaningful information, 
+                                      * since the element has no degrees of
+                                      * freedom.
+                                      */
     virtual
     void
     fill_fe_subface_values (const Mapping<dim> & mapping,
@@ -104,11 +233,50 @@ class FE_Nothing : public FiniteElement<dim>
                            typename Mapping<dim>::InternalDataBase & fedata,
                            FEValuesData<dim,dim> & data) const;
 
+                                     /**
+                                      * Prepare internal data
+                                      * structures and fill in values
+                                      * independent of the
+                                      * cell. Returns a pointer to an
+                                      * object of which the caller of
+                                      * this function then has to
+                                      * assume ownership (which
+                                      * includes destruction when it
+                                      * is no more needed).
+                                      *
+                                      * In the current case, this function 
+                                      * just returns a default pointer, since
+                                      * no meaningful data exists for this 
+                                      * element.
+                                      */
     virtual
     typename Mapping<dim>::InternalDataBase *
     get_data (const UpdateFlags     update_flags,
              const Mapping<dim>    & mapping,
              const Quadrature<dim> & quadrature) const;
+        
+                                     /**
+                                      * Return whether this element dominates
+                                      * the one given as argument when they
+                                      * meet at a common face,
+                                      * whether it is the other way around,
+                                      * whether neither dominates, or if
+                                      * either could dominate.
+                                      *
+                                      * For a definition of domination, see
+                                      * FiniteElementBase::Domination and in
+                                      * particular the @ref hp_paper "hp paper".
+                                      *
+                                      * In the current case, the other element
+                                      * is always assumed to dominate, unless 
+                                      * it is also of type FE_Nothing().  In
+                                      * that situation, either element can
+                                      * dominate.
+                                      */
+    virtual
+    FiniteElementDomination::Domination
+    compare_for_face_domination (const FiniteElement<dim> & fe_other) const;
+    
 };
 
 
index d4d7555ae963ad7365070e0b5e27a612562a49ff..eb86ec5cf148b47a9950fc0a60421269c436c0de 100644 (file)
@@ -120,8 +120,8 @@ FE_Nothing<dim>::get_data (const UpdateFlags  /*flags*/,
                           const Mapping<dim> & /*mapping*/,
                           const Quadrature<dim> & /*quadrature*/) const
 {
-  Assert(false,ExcMessage(zero_dof_message));
-  return NULL;
+  typename Mapping<dim>::InternalDataBase* data = new typename Mapping<dim>::InternalDataBase();
+  return data;
 }
 
 
@@ -129,7 +129,7 @@ FE_Nothing<dim>::get_data (const UpdateFlags  /*flags*/,
 template <int dim>
 void
 FE_Nothing<dim>::
-fill_fe_values (const Mapping<dim>                      & /*mapping*/,
+fill_fe_values (const Mapping<dim> & /*mapping*/,
                const typename Triangulation<dim>::cell_iterator & /*cell*/,
                const Quadrature<dim> & /*quadrature*/,
                typename Mapping<dim>::InternalDataBase & /*mapping_data*/,
@@ -174,6 +174,17 @@ fill_fe_subface_values (const Mapping<dim> & /*mapping*/,
 }
 
 
+template <int dim>
+FiniteElementDomination::Domination
+FE_Nothing<dim> ::
+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;
+  else
+    return FiniteElementDomination::other_element_dominates;
+}
+    
 
 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.