From: oliver Date: Mon, 26 Dec 2005 21:43:37 +0000 (+0000) Subject: Modified the hpFeValues classes to support quadrature rules which differ X-Git-Tag: v8.0.0~12735 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bbc7402889e35fe64b329a10ef3ec1210b34b815;p=dealii.git Modified the hpFeValues classes to support quadrature rules which differ 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 --- diff --git a/deal.II/deal.II/include/fe/hp_fe_values.h b/deal.II/deal.II/include/fe/hp_fe_values.h index 4d40559961..150df7c745 100644 --- a/deal.II/deal.II/include/fe/hp_fe_values.h +++ b/deal.II/deal.II/include/fe/hp_fe_values.h @@ -75,7 +75,8 @@ namespace internal * to derived classes. */ virtual - FEValues * create_fe_values (const FiniteElement &fe) const = 0; + FEValues * create_fe_values (const FiniteElement &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 &fe); - + select_fe_values (const FiniteElement &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 >, + std::map >, unsigned int>, boost::shared_ptr > fe_to_fe_values_map; /** @@ -284,7 +286,8 @@ class hpFEValues : public internal::FEValuesMap >, */ virtual FEValues * - create_fe_values (const FiniteElement &fe) const; + create_fe_values (const FiniteElement &fe, + const unsigned int active_fe_index) const; }; @@ -360,6 +363,30 @@ class hpFEFaceValues : public internal::FEValuesMap >, reinit (const typename hpDoFHandler::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::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 >, */ virtual FEFaceValues * - create_fe_values (const FiniteElement &fe) const; + create_fe_values (const FiniteElement &fe, + const unsigned int active_fe_index) const; }; @@ -445,6 +473,32 @@ class hpFESubfaceValues : public internal::FEValuesMap 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::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 */ virtual FESubfaceValues * - create_fe_values (const FiniteElement &fe) const; + create_fe_values (const FiniteElement &fe, + const unsigned int active_fe_index) const; }; diff --git a/deal.II/deal.II/source/fe/hp_fe_values.cc b/deal.II/deal.II/source/fe/hp_fe_values.cc index e5cae4e6a2..5dbf7d07de 100644 --- a/deal.II/deal.II/source/fe/hp_fe_values.cc +++ b/deal.II/deal.II/source/fe/hp_fe_values.cc @@ -27,21 +27,22 @@ namespace internal template FEValues & - FEValuesMap::select_fe_values (const FiniteElement &fe) + FEValuesMap::select_fe_values (const FiniteElement &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 (create_fe_values (fe)); + fe_to_fe_values_map[std::make_pair(&fe, active_fe_index)] + = boost::shared_ptr (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::reinit (const typename hpDoFHandler::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 FEValues * -hpFEValues::create_fe_values (const FiniteElement &fe) const +hpFEValues::create_fe_values (const FiniteElement &fe, + const unsigned int active_fe_index) const { return new FEValues ( - 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::reinit (const typename hpDoFHandler::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 +void +hpFEFaceValues::reinit (const typename hpDoFHandler::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 FEFaceValues * -hpFEFaceValues::create_fe_values (const FiniteElement &fe) const +hpFEFaceValues::create_fe_values (const FiniteElement &fe, + const unsigned int active_fe_index) const { return new FEFaceValues ( - 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::reinit (const typename hpDoFHandler::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 +void +hpFESubfaceValues::reinit (const typename hpDoFHandler::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 FESubfaceValues * -hpFESubfaceValues::create_fe_values (const FiniteElement &fe) const +hpFESubfaceValues::create_fe_values (const FiniteElement &fe, + const unsigned int active_fe_index) const { return new FESubfaceValues ( - 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); }