From 8fed6a501d48d8df1ebe4e9b7e287fc5540c2da4 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Fri, 4 Apr 2014 12:38:36 +0000 Subject: [PATCH] Forgot one directory git-svn-id: https://svn.dealii.org/trunk@32715 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/source/fe/fe.cc | 6 ++-- deal.II/source/fe/fe_dgp.cc | 5 +-- deal.II/source/fe/fe_dgq.cc | 8 ++--- deal.II/source/fe/fe_face.cc | 15 +++++---- deal.II/source/fe/fe_nedelec.cc | 8 +++-- deal.II/source/fe/fe_q_base.cc | 11 +++--- deal.II/source/fe/fe_q_dg0.cc | 20 +++++++++++ deal.II/source/fe/fe_q_hierarchical.cc | 5 +-- deal.II/source/fe/fe_raviart_thomas.cc | 8 +++-- deal.II/source/fe/fe_system.cc | 46 ++++++++++++++++++++------ 10 files changed, 97 insertions(+), 35 deletions(-) diff --git a/deal.II/source/fe/fe.cc b/deal.II/source/fe/fe.cc index 8237022f8a..7e9d79bf4d 100644 --- a/deal.II/source/fe/fe.cc +++ b/deal.II/source/fe/fe.cc @@ -1104,11 +1104,13 @@ FiniteElement::has_support_on_face ( template -Table<2,bool> +std::pair, std::vector > FiniteElement::get_constant_modes () const { Assert (false, ExcNotImplemented()); - return Table<2,bool>(this->n_components(), this->dofs_per_cell); + return std::pair, std::vector > + (Table<2,bool>(this->n_components(), this->dofs_per_cell), + std::vector(this->n_components())); } diff --git a/deal.II/source/fe/fe_dgp.cc b/deal.II/source/fe/fe_dgp.cc index ca52e4f149..daa2b51d3a 100644 --- a/deal.II/source/fe/fe_dgp.cc +++ b/deal.II/source/fe/fe_dgp.cc @@ -238,12 +238,13 @@ FE_DGP::has_support_on_face (const unsigned int, template -Table<2,bool> +std::pair, std::vector > FE_DGP::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); constant_modes(0,0) = true; - return constant_modes; + return std::pair, std::vector > + (constant_modes, std::vector(1, 0)); } diff --git a/deal.II/source/fe/fe_dgq.cc b/deal.II/source/fe/fe_dgq.cc index 9f8a69ce29..8fe34b6c60 100644 --- a/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/source/fe/fe_dgq.cc @@ -765,13 +765,13 @@ FE_DGQ::has_support_on_face (const unsigned int shape_index, template -Table<2,bool> +std::pair, std::vector > FE_DGQ::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); - for (unsigned int i=0; idofs_per_cell; ++i) - constant_modes(0,i) = true; - return constant_modes; + constant_modes.fill(true); + return std::pair, std::vector > + (constant_modes, std::vector(1, 0)); } diff --git a/deal.II/source/fe/fe_face.cc b/deal.II/source/fe/fe_face.cc index 9d62ff19d5..3a210bd849 100644 --- a/deal.II/source/fe/fe_face.cc +++ b/deal.II/source/fe/fe_face.cc @@ -289,13 +289,14 @@ compare_for_face_domination (const FiniteElement &fe_other) const template -Table<2,bool> +std::pair, std::vector > FE_FaceQ::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); for (unsigned int i=0; idofs_per_cell; ++i) constant_modes(0,i) = true; - return constant_modes; + return std::pair, std::vector > + (constant_modes, std::vector(1, 0)); } @@ -416,13 +417,14 @@ compare_for_face_domination (const FiniteElement<1,spacedim> &fe_other) const template -Table<2,bool> +std::pair, std::vector > FE_FaceQ<1,spacedim>::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); for (unsigned int i=0; idofs_per_cell; ++i) constant_modes(0,i) = true; - return constant_modes; + return std::pair, std::vector > + (constant_modes, std::vector(1,0)); } @@ -780,13 +782,14 @@ get_subface_interpolation_matrix (const FiniteElement &x_source_fe template -Table<2,bool> +std::pair, std::vector > FE_FaceP::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); for (unsigned int face=0; face::faces_per_cell; ++face) constant_modes(0, face*this->dofs_per_face) = true; - return constant_modes; + return std::pair, std::vector > + (constant_modes, std::vector(1, 0)); } diff --git a/deal.II/source/fe/fe_nedelec.cc b/deal.II/source/fe/fe_nedelec.cc index 29c15ced94..40c14953f4 100644 --- a/deal.II/source/fe/fe_nedelec.cc +++ b/deal.II/source/fe/fe_nedelec.cc @@ -5466,14 +5466,18 @@ const template -Table<2,bool> +std::pair, std::vector > FE_Nedelec::get_constant_modes() const { Table<2,bool> constant_modes(dim, this->dofs_per_cell); for (unsigned int d=0; ddofs_per_cell; ++i) constant_modes(d,i) = true; - return constant_modes; + std::vector components; + for (unsigned int d=0; d, std::vector > + (constant_modes, components); } diff --git a/deal.II/source/fe/fe_q_base.cc b/deal.II/source/fe/fe_q_base.cc index 7212212063..0e9b94ec33 100644 --- a/deal.II/source/fe/fe_q_base.cc +++ b/deal.II/source/fe/fe_q_base.cc @@ -1495,14 +1495,15 @@ FE_Q_Base::has_support_on_face (const unsigned int shape_inde template -Table<2,bool> +std::pair, std::vector > FE_Q_Base::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); - // leave out last component - for (unsigned int i=0; i(this->degree+1); ++i) - constant_modes(0, i) = true; - return constant_modes; + // FE_Q_DG0 should not use this function... + AssertDimension(this->dofs_per_cell, Utilities::fixed_power(this->degree+1)); + constant_modes.fill(true); + return std::pair, std::vector > + (constant_modes, std::vector(1, 0)); } diff --git a/deal.II/source/fe/fe_q_dg0.cc b/deal.II/source/fe/fe_q_dg0.cc index 8ecc7df865..99ba0e946b 100644 --- a/deal.II/source/fe/fe_q_dg0.cc +++ b/deal.II/source/fe/fe_q_dg0.cc @@ -307,6 +307,26 @@ FE_Q_DG0::has_support_on_face (const unsigned int shape_index, } + +template +std::pair, std::vector > +FE_Q_DG0::get_constant_modes () const +{ + Table<2,bool> constant_modes(2, this->dofs_per_cell); + + // 1 represented by FE_Q part + for (unsigned int i=0; idofs_per_cell-1; ++i) + constant_modes(0, i) = true; + + // 1 represented by DG0 part + constant_modes(1, this->dofs_per_cell-1) = true; + + return std::pair, std::vector > + (constant_modes, std::vector (2, 0)); +} + + + // explicit instantiations #include "fe_q_dg0.inst" diff --git a/deal.II/source/fe/fe_q_hierarchical.cc b/deal.II/source/fe/fe_q_hierarchical.cc index b1800b6311..8d2d5bf0c6 100644 --- a/deal.II/source/fe/fe_q_hierarchical.cc +++ b/deal.II/source/fe/fe_q_hierarchical.cc @@ -1862,7 +1862,7 @@ FE_Q_Hierarchical::get_embedding_dofs (const unsigned int sub_degree) const template -Table<2,bool> +std::pair, std::vector > FE_Q_Hierarchical::get_constant_modes () const { Table<2,bool> constant_modes(1, this->dofs_per_cell); @@ -1870,7 +1870,8 @@ FE_Q_Hierarchical::get_constant_modes () const constant_modes(0,i) = true; for (unsigned int i=GeometryInfo::vertices_per_cell; idofs_per_cell; ++i) constant_modes(0,i) = false; - return constant_modes; + return std::pair, std::vector > + (constant_modes, std::vector(1, 0)); } diff --git a/deal.II/source/fe/fe_raviart_thomas.cc b/deal.II/source/fe/fe_raviart_thomas.cc index 22208bae95..c799221bc7 100644 --- a/deal.II/source/fe/fe_raviart_thomas.cc +++ b/deal.II/source/fe/fe_raviart_thomas.cc @@ -404,14 +404,18 @@ FE_RaviartThomas::get_dpo_vector (const unsigned int deg) template -Table<2,bool> +std::pair, std::vector > FE_RaviartThomas::get_constant_modes() const { Table<2,bool> constant_modes(dim, this->dofs_per_cell); for (unsigned int d=0; ddofs_per_cell; ++i) constant_modes(d,i) = true; - return constant_modes; + std::vector components; + for (unsigned int d=0; d, std::vector > + (constant_modes, components); } diff --git a/deal.II/source/fe/fe_system.cc b/deal.II/source/fe/fe_system.cc index 24cc2213c5..bf072c4a01 100644 --- a/deal.II/source/fe/fe_system.cc +++ b/deal.II/source/fe/fe_system.cc @@ -2983,28 +2983,54 @@ FESystem::unit_face_support_point (const unsigned int index) const template -Table<2,bool> +std::pair, std::vector > FESystem::get_constant_modes () const { + // Note that this->n_components() is actually only an estimate of how many + // constant modes we will need. There might be more than one such mode + // (e.g. FE_Q_DG0). Table<2,bool> constant_modes(this->n_components(), this->dofs_per_cell); - unsigned int comp=0; + std::vector components; for (unsigned int i=0; i base_table = base_elements[i].first->get_constant_modes(); - const unsigned int n_base_components = base_elements[i].first->n_components(); + std::pair, std::vector > + base_table = base_elements[i].first->get_constant_modes(); + AssertDimension(base_table.first.n_rows(), base_table.second.size()); + const unsigned int element_multiplicity = this->element_multiplicity(i); + + // there might be more than one constant mode for some scalar elements, + // so make sure the table actually fits: Create a new table with more + // rows + const unsigned int comp = components.size(); + if (constant_modes.n_rows() < comp+base_table.first.n_rows()*element_multiplicity) + { + Table<2,bool> new_constant_modes(comp+base_table.first.n_rows()* + element_multiplicity, + constant_modes.n_cols()); + std::copy(&constant_modes(0,0), &constant_modes(0,0)+this->dofs_per_cell*comp, + &new_constant_modes(0,0)); + constant_modes.swap(new_constant_modes); + } + + // next, fill the constant modes from the individual components as well + // as the component numbers corresponding to the constant mode rows for (unsigned int k=0; kdofs_per_cell; ++k) { std::pair, unsigned int> ind = this->system_to_base_index(k); if (ind.first.first == i) - for (unsigned int c=0; cbase_elements[i].first->n_components() + +base_table.second[c]); } - AssertDimension(comp, this->n_components()); - return constant_modes; + AssertDimension(components.size(), constant_modes.n_rows()); + return std::pair, std::vector >(constant_modes, + components); } -- 2.39.5