]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
move primitivity to FiniteElementData
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 17 Sep 2010 16:21:31 +0000 (16:21 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 17 Sep 2010 16:21:31 +0000 (16:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@22017 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_base.h
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_data.cc

index 0e448f8d11ba70890b90b35ae807420f0531b061..0c89f0ab561d253c6346730c253848081ddadf06 100644 (file)
@@ -1397,24 +1397,7 @@ class FiniteElement : public Subscriptor,
     bool
     is_primitive (const unsigned int i) const;
 
-                                    /**
-                                     * Return whether the entire
-                                     * finite element is primitive,
-                                     * in the sense that all its
-                                     * shape functions are
-                                     * primitive. If the finite
-                                     * element is scalar, then this
-                                     * is always the case.
-                                     *
-                                     * Since this is an extremely
-                                     * common operation, the result
-                                     * is cached in the
-                                     * #cached_primitivity
-                                     * variable which is computed in
-                                     * the constructor.
-                                     */
-    bool
-    is_primitive () const;
+    FiniteElementData<dim>::is_primitive;
     
                                     /**
                                      * Number of base elements in a
@@ -2000,18 +1983,6 @@ class FiniteElement : public Subscriptor,
     void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false,
                                                      const bool isotropic_prolongation_only=false);
     
-    
-                                    /**
-                                     * Store whether all shape
-                                     * functions are primitive. Since
-                                     * finding this out is a very
-                                     * common operation, we cache the
-                                     * result, i.e. compute the value
-                                     * in the constructor for simpler
-                                     * access.
-                                     */
-    const bool cached_primitivity;
-
                                     /**
                                      * Vector of projection
                                      * matrices. See
@@ -2862,22 +2833,12 @@ FiniteElement<dim,spacedim>::is_primitive (const unsigned int i) const
                                   //
                                   // for good measure, short circuit the test
                                   // if the entire FE is primitive
-  return ((cached_primitivity == true)
-         ||
+  return (is_primitive() ||
          (n_nonzero_components_table[i] == 1));
 }
 
 
 
-template <int dim, int spacedim>  
-inline
-bool
-FiniteElement<dim,spacedim>::is_primitive () const
-{
-  return cached_primitivity;
-}
-
-
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index dae7c287e52dbc761dd27a25967b498543f6e9d9..916467b88ca07f0586db686b68a78e9966b57dc1 100644 (file)
@@ -517,6 +517,24 @@ class FiniteElementData
                                      */
     unsigned int n_blocks () const;
 
+                                    /**
+                                     * Return whether the entire
+                                     * finite element is primitive,
+                                     * in the sense that all its
+                                     * shape functions are
+                                     * primitive. If the finite
+                                     * element is scalar, then this
+                                     * is always the case.
+                                     *
+                                     * Since this is an extremely
+                                     * common operation, the result
+                                     * is cached in the
+                                     * #cached_primitivity
+                                     * variable which is computed in
+                                     * the constructor.
+                                     */
+    bool is_primitive () const;
+    
                                     /**
                                      * Maximal polynomial degree of a
                                      * shape function in a single
@@ -598,6 +616,30 @@ class FiniteElementData
                                      * Comparison operator.
                                      */
     bool operator == (const FiniteElementData &) const;
+
+  protected:
+
+                                    /**
+                                     * Set the primitivity of the
+                                     * element. This is usually done
+                                     * by the constructor of a
+                                     * derived class.
+                                     * See @ref GlossPrimitive "primitive"
+                                     * for details.
+                                     */
+    void set_primitivity(const bool value);
+    
+  private:
+                                    /**
+                                     * Store whether all shape
+                                     * functions are primitive. Since
+                                     * finding this out is a very
+                                     * common operation, we cache the
+                                     * result, i.e. compute the value
+                                     * in the constructor for simpler
+                                     * access.
+                                     */
+    bool cached_primitivity;
 };
 
 
@@ -750,6 +792,24 @@ FiniteElementData<dim>::n_components () const
 
 
 
+template <int dim>  
+inline
+bool
+FiniteElementData<dim>::is_primitive () const
+{
+  return cached_primitivity;
+}
+
+
+template <int dim>  
+inline
+void
+FiniteElementData<dim>::set_primitivity (const bool value)
+{
+  cached_primitivity = value;
+}
+
+
 template <int dim>
 inline
 unsigned int 
index 3a44e0290aa5f6e1e61b24e406b5a39bef916d56..32784952bcc65c86676438f0919275b4baa4bda5 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -113,7 +113,6 @@ FiniteElement<dim,spacedim>::FiniteElement (
   const std::vector<std::vector<bool> > &nonzero_c)
                :
                FiniteElementData<dim> (fe_data),
-               cached_primitivity(false),
                adjust_quad_dof_index_for_face_orientation_table (dim == 3 ?
                                                                  this->dofs_per_quad : 0 ,
                                                                  dim==3 ? 8 : 0),
@@ -156,12 +155,11 @@ FiniteElement<dim,spacedim>::FiniteElement (
                                    // nonzero_components vector.
    const_cast<std::vector<unsigned int>&>
    (n_nonzero_components_table) = compute_n_nonzero_components(nonzero_components);
-   const_cast<bool&>
-   (cached_primitivity) = std::find_if (n_nonzero_components_table.begin(),
-                                     n_nonzero_components_table.end(),
-                                     std::bind2nd(std::not_equal_to<unsigned int>(),
-                                                  1U))
-                     == n_nonzero_components_table.end();
+   this->set_primitivity(std::find_if (n_nonzero_components_table.begin(),
+                                      n_nonzero_components_table.end(),
+                                      std::bind2nd(std::not_equal_to<unsigned int>(),
+                                                   1U))
+                        == n_nonzero_components_table.end());
    
    
   Assert (restriction_is_additive_flags.size() == this->dofs_per_cell,
@@ -189,7 +187,7 @@ FiniteElement<dim,spacedim>::FiniteElement (
                                   // only one (vector-)component; if
                                   // the element is not primitive,
                                   // leave these tables empty.
-  if (cached_primitivity)
+  if (this->is_primitive())
     {
       system_to_component_table.resize(this->dofs_per_cell);
       face_system_to_component_table.resize(this->dofs_per_face);
index 7b7588d9e978bafd32d0b82e1a1497d0a64dd281..fe483d21397bb6064f332d341ceb8e011aa36e53 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -34,7 +34,8 @@ FiniteElementData<dim>::FiniteElementData ()
                components(0),
                blocks(0),
                degree(0),
-               conforming_space(unknown)
+               conforming_space(unknown),
+               cached_primitivity(false)
 {}
 
 

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.