]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Modified the hpFeValues classes to support quadrature rules which differ
authoroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 26 Dec 2005 21:43:37 +0000 (21:43 +0000)
committeroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 26 Dec 2005 21:43:37 +0000 (21:43 +0000)
from the default quadrature rule on a face. This leads to about 10-15%
performance gain during matrix assembly.

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

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

index 4d40559961c91a5d44ca1172342628b7bb5b44cd..150df7c745be73d8e2b4667520c7dbcc76b218ed 100644 (file)
@@ -75,7 +75,8 @@ namespace internal
                                         * to derived classes.
                                         */
       virtual
-      FEValues * create_fe_values (const FiniteElement<dim> &fe) const = 0;
+      FEValues * create_fe_values (const FiniteElement<dim> &fe,
+                                  const unsigned int active_fe_index) const = 0;
       
                                        /**
                                         * Select the @p FEValues
@@ -92,8 +93,9 @@ namespace internal
                                         * this object is returned.
                                         */
       FEValues &
-      select_fe_values (const FiniteElement<dim> &fe);
-      
+      select_fe_values (const FiniteElement<dim> &fe,
+                       const unsigned int active_fe_index);
+
 
                                        /**
                                         * This field remembers the
@@ -116,7 +118,7 @@ namespace internal
                                         * destructor of this class is
                                         * run.
                                         */
-      std::map<SmartPointer<const FiniteElement<dim> >,
+      std::map<std::pair<SmartPointer<const FiniteElement<dim> >, unsigned int>,
                boost::shared_ptr<FEValues> > fe_to_fe_values_map;
 
                                        /**
@@ -284,7 +286,8 @@ class hpFEValues : public internal::FEValuesMap<dim,FEValues<dim> >,
                                       */
     virtual
     FEValues<dim> *
-    create_fe_values (const FiniteElement<dim> &fe) const;
+    create_fe_values (const FiniteElement<dim> &fe,
+                     const unsigned int active_fe_index) const;
 };
 
 
@@ -360,6 +363,30 @@ class hpFEFaceValues : public internal::FEValuesMap<dim,FEFaceValues<dim> >,
     reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
             const unsigned int face_no);
 
+
+                                     /**
+                                      * Reinitialize the object for
+                                      * the given cell. In this
+                                     * method, the user can specify
+                                     * which active_fe_index
+                                     * should be used to initialise
+                                     * the underlying FEValues
+                                     * object. This functionality
+                                     * is required, if the face terms
+                                     * between two cells with different
+                                     * polynomial degree should be
+                                     * assembled. In this case the
+                                     * values on the face of the
+                                     * lower order element have to
+                                     * be evaluated at the quadratrure
+                                     * points of the higher order
+                                     * element.
+                                      */
+    void
+    reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
+            const unsigned int face_no,
+           const unsigned int active_fe_index);
+
   protected:
                                      /**
                                       * Create an object of type
@@ -368,7 +395,8 @@ class hpFEFaceValues : public internal::FEValuesMap<dim,FEFaceValues<dim> >,
                                       */
     virtual
     FEFaceValues<dim> *
-    create_fe_values (const FiniteElement<dim> &fe) const;
+    create_fe_values (const FiniteElement<dim> &fe,
+                     const unsigned int active_fe_index) const;
 };
 
 
@@ -445,6 +473,32 @@ class hpFESubfaceValues : public internal::FEValuesMap<dim,FESubfaceValues<dim>
             const unsigned int face_no,
             const unsigned int subface_no);
 
+
+                                     /**
+                                      * Reinitialize the object for
+                                      * the given cell. In this
+                                     * method, the user can specify
+                                     * which active_fe_index
+                                     * should be used to initialise
+                                     * the underlying FEValues
+                                     * object. This functionality
+                                     * is required, if the face terms
+                                     * between two cells with different
+                                     * polynomial degree should be
+                                     * assembled. In this case the
+                                     * values on the face of the
+                                     * lower order element have to
+                                     * be evaluated at the quadratrure
+                                     * points of the higher order
+                                     * element.
+                                      */
+    void
+    reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
+            const unsigned int face_no,
+            const unsigned int subface_no,
+           const unsigned int active_fe_index);
+
+
   protected:
                                      /**
                                       * Create an object of type
@@ -453,7 +507,8 @@ class hpFESubfaceValues : public internal::FEValuesMap<dim,FESubfaceValues<dim>
                                       */
     virtual
     FESubfaceValues<dim> *
-    create_fe_values (const FiniteElement<dim> &fe) const;
+    create_fe_values (const FiniteElement<dim> &fe,
+                     const unsigned int active_fe_index) const;
 };
 
 
index e5cae4e6a299224ab302e8c7ff0ca7fcca22f815..5dbf7d07de0c36a5c243810620895a8814ac6e73 100644 (file)
@@ -27,21 +27,22 @@ namespace internal
   
   template <int dim, class FEValues>
   FEValues &
-  FEValuesMap<dim,FEValues>::select_fe_values (const FiniteElement<dim> &fe)
+  FEValuesMap<dim,FEValues>::select_fe_values (const FiniteElement<dim> &fe,
+                                              const unsigned int active_fe_index)
   {
                                      // check if the finite element
                                      // does not exist as a key in the
                                      // map
-    if (fe_to_fe_values_map.find (&fe) ==
+    if (fe_to_fe_values_map.find (std::make_pair(&fe, active_fe_index)) ==
         fe_to_fe_values_map.end())
                                        // a-ha! doesn't yet, so let's
                                        // make it up
-      fe_to_fe_values_map[&fe]
-        = boost::shared_ptr<FEValues> (create_fe_values (fe));
+      fe_to_fe_values_map[std::make_pair(&fe, active_fe_index)]
+        = boost::shared_ptr<FEValues> (create_fe_values (fe, active_fe_index));
 
 
                                      // now there definitely is one!
-    present_fe_values = fe_to_fe_values_map[&fe];
+    present_fe_values = fe_to_fe_values_map[std::make_pair (&fe, active_fe_index)];
 
     return *present_fe_values;
   }
@@ -109,18 +110,19 @@ void
 hpFEValues<dim>::reinit (const typename hpDoFHandler<dim>::cell_iterator &cell)
 {
   this->present_fe_index = cell->active_fe_index ();
-  this->select_fe_values (cell->get_fe()).reinit (cell);
+  this->select_fe_values (cell->get_fe(), this->present_fe_index).reinit (cell);
 }
 
 
 
 template <int dim>
 FEValues<dim> *
-hpFEValues<dim>::create_fe_values (const FiniteElement<dim> &fe) const
+hpFEValues<dim>::create_fe_values (const FiniteElement<dim> &fe,
+                                  const unsigned int active_fe_index) const
 {
   return new FEValues<dim> (
-      this->mapping_collection.get_mapping (this->present_fe_index), fe,
-      this->qcollection.get_quadrature (this->present_fe_index),
+      this->mapping_collection.get_mapping (active_fe_index), fe,
+      this->qcollection.get_quadrature (active_fe_index),
       this->update_flags);
 }
 
@@ -156,18 +158,29 @@ hpFEFaceValues<dim>::reinit (const typename hpDoFHandler<dim>::cell_iterator &ce
                              const unsigned int face_no)
 {
   this->present_fe_index = cell->active_fe_index ();
-  this->select_fe_values (cell->get_fe()).reinit (cell, face_no);
+  this->select_fe_values (cell->get_fe(), this->present_fe_index).reinit (cell, face_no);
 }
 
 
+template <int dim>
+void
+hpFEFaceValues<dim>::reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
+                            const unsigned int face_no,
+                            const unsigned int active_fe_index)
+{
+  this->present_fe_index = active_fe_index;
+  this->select_fe_values (cell->get_fe(), active_fe_index).reinit (cell, face_no);
+}
+
 
 template <int dim>
 FEFaceValues<dim> *
-hpFEFaceValues<dim>::create_fe_values (const FiniteElement<dim> &fe) const
+hpFEFaceValues<dim>::create_fe_values (const FiniteElement<dim> &fe,
+                                      const unsigned int active_fe_index) const
 {
   return new FEFaceValues<dim> (
-      this->mapping_collection.get_mapping (this->present_fe_index), fe,
-      this->qcollection.get_quadrature (this->present_fe_index),
+      this->mapping_collection.get_mapping (active_fe_index), fe,
+      this->qcollection.get_quadrature (active_fe_index),
       this->update_flags);
 }
 
@@ -204,18 +217,30 @@ hpFESubfaceValues<dim>::reinit (const typename hpDoFHandler<dim>::cell_iterator
                                 const unsigned int subface_no)
 {
   this->present_fe_index = cell->active_fe_index ();
-  this->select_fe_values (cell->get_fe()).reinit (cell, face_no, subface_no);
+  this->select_fe_values (cell->get_fe(), this->present_fe_index).reinit (cell, face_no, subface_no);
 }
 
 
+template <int dim>
+void
+hpFESubfaceValues<dim>::reinit (const typename hpDoFHandler<dim>::cell_iterator &cell,
+                                const unsigned int face_no,
+                                const unsigned int subface_no,
+                               const unsigned int active_fe_index)
+{
+  this->present_fe_index = active_fe_index;
+  this->select_fe_values (cell->get_fe(), active_fe_index).reinit (cell, face_no, subface_no);
+}
+
 
 template <int dim>
 FESubfaceValues<dim> *
-hpFESubfaceValues<dim>::create_fe_values (const FiniteElement<dim> &fe) const
+hpFESubfaceValues<dim>::create_fe_values (const FiniteElement<dim> &fe,
+                                         const unsigned int active_fe_index) const
 {
   return new FESubfaceValues<dim> (
-      this->mapping_collection.get_mapping (this->present_fe_index), fe,
-      this->qcollection.get_quadrature (this->present_fe_index),
+      this->mapping_collection.get_mapping (active_fe_index), fe,
+      this->qcollection.get_quadrature (active_fe_index),
       this->update_flags);
 }
 

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.