]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Added a method "hp_constraints_are_implemented" which returns a flag
authorkayser-herold <kayser-herold@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 28 Jul 2006 16:46:00 +0000 (16:46 +0000)
committerkayser-herold <kayser-herold@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 28 Jul 2006 16:46:00 +0000 (16:46 +0000)
indicating wheter the FiniteElement already supports the new style
of hanging node matrices.

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

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

index 37aaef223e4dffc0089aa50a876222f94e2aefd1..86bd91386a03637d231b7896a28affc2f2efccd2 100644 (file)
@@ -874,6 +874,47 @@ class FiniteElement : public Subscriptor,
                                       */
     bool constraints_are_implemented () const;
 
+
+                                     /**
+                                      * Return whether this element
+                                      * implements its hanging node
+                                      * constraints in the new way,
+                                     * which has to be used to make
+                                     * elements "hp compatible".
+                                     * That means, the element properly
+                                     * implements the
+                                     * get_face_interpolation_matrix
+                                     * and get_subface_interpolation_matrix
+                                     * methods. Therefore the return
+                                     * value also indicates whether a call
+                                      * to the @p get_face_interpolation_matrix
+                                     * method and the @p get_subface_interpolation_matrix
+                                     * method will generate an error or not.
+                                      *
+                                     * Currently the main purpose
+                                     * of this function is to
+                                     * allow the make_hanging_node_constraints
+                                     * method to decide wheter the
+                                     * new procedures, which are supposed
+                                     * to work in the hp framework
+                                     * can be used, or if the old
+                                     * well verified but not hp capable
+                                     * functions should be used.
+                                     * Once, the transition to the new
+                                     * scheme for computing the
+                                     * interface constraints is completed,
+                                     * this function will probably mostly
+                                     * superfluous.
+                                     *
+                                     * Derived classes should implement
+                                     * this function accordingly. The
+                                     * default assumption is that a
+                                     * finite element does not provide
+                                     * hp capable face interpolation.
+                                      */
+    virtual bool hp_constraints_are_implemented () const;
+
+
                                     /**
                                      * Return the matrix
                                      * interpolating from the given
@@ -929,7 +970,6 @@ class FiniteElement : public Subscriptor,
     virtual void
     get_face_interpolation_matrix (const FiniteElement<dim> &source,
                                   FullMatrix<double>       &matrix) const;
-                                    //@}
     
 
                                     /**
index a3e6915032cfd49ee43558e97fa6c3b2d93c07f0..92b0873fad38c055069f54b91c0e7d9180f2f991 100644 (file)
@@ -327,8 +327,24 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      const unsigned int        subface,
                                      FullMatrix<double>       &matrix) const;
                                     //@}
+
     
+                                     /**
+                                      * Return whether this element
+                                      * implements its hanging node
+                                      * constraints in the new way,
+                                     * which has to be used to make
+                                     * elements "hp compatible".
+                                      *
+                                     * For the FE_Q class the
+                                     * result is always true, as
+                                     * it implements the complete
+                                     * set of functions necessary
+                                     * for hp capability.
+                                      */
+    virtual bool hp_constraints_are_implemented () const;
 
+    
                                     /**
                                      * Check for non-zero values on a face.
                                      *
index dcc9cc824b4b316d7eac60d3cc85517abc0d0f82..fc87a263c12589cf760b00abda0719930b81d612 100644 (file)
@@ -383,6 +383,15 @@ FiniteElement<dim>::constraints_are_implemented () const
 
 
 
+template <int dim>
+bool
+FiniteElement<dim>::hp_constraints_are_implemented () const
+{
+  return (this->dofs_per_face  == 0);
+}
+
+
+
 template <int dim>
 const FullMatrix<double> &
 FiniteElement<dim>::constraints () const
index fb7606b83b886d4f3e0198ad6c735c5598bf4c6b..34bfa21f0dc1fc15c1e7f676172d2e9b6038d39d 100644 (file)
@@ -528,6 +528,15 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
 
 #endif
 
+
+template <int dim>
+bool
+FE_Q<dim>::hp_constraints_are_implemented () const
+{
+  return true;
+}
+
+
 template <int dim>
 std::vector<std::pair<unsigned int, unsigned int> >
 FE_Q<dim>::

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.