]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Removed polynomial degree.
authorLuca Heltai <luca.heltai@sissa.it>
Mon, 29 Feb 2016 16:46:00 +0000 (17:46 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 6 Apr 2016 11:00:14 +0000 (13:00 +0200)
include/deal.II/fe/mapping_manifold.h
source/fe/mapping_manifold.cc

index 2b072b9615740b20f9621225f61b60ead4426a8a..7819542107fb1d174e2a187a10f24f298d5084f1 100644 (file)
@@ -39,61 +39,28 @@ template <int,int> class MappingQ;
 
 
 /**
- * This class implements the functionality for polynomial mappings $Q_p$ of
- * polynomial degree $p$ that will be used on all cells of the mesh. The
- * MappingQ1 and MappingQ classes specialize this behavior slightly.
+ * This class implements the functionality for Manifold conforming
+ * mappings. This Mapping computes the transformation between the
+ * reference and real cell by exploiting the geometrical information
+ * coming from the underlying Manifold object.
  *
- * The class is poorly named. It should really have been called MappingQ
- * because it consistently uses $Q_p$ mappings on all cells of a
- * triangulation. However, the name MappingQ was already taken when we rewrote
- * the entire class hierarchy for mappings. One might argue that one should
- * always use MappingQGeneric over the existing class MappingQ (which, unless
- * explicitly specified during the construction of the object, only uses
- * mappings of degree $p$ <i>on cells at the boundary of the domain</i>). On
- * the other hand, there are good reasons to use MappingQ in many situations:
- * in many situations, curved domains are only provided with information about
- * how exactly edges at the boundary are shaped, but we do not know anything
- * about internal edges. Thus, in the absence of other information, we can
- * only assume that internal edges are straight lines, and in that case
- * internal cells may as well be treated is bilinear quadrilaterals or
- * trilinear hexahedra. (An example of how such meshes look is shown in step-1
- * already, but it is also discussed in the "Results" section of step-6.)
- * Because bi-/trilinear mappings are significantly cheaper to compute than
- * higher order mappings, it is advantageous in such situations to use the
- * higher order mapping only on cells at the boundary of the domain -- i.e.,
- * the behavior of MappingQ. Of course, MappingQGeneric also uses bilinear
- * mappings for interior cells as long as it has no knowledge about curvature
- * of interior edges, but it implements this the expensive way: as a general
- * $Q_p$ mapping where the mapping support points just <i>happen</i> to be
- * arranged along linear or bilinear edges or faces.
+ * Quadrature points computed using this mapping lye on the exact
+ * geometrical objects, and tangent and normal vectors computed using
+ * this class are normal and tangent to the underlying geometry. This
+ * is in constrast with the MappingQ class, which approximates the
+ * geometry using a polynomial of some order, and then computes the
+ * normals and tangents using the approximated surface.
  *
- * There are a number of special cases worth considering:
- * - If you really want to use a higher order mapping for all cells,
- * you can do this using the current class, but this only makes sense if you
- * can actually provide information about how interior edges and faces of the
- * mesh should be curved. This is typically done by associating a Manifold
- * with interior cells and edges. A simple example of this is discussed in the
- * "Results" section of step-6; a full discussion of manifolds is provided in
- * step-53.
- * - If you are working on meshes that describe a (curved) manifold
- * embedded in higher space dimensions, i.e., if dim!=spacedim, then every
- * cell is at the boundary of the domain you will likely already have attached
- * a manifold object to all cells that can then also be used by the mapping
- * classes for higher order mappings.
- *
- *
- * @author Wolfgang Bangerth, 2015
+ * @author Luca Heltai, Wolfgang Bangerth, Alberto Sartori 2016
  */
 template <int dim, int spacedim=dim>
 class MappingManifold : public Mapping<dim,spacedim>
 {
 public:
   /**
-   * Constructor.  @p polynomial_degree denotes the polynomial degree of the
-   * polynomials that are used to map cells from the reference to the real
-   * cell.
+   * Constructor.
    */
-  MappingManifold (const unsigned int polynomial_degree);
+  MappingManifold ();
 
   /**
    * Copy constructor.
@@ -104,12 +71,6 @@ public:
   virtual
   Mapping<dim,spacedim> *clone () const;
 
-  /**
-   * Return the degree of the mapping, i.e. the value which was passed to the
-   * constructor.
-   */
-  unsigned int get_degree () const;
-
   /**
    * Always returns @p true because the default implementation of functions in
    * this class preserves vertex locations.
@@ -208,10 +169,9 @@ public:
   {
   public:
     /**
-     * Constructor. The argument denotes the polynomial degree of the mapping
-     * to which this object will correspond.
+     * Constructor.
      */
-    InternalData(const unsigned int polynomial_degree);
+    InternalData();
 
     /**
      * Initialize the object's member variables related to cell data based on
index e777c0e301018a44612c38acac88c7ee61822598..5aa00c13a4de7b996f3e0ed0c8e3a74517c3bfec 100644 (file)
@@ -662,9 +662,9 @@ namespace internal
 
 
 template<int dim, int spacedim>
-MappingManifold<dim,spacedim>::InternalData::InternalData (const unsigned int polynomial_degree)
+MappingManifold<dim,spacedim>::InternalData::InternalData ()
   :
-  polynomial_degree (polynomial_degree),
+  polynomial_degree (1),
   n_shape_functions (Utilities::fixed_power<dim>(polynomial_degree+1))
 {}
 
@@ -821,8 +821,9 @@ namespace
 {
   template <int dim>
   std::vector<unsigned int>
-  get_dpo_vector (const unsigned int degree)
+  get_dpo_vector ()
   {
+    unsigned int degree = 1;
     std::vector<unsigned int> dpo(dim+1, 1U);
     for (unsigned int i=1; i<dpo.size(); ++i)
       dpo[i]=dpo[i-1]*(degree-1);
@@ -861,7 +862,7 @@ compute_shape_function_values (const std::vector<Point<dim> > &unit_points)
       const std::vector<unsigned int>
       renumber (FETools::
                 lexicographic_to_hierarchic_numbering (
-                  FiniteElementData<dim> (get_dpo_vector<dim>(polynomial_degree), 1,
+                  FiniteElementData<dim> (get_dpo_vector<dim>(), 1,
                                           polynomial_degree)));
 
       std::vector<double> values;
@@ -1155,16 +1156,14 @@ namespace
 
 
 template<int dim, int spacedim>
-MappingManifold<dim,spacedim>::MappingManifold (const unsigned int p)
+MappingManifold<dim,spacedim>::MappingManifold ()
   :
-  polynomial_degree(p),
+  polynomial_degree(1),
   line_support_points(this->polynomial_degree+1),
   fe_q(dim == 3 ? new FE_Q<dim>(this->polynomial_degree) : 0),
   support_point_weights_on_quad (compute_support_point_weights_on_quad<dim>(this->polynomial_degree)),
   support_point_weights_on_hex (compute_support_point_weights_on_hex<dim>(this->polynomial_degree))
 {
-  Assert (p >= 1, ExcMessage ("It only makes sense to create polynomial mappings "
-                              "with a polynomial degree greater or equal to one."));
 }
 
 
@@ -1191,16 +1190,6 @@ MappingManifold<dim,spacedim>::clone () const
 
 
 
-
-template<int dim, int spacedim>
-unsigned int
-MappingManifold<dim,spacedim>::get_degree() const
-{
-  return polynomial_degree;
-}
-
-
-
 template<int dim, int spacedim>
 Point<spacedim>
 MappingManifold<dim,spacedim>::
@@ -1932,7 +1921,7 @@ typename MappingManifold<dim,spacedim>::InternalData *
 MappingManifold<dim,spacedim>::get_data (const UpdateFlags update_flags,
                                          const Quadrature<dim> &q) const
 {
-  InternalData *data = new InternalData(polynomial_degree);
+  InternalData *data = new InternalData();
   data->initialize (this->requires_update_flags(update_flags), q, q.size());
 
   return data;

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.