From: guido Date: Wed, 29 Jun 2005 15:41:20 +0000 (+0000) Subject: simplified constructor for FiniteElementBase X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=14cde46fcbd87b97f978c0da382151f59d6b0eac;p=dealii-svn.git simplified constructor for FiniteElementBase git-svn-id: https://svn.dealii.org/trunk@10970 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_base.h b/deal.II/deal.II/include/fe/fe_base.h index 9f39c08923..82dae1329b 100644 --- a/deal.II/deal.II/include/fe/fe_base.h +++ b/deal.II/deal.II/include/fe/fe_base.h @@ -362,7 +362,7 @@ class FiniteElementData * transfer, respectively. * * - * @sect3{Support points} + *

Support points

* * 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} + *

Finite elements in one dimension

* * 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} + *

Finite elements in two dimensions

* * 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} + *

Finite elements in three dimensions

* * 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 + * dofs_per_cell. 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 &fe_data, const std::vector &restriction_is_additive_flags, diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 7b973f5dfe..994005544a 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -106,10 +106,10 @@ FiniteElementBase::InternalDataBase::~InternalDataBase () template -FiniteElementBase:: -FiniteElementBase (const FiniteElementData &fe_data, - const std::vector &restriction_is_additive_flags, - const std::vector > &nonzero_components) +FiniteElementBase::FiniteElementBase ( + const FiniteElementData &fe_data, + const std::vector &r_i_a_f, + const std::vector > &nonzero_c) : FiniteElementData (fe_data), system_to_component_table (this->dofs_per_cell), @@ -118,30 +118,47 @@ FiniteElementBase (const FiniteElementData &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(), - 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& riaf = const_cast&>(restriction_is_additive_flags); - riaf.resize(ndofs, restriction_is_additive_flags[0]); + std::vector& aux + = const_cast&> (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 >& aux + = const_cast >&> (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&> + (n_nonzero_components_table) = compute_n_nonzero_components(nonzero_components); + const_cast + (cached_primitivity) = std::find_if (n_nonzero_components_table.begin(), + n_nonzero_components_table.end(), + std::bind2nd(std::not_equal_to(), + 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 &mapping, template std::vector -FiniteElementBase:: -compute_n_nonzero_components (const std::vector > &nonzero_components) +FiniteElementBase::compute_n_nonzero_components ( + const std::vector > &nonzero_components) { std::vector retval (nonzero_components.size()); for (unsigned int i=0; i::FE_Q (const unsigned int degree) : FE_Poly, dim> ( TensorProductPolynomials(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)), - FiniteElementData(get_dpo_vector(degree),1, degree, FiniteElementData::H1), - std::vector (FiniteElementData( - get_dpo_vector(degree),1, degree).dofs_per_cell, false), - std::vector >(FiniteElementData( - get_dpo_vector(degree),1, degree).dofs_per_cell, std::vector(1,true))), - face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (degree))) + FiniteElementData(get_dpo_vector(degree), + 1, degree, + FiniteElementData::H1), + std::vector (1, false), + std::vector >(1, std::vector(1,true))), + face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (degree))) { std::vector renumber (this->dofs_per_cell); FETools::hierarchic_to_lexicographic_numbering (*this, renumber);