]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
simplified constructor for FiniteElementBase
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 29 Jun 2005 15:41:20 +0000 (15:41 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 29 Jun 2005 15:41:20 +0000 (15:41 +0000)
git-svn-id: https://svn.dealii.org/trunk@10970 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 9f39c089233ab0d2f13609c4470fe5f1485f4dd1..82dae1329bda5c3c889945cd9be1433401e59734 100644 (file)
@@ -362,7 +362,7 @@ class FiniteElementData
  * transfer, respectively.
  *
  *
- * @sect3{Support points}
+ * <h3>Support points</h3>
  *
  * Since a FiniteElement does not have information on the actual
  * grid cell, it can only provide support points on the unit
@@ -391,7 +391,7 @@ class FiniteElementData
  * If the mapping of all support points is needed, the first variant should
  * be preferred for efficiency.
  *
- * @sect3{Finite elements in one dimension}
+ * <h3>Finite elements in one dimension</h3>
  *
  * Finite elements in one dimension need only set the #restriction
  * and #prolongation matrices. The constructor of this class in one
@@ -399,7 +399,7 @@ class FiniteElementData
  * dimension zero. Changing this behaviour in derived classes is
  * generally not a reasonable idea and you risk getting into trouble.
  * 
- * @sect3{Finite elements in two dimensions}
+ * <h3>Finite elements in two dimensions</h3>
  * 
  * In addition to the fields already present in 1D, a constraint
  * matrix is needed, if the finite element has node values located on
@@ -433,7 +433,7 @@ class FiniteElementData
  * at the time of this writing whether this is a constraint itself.
  *
  * 
- * @sect3{Finite elements in three dimensions}
+ * <h3>Finite elements in three dimensions</h3>
  *
  * For the interface constraints, almost the same holds as for the 2D case.
  * The numbering for the indices $n$ on the mother face is obvious and keeps
@@ -580,9 +580,18 @@ class FiniteElementBase : public Subscriptor,
                                      * existing finite element
                                      * classes. For the second and
                                      * third parameter of this
-                                     * constructor, see the document
+                                     * constructor, see the documentation
                                      * of the respective member
                                      * variables.
+                                     *
+                                     * @note Both vector parameters
+                                     * should have length
+                                     * <tt>dofs_per_cell</tt>. Nevertheless,
+                                     * it is allowed to use vectors
+                                     * of length one. In this case,
+                                     * the vector is resized to the
+                                     * correct length and filled with
+                                     * the entry value.
                                      */
     FiniteElementBase (const FiniteElementData<dim> &fe_data,
                       const std::vector<bool> &restriction_is_additive_flags,
index 7b973f5dfe3ea884ffcbce315060f5196625d74d..994005544a2cb7e18483ccd27f3e71cc6ce1ffbd 100644 (file)
@@ -106,10 +106,10 @@ FiniteElementBase<dim>::InternalDataBase::~InternalDataBase ()
 
 
 template <int dim>
-FiniteElementBase<dim>::
-FiniteElementBase (const FiniteElementData<dim> &fe_data,
-                   const std::vector<bool> &restriction_is_additive_flags,
-                   const std::vector<std::vector<bool> > &nonzero_components)
+FiniteElementBase<dim>::FiniteElementBase (
+  const FiniteElementData<dim> &fe_data,
+  const std::vector<bool> &r_i_a_f,
+  const std::vector<std::vector<bool> > &nonzero_c)
                :
                FiniteElementData<dim> (fe_data),
                 system_to_component_table (this->dofs_per_cell),
@@ -118,30 +118,47 @@ FiniteElementBase (const FiniteElementData<dim> &fe_data,
                 face_system_to_base_table(this->dofs_per_face),                
                 component_to_base_table (this->components,
                                          std::make_pair(0U, 0U)),
-                restriction_is_additive_flags(restriction_is_additive_flags),
-               nonzero_components (nonzero_components),
-               n_nonzero_components_table (compute_n_nonzero_components(nonzero_components)),
-               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())
+                restriction_is_additive_flags(r_i_a_f),
+               nonzero_components (nonzero_c)
 {
                                   // Special handling of vectors of
                                   // length one: in this case, we
                                   // assume that all entries were
                                   // supposed to be equal.
+
+                                  // Normally, we should be careful
+                                  // with const_cast, but since this
+                                  // is the constructor and we do it
+                                  // here only, we are fine.
    unsigned int ndofs = this->dofs_per_cell;
-   if (restriction_is_additive_flags.size() == 1 && this->dofs_per_cell > 1)
+   if (restriction_is_additive_flags.size() == 1 && ndofs > 1)
      {
-       std::vector<bool>& riaf = const_cast<std::vector<bool>&>(restriction_is_additive_flags);
-       riaf.resize(ndofs, restriction_is_additive_flags[0]);
+       std::vector<bool>& aux
+        = const_cast<std::vector<bool>&> (restriction_is_additive_flags);
+       aux.resize(ndofs, restriction_is_additive_flags[0]);
      }
    
-//   if (nonzero_components.size() == 1 && this->dofs_per_cell > 1)
-//     nonzero_components.resize(ndofs, nonzero_components[0]);
-  
+   if (nonzero_components.size() == 1 && ndofs > 1)
+     {
+       std::vector<std::vector<bool> >& aux
+        = const_cast<std::vector<std::vector<bool> >&> (nonzero_components);
+       aux.resize(ndofs, nonzero_components[0]);
+     }
+
+                                   // These used to be initialized in
+                                   // the constructor, but here we
+                                   // have the possibly corrected
+                                   // 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();
+   
+   
   Assert (restriction_is_additive_flags.size() == this->dofs_per_cell,
          ExcDimensionMismatch(restriction_is_additive_flags.size(),
                               this->dofs_per_cell));
@@ -639,8 +656,8 @@ compute_2nd (const Mapping<dim>                   &mapping,
 
 template <int dim>
 std::vector<unsigned int>
-FiniteElementBase<dim>::
-compute_n_nonzero_components (const std::vector<std::vector<bool> > &nonzero_components)
+FiniteElementBase<dim>::compute_n_nonzero_components (
+  const std::vector<std::vector<bool> > &nonzero_components)
 {
   std::vector<unsigned int> retval (nonzero_components.size());
   for (unsigned int i=0; i<nonzero_components.size(); ++i)
index 5045dcf58fe1f989786149199c012bb1a9b98eec..d534e4e83f981143731dcdd754e38afe38d43a36 100644 (file)
@@ -185,12 +185,12 @@ FE_Q<dim>::FE_Q (const unsigned int degree)
                :
                FE_Poly<TensorProductPolynomials<dim>, dim> (
                  TensorProductPolynomials<dim>(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)),
-                 FiniteElementData<dim>(get_dpo_vector(degree),1, degree, FiniteElementData<dim>::H1),
-                 std::vector<bool> (FiniteElementData<dim>(
-                   get_dpo_vector(degree),1, degree).dofs_per_cell, false),
-                 std::vector<std::vector<bool> >(FiniteElementData<dim>(
-                   get_dpo_vector(degree),1, degree).dofs_per_cell, std::vector<bool>(1,true))),
-                                                   face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (degree)))
+                 FiniteElementData<dim>(get_dpo_vector(degree),
+                                        1, degree,
+                                        FiniteElementData<dim>::H1),
+                 std::vector<bool> (1, false),
+                 std::vector<std::vector<bool> >(1, std::vector<bool>(1,true))),
+               face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (degree)))
 {
   std::vector<unsigned int> renumber (this->dofs_per_cell);
   FETools::hierarchic_to_lexicographic_numbering (*this, renumber);

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.