From: wolf Date: Mon, 16 Sep 2002 22:12:01 +0000 (+0000) Subject: Get rid of component_to_system_table, and fix a couple of other problems with non... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=42eede34a9c1124f2fd7f0d6c9327e793f340643;p=dealii-svn.git Get rid of component_to_system_table, and fix a couple of other problems with non-primitive elements as arguments to FESystem. git-svn-id: https://svn.dealii.org/trunk@6438 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 220f4ccb40..bbdbf6d1ff 100644 --- a/deal.II/deal.II/include/fe/fe_base.h +++ b/deal.II/deal.II/include/fe/fe_base.h @@ -696,45 +696,6 @@ class FiniteElementBase : public Subscriptor, */ bool operator == (const FiniteElementBase &) const; - /** - * Given a vector component and - * an index of a shape function - * within the shape functions - * corresponding to this vector - * component, return the index of - * this shape function within the - * shape functions of this - * element. If this is a scalar - * element, then the given - * component may only be zero, - * and the given component index - * is also the return value. - * - * If the finite element is - * vector-valued and has - * non-primitive shape functions, - * i.e. some of its shape - * functions are non-zero in more - * than just one vector - * component, then this function - * cannot be used since shape - * functions are no more - * associated with individual - * vector components, and an - * exception of type - * @p{ExcFENotPrimitive} is - * thrown. - */ - unsigned int component_to_system_index (const unsigned int component, - const unsigned int component_index) const; - - /** - * Same as above, but compute the - * data from the index of a shape - * function on a face. - */ - unsigned int face_component_to_system_index (const unsigned int component, - const unsigned int component_index) const; /** * Compute vector component and @@ -1213,37 +1174,6 @@ class FiniteElementBase : public Subscriptor, std::vector,unsigned int> > face_system_to_base_table; - /** - * Map between component and - * linear dofs: For each pair of - * vector component and index - * within this component, store - * the global dof number in the - * composed element. If the - * element is scalar, then the - * outer (component) index can - * only be zero, and the inner - * index is equal to the stored - * value. - * - * If the element is not - * primitive, i.e. there are - * shape functions that are - * non-zero in more than one - * vector-component, then this - * function is obviously useless, - * and all entries will be - * invalid. - */ - std::vector< std::vector > component_to_system_table; - - /** - * Map between component and - * linear dofs on a face. Same - * applies as above. - */ - std::vector< std::vector > face_component_to_system_table; - /** * The base element establishing * a component. @@ -1542,23 +1472,6 @@ FiniteElementData::n_components () const }; -template -inline -unsigned int -FiniteElementBase::component_to_system_index (const unsigned int component, - const unsigned int component_index) const -{ - Assert(componentn_components(), - ExcIndexRange(component, 0, this->n_components())); - Assert(component_index::ExcFENotPrimitive()); - - return component_to_system_table[component][component_index]; -} - template inline @@ -1573,22 +1486,6 @@ FiniteElementBase::system_to_component_index (const unsigned int index) con } -template -inline -unsigned int -FiniteElementBase::face_component_to_system_index (const unsigned int component, - const unsigned int component_index) const -{ - Assert(componentn_components(), - ExcIndexRange(component, 0, this->n_components())); - Assert(component_index::ExcFENotPrimitive()); - - return face_component_to_system_table[component][component_index]; -} template diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 9c88dfd678..cea40f5b4a 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -1865,8 +1865,14 @@ DoFTools::compute_intergrid_weights_1 (const DoFHandler &coars for (unsigned int local_coarse_dof=0; local_coarse_dof &fe_data, face_system_to_component_table(this->dofs_per_face), system_to_base_table(this->dofs_per_cell), face_system_to_base_table(this->dofs_per_face), - component_to_system_table(this->components, - std::vector(this->dofs_per_cell)), - face_component_to_system_table(this->components, - std::vector(this->dofs_per_face)), component_to_base_table (this->components, std::make_pair(0U, 0U)), restriction_is_additive_flags(restriction_is_additive_flags), @@ -198,13 +194,11 @@ FiniteElementBase (const FiniteElementData &fe_data, { system_to_component_table[j] = std::pair(0,j); system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j); - component_to_system_table[0][j] = j; } for (unsigned int j=0 ; jdofs_per_face ; ++j) { face_system_to_component_table[j] = std::pair(0,j); face_system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j); - face_component_to_system_table[0][j] = j; } }; @@ -394,8 +388,6 @@ FiniteElementBase::memory_consumption () const MemoryConsumption::memory_consumption (face_system_to_component_table) + MemoryConsumption::memory_consumption (system_to_base_table) + MemoryConsumption::memory_consumption (face_system_to_base_table) + - MemoryConsumption::memory_consumption (component_to_system_table) + - MemoryConsumption::memory_consumption (face_component_to_system_table) + MemoryConsumption::memory_consumption (component_to_base_table) + MemoryConsumption::memory_consumption (restriction_is_additive_flags) + MemoryConsumption::memory_consumption (nonzero_components) + diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 67a339c228..e2a43b3c29 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -655,42 +655,17 @@ FESystem::compute_fill (const Mapping &mapping, Assert (face_quadrature != 0, ExcInternalError()); } - - for (unsigned int base_no=0, comp=0; base_no &base_fe=base_element(base_no); + const FiniteElement & + base_fe = base_element(base_no); typename FiniteElementBase::InternalDataBase & base_fe_data = fe_data.get_fe_data(base_no); - - // Make sure that in the - // case of fill_fe_values - // the data is only copied - // from base_data to data - // if base_data is - // changed. therefore use - // fe_fe_data.current_update_flags() - - // for the case of - // fill_fe_(sub)face_values - // the data needs to be - // copied from base_data to - // data on each face, - // therefore use - // base_fe_data.update_flags. - - // Store these flags into - // base_flags before - // calling - // base_fe.fill_fe_([sub]face_)values - // as the latter changes - // the return value of - // base_fe_data.current_update_flags() - const UpdateFlags base_flags(dim_1==dim ? - base_fe_data.current_update_flags() : - base_fe_data.update_flags); - - FEValuesData &base_data=fe_data.get_fe_values_data(base_no); - + FEValuesData & + base_data = fe_data.get_fe_values_data(base_no); + if (face_no==invalid_face_number) base_fe.fill_fe_values(mapping, cell, *cell_quadrature, mapping_data, base_fe_data, base_data); @@ -700,25 +675,133 @@ FESystem::compute_fill (const Mapping &mapping, else base_fe.fill_fe_subface_values(mapping, cell, face_no, sub_no, *face_quadrature, mapping_data, base_fe_data, base_data); + }; - for (unsigned int m=0; mcomponent_to_system_index(comp,k); - if (base_flags & update_values) - data.shape_values(system_index, point)= - base_data.shape_values(k,point); - - if (base_flags & update_gradients) - data.shape_gradients[system_index][point]= - base_data.shape_gradients[k][point]; - - if (base_flags & update_second_derivatives) - data.shape_2nd_derivatives[system_index][point]= - base_data.shape_2nd_derivatives[k][point]; - } - } + // now data has been generated, + // so copy it. we used to work + // by looping over all base + // elements, then over + // multiplicity, then over the + // shape functions from that + // base element, but that + // requires that we can infer + // the global number of a shape + // function from its number in + // the base element. for that + // we had the + // component_to_system_table. + // + // however, this does of course + // no longer work since we have + // non-primitive elements. so + // we go the other way round: + // loop over all shape + // functions of the composed + // element, and find from where + // to copy the data. to be able + // to cache some data, make + // things a little bit more + // complicated: loop over all + // base elements, then over all + // shape functions of the + // composed element, and only + // treat those shape functions + // that belong to a given base + // element + for (unsigned int base_no=0; base_nodofs_per_cell; + ++system_index) + if (this->system_to_base_table[system_index].first.first == base_no) + { + const FiniteElement & + base_fe = base_element(base_no); + typename FiniteElementBase::InternalDataBase & + base_fe_data = fe_data.get_fe_data(base_no); + FEValuesData & + base_data = fe_data.get_fe_values_data(base_no); + + const unsigned int n_q_points = quadrature.n_quadrature_points; + + // Make sure that in + // the case of + // fill_fe_values the + // data is only copied + // from base_data to + // data if base_data is + // changed. therefore + // use + // fe_fe_data.current_update_flags() + // + // for the case of + // fill_fe_(sub)face_values + // the data needs to be + // copied from + // base_data to data on + // each face, therefore + // use + // base_fe_data.update_flags. + // + // Store these flags into + // base_flags before + // calling + // base_fe.fill_fe_([sub]face_)values + // as the latter changes + // the return value of + // base_fe_data.current_update_flags() + const UpdateFlags base_flags(dim_1==dim ? + base_fe_data.current_update_flags() : + base_fe_data.update_flags); + + const unsigned int k = system_to_base_table[system_index].second; + Assert (kn_nonzero_components(i); + unsigned int in_index = 0; + for (unsigned int i=0; in_nonzero_components(system_index) == + base_fe.n_nonzero_components(k), + ExcInternalError()); + for (unsigned int s=0; sn_nonzero_components(system_index); ++s) + { + if (base_flags & update_values) + for (unsigned int q=0; q::build_cell_tables() this->system_to_component_table[total_index] = non_primitive_index; } } - - // Initialize mapping from - // components to linear - // index. Fortunately, this is the - // inverse of what we just did. - // - // note that we can only do this if - // all sub-elements are primiive, - // otherwise fill the table with - // invalid values - if (this->is_primitive()) - { - std::vector dofs_per_component (this->n_components(), 0); - for (unsigned int sys=0; sysdofs_per_cell; ++sys) - ++dofs_per_component[this->system_to_component_index(sys).first]; - for (unsigned int component=0; componentn_components(); ++component) - this->component_to_system_table[component].resize(dofs_per_component[component]); - - // then go the reverse way to fill the array - for (unsigned int sys=0; sysdofs_per_cell; ++sys) - { - const unsigned int - comp = this->system_to_component_index(sys).first, - index_in_comp = this->system_to_component_index(sys).second; - - Assert (comp < this->component_to_system_table.size(), - ExcInternalError()); - Assert (index_in_comp < this->component_to_system_table[comp].size(), - ExcInternalError()); - this->component_to_system_table[comp][index_in_comp] = sys; - }; - } - else - // element is not primitive, so - // fill elements with nonsense - { - std::vector< std::vector > tmp; - this->component_to_system_table.swap (tmp); - }; }; @@ -1159,40 +1203,6 @@ FESystem::build_face_tables() ExcInternalError()); Assert (total_index == this->face_system_to_base_table.size(), ExcInternalError()); - - // finally, initialize reverse - // mapping. same applies - // w.r.t. primitivity as above: - if (this->is_primitive()) - { - std::vector dofs_per_component (this->n_components(), 0); - for (unsigned int sys=0; sysdofs_per_face; ++sys) - ++dofs_per_component[this->face_system_to_component_index(sys).first]; - for (unsigned int component=0; componentn_components(); ++component) - this->face_component_to_system_table[component].resize(dofs_per_component[component]); - - // then go the reverse way to fill - // the array - for (unsigned int sys=0; sysdofs_per_face; ++sys) - { - const unsigned int - comp = this->face_system_to_component_index(sys).first, - index_in_comp = this->face_system_to_component_index(sys).second; - - Assert (comp < this->face_component_to_system_table.size(), - ExcInternalError()); - Assert (index_in_comp < this->face_component_to_system_table[comp].size(), - ExcInternalError()); - this->face_component_to_system_table[comp][index_in_comp] = sys; - }; - } - else - // element is not primitive, so - // fill elements with nonsense - { - std::vector< std::vector > tmp; - this->face_component_to_system_table.swap (tmp); - }; }; diff --git a/tests/fe/up_and_down.cc b/tests/fe/up_and_down.cc index e2d2654b99..e3cc52ea27 100644 --- a/tests/fe/up_and_down.cc +++ b/tests/fe/up_and_down.cc @@ -206,19 +206,19 @@ void test () // with scalar and Nedelec // elements, to make things // really difficult -// (dim != 1 ? -// new FESystem(FE_Nedelec(1), 2) -// : 0), -// (dim != 1 ? -// new FESystem(FE_Nedelec(1), 2, -// FE_Q (2), 2) -// : 0), -// (dim != 1 ? -// new FESystem(FE_Nedelec(1), 2, -// FE_DGQ (2), 2, -// FESystem(FE_Nedelec(1), 2, -// FE_Q (2), 2), 2) -// : 0), + (dim != 1 ? + new FESystem(FE_Nedelec(1), 2) + : 0), + (dim != 1 ? + new FESystem(FE_Nedelec(1), 2, + FE_Q (2), 2) + : 0), + (dim != 1 ? + new FESystem(FE_Nedelec(1), 2, + FE_DGQ (2), 2, + FESystem(FE_Nedelec(1), 2, + FE_Q (2), 2), 2) + : 0), }; for (unsigned int i=0; i