From c1af58fafd4e294d883547f6f4d34d83a622d2ab Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 7 Mar 2001 17:53:20 +0000 Subject: [PATCH] Further clean-ups. git-svn-id: https://svn.dealii.org/trunk@4157 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_system.h | 38 +++-- deal.II/deal.II/source/fe/fe_system.cc | 212 +++++++++++++++---------- 2 files changed, 149 insertions(+), 101 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_system.h b/deal.II/deal.II/include/fe/fe_system.h index ca4a9249f1..c5eef93c96 100644 --- a/deal.II/deal.II/include/fe/fe_system.h +++ b/deal.II/deal.II/include/fe/fe_system.h @@ -202,11 +202,13 @@ class FESystem : public FiniteElement * the @p{cell->get_dof_indices} * function, but: * - * If the shape functions of one - * of the base elements are not - * Lagrangian interpolants at some - * points, the size of the array - * will be zero. + * If one of the base elements + * has no support points, then it + * makes no sense to define + * support points for the + * composed element, so return an + * empty array to demonstrate + * that fact. */ virtual void get_unit_face_support_points (typename std::vector > &) const; @@ -311,6 +313,13 @@ class FESystem : public FiniteElement private: + /** + * Value to indicate that a given + * face or subface number is + * invalid. + */ + static const unsigned int invalid_face_number = static_cast(-1); + /** * Pairs of multiplicity and * element type. @@ -404,17 +413,14 @@ class FESystem : public FiniteElement const unsigned int N3); /** - * This function is simply singled out of - * the constructor. It sets up the - * index table for the system as well as - * @p{restriction} and @p{prolongation} - * matrices. Since the operation of this - * function can be done without explicit - * knowledge of the data type of the - * underlying finite element class, we - * don't want to have this function in - * the general template definition in - * the @p{.h} file. + * This function is simply + * singled out of the + * constructors since there are + * several of them. It sets up + * the index table for the system + * as well as @p{restriction} and + * @p{prolongation} + * matrices. */ void initialize(); diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 8dec1d7168..140c00d0ab 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -31,6 +31,36 @@ using namespace std; /* ----------------------- FESystem::InternalData ------------------- */ +template +FESystem::InternalData::InternalData(const unsigned int n_base_elements): + base_fe_datas(n_base_elements), + base_fe_values_datas(n_base_elements) +{} + + + +template +FESystem::InternalData::~InternalData() +{ + // delete pointers and set them to + // zero to avoid inadvertant use + for (unsigned int i=0; i inline FiniteElementBase::InternalDataBase & FESystem:: @@ -99,6 +129,10 @@ InternalData::delete_fe_values_data (const unsigned int base_no) /* ---------------------------------- FESystem ------------------- */ +template +const unsigned int FESystem::invalid_face_number; + + template FESystem::FESystem (const FiniteElement &fe, const unsigned int n_elements) : FiniteElement (multiply_dof_numbers(fe, n_elements), @@ -214,6 +248,7 @@ get_unit_support_points (typename std::vector > &unit_support_points) base_element_dofs_per_cell = base_element(base_el).dofs_per_cell; base_element(base_el).get_unit_support_points (base_unit_support_points); + // if one of the base elements // has no support points, then // it makes no sense to define @@ -228,9 +263,10 @@ get_unit_support_points (typename std::vector > &unit_support_points) } // otherwise distribute the - // support of this base element - // to all degrees of freedom - // contributed by it + // support points of this base + // element to all degrees of + // freedom contributed by this + // base element Assert(base_unit_support_points.size()==base_element_dofs_per_cell, ExcInternalError()); for (unsigned int n=0; n > &unit_support_points) template -void FESystem::get_unit_face_support_points ( - typename std::vector > &unit_support_points) const +void +FESystem:: +get_unit_face_support_points (typename std::vector > &unit_support_points) const { + // generate unit face support points + // from unit support points of sub + // elements unit_support_points.resize(dofs_per_face); - typename std::vector > base_unit_support_points (base_element(0).dofs_per_cell); + typename std::vector > + base_unit_support_points (base_element(0).dofs_per_cell); + unsigned int comp = 0; for (unsigned int base_el=0 ; base_el::get_unit_face_support_points ( } + template UpdateFlags -FESystem::update_once (UpdateFlags flags) const +FESystem::update_once (const UpdateFlags flags) const { UpdateFlags out = update_default; + // generate maximal set of flags + // that are necessary for (unsigned int base_no=0; base_no UpdateFlags -FESystem::update_each (UpdateFlags flags) const +FESystem::update_each (const UpdateFlags flags) const { UpdateFlags out = update_default; + // generate maximal set of flags + // that are necessary for (unsigned int base_no=0; base_no Mapping::InternalDataBase* -FESystem::get_data ( UpdateFlags flags, +FESystem::get_data (UpdateFlags flags, const Mapping &mapping, const Quadrature &quadrature) const { @@ -318,7 +377,9 @@ FESystem::get_data ( UpdateFlags flags, data->initialize_2nd (this, mapping, quadrature); } - + + // get data objects from each of + // the base elements and store them for (unsigned int base_no=0; base_no::InternalDataBase *base_fe_data_base = @@ -363,9 +424,11 @@ FESystem::get_data ( UpdateFlags flags, FEValuesData *base_data = new FEValuesData(); data->set_fe_values_data(base_no, base_data); } + data->update_flags=data->update_once | data->update_each; Assert(data->update_once==update_once(flags), ExcInternalError()); Assert(data->update_each==update_each(flags), ExcInternalError()); + return data; } @@ -376,14 +439,14 @@ void FESystem::fill_fe_values (const Mapping &mapping, const DoFHandler::cell_iterator &cell, const Quadrature &quadrature, - Mapping::InternalDataBase &mapping_data, - Mapping::InternalDataBase &fe_data, + Mapping::InternalDataBase &mapping_data, + Mapping::InternalDataBase &fe_data, FEValuesData &data) const { - const unsigned int minus_1=static_cast (-1); - compute_fill(mapping, cell, minus_1, minus_1, + compute_fill(mapping, cell, invalid_face_number, invalid_face_number, quadrature, mapping_data, fe_data, data); -} +}; + template @@ -392,15 +455,14 @@ FESystem::fill_fe_face_values (const Mapping &mappin const DoFHandler::cell_iterator &cell, const unsigned int face_no, const Quadrature &quadrature, - Mapping::InternalDataBase &mapping_data, - Mapping::InternalDataBase &fe_data, + Mapping::InternalDataBase &mapping_data, + Mapping::InternalDataBase &fe_data, FEValuesData &data) const { - const unsigned int minus_1=static_cast (-1); - compute_fill(mapping, cell, face_no, minus_1, + compute_fill(mapping, cell, face_no, invalid_face_number, quadrature, mapping_data, fe_data, data); +}; -} @@ -411,8 +473,8 @@ FESystem::fill_fe_subface_values (const Mapping &map const unsigned int face_no, const unsigned int sub_no, const Quadrature &quadrature, - Mapping::InternalDataBase &mapping_data, - Mapping::InternalDataBase &fe_data, + Mapping::InternalDataBase &mapping_data, + Mapping::InternalDataBase &fe_data, FEValuesData &data) const { compute_fill(mapping, cell, face_no, sub_no, @@ -429,11 +491,15 @@ FESystem::compute_fill (const Mapping &mapping, const unsigned int face_no, const unsigned int sub_no, const Quadrature &quadrature, - Mapping::InternalDataBase &mapping_data, - Mapping::InternalDataBase &fedata, + Mapping::InternalDataBase &mapping_data, + Mapping::InternalDataBase &fedata, FEValuesData &data) const { - InternalData& fe_data = dynamic_cast (fedata); + // convert data object to internal + // data for this class. fails with + // an exception if that is not + // possible + InternalData & fe_data = dynamic_cast (fedata); // Either dim_1==dim (fill_fe_values) // or dim_1==dim-1 (fill_fe_(sub)face_values) @@ -494,8 +560,7 @@ FESystem::compute_fill (const Mapping &mapping, // Quadrature: const Subscriptor* quadrature_base_pointer = &quadrature; - const unsigned int minus_1=static_cast (-1); - if (face_no==minus_1) + if (face_no==invalid_face_number) { Assert(dim_1==dim, ExcDimensionMismatch(dim_1,dim)); cell_quadrature=dynamic_cast *>(quadrature_base_pointer); @@ -544,10 +609,10 @@ FESystem::compute_fill (const Mapping &mapping, FEValuesData &base_data=fe_data.get_fe_values_data(base_no); - if (face_no==minus_1) + if (face_no==invalid_face_number) base_fe.fill_fe_values(mapping, cell, *cell_quadrature, mapping_data, base_fe_data, base_data); - else if (sub_no==minus_1) + else if (sub_no==invalid_face_number) base_fe.fill_fe_face_values(mapping, cell, face_no, *face_quadrature, mapping_data, base_fe_data, base_data); else @@ -601,11 +666,12 @@ FESystem::compute_fill (const Mapping &mapping, } } } + if (fe_data.compute_second_derivatives) { unsigned int offset = 0; - if (face_no != static_cast (-1)) - offset = (sub_no == static_cast (-1)) + if (face_no != invalid_face_number) + offset = (sub_no == invalid_face_number) ? face_no * quadrature.n_quadrature_points :(face_no * GeometryInfo::subfaces_per_face + sub_no) * quadrature.n_quadrature_points; @@ -732,7 +798,7 @@ FESystem::build_cell_table() } } - // Inintialize mapping from component + // Initialize mapping from component // to base element // Initialize mapping from components to // linear index. Fortunately, this is @@ -747,6 +813,7 @@ FESystem::build_cell_table() } + template void FESystem::build_face_table() @@ -843,7 +910,7 @@ FESystem::build_face_table() Assert (total_index <= face_system_to_component_table.size(), ExcInternalError()); - // Inintialize mapping from component + // Initialize mapping from component // to base element // Initialize mapping from components to // linear index. Fortunately, this is @@ -858,6 +925,7 @@ FESystem::build_face_table() } + template void FESystem::build_interface_constraints () { @@ -1052,18 +1120,19 @@ void FESystem::build_interface_constraints () }; + template void FESystem::initialize () { build_cell_table(); build_face_table(); - // Check if matrices are void. - + // Check if some of the matrices of + // the base elements are void. bool do_restriction = true; bool do_prolongation = true; - for (unsigned int i=0;i::initialize () do_prolongation = false; } - // Void matrices if not defined for all base elements - + // if we encountered void matrices, + // disable them for the composite + // element as well if (!do_restriction) for (unsigned int i=0;i::children_per_cell;++i) restriction[i].reinit(0,0); @@ -1107,10 +1177,9 @@ void FESystem::initialize () // this is kind'o hairy, so don't try // to do it dimension independent - // TODO: there's an assertion thrown for - // dim=3 and for FESystem(FE_Q (3), 2) and for - // FESystem(FE_Q (1), 2, FE_Q (3), 1) - // and for FESystem(FE_Q (4), 2)) +//TODO: there's an assertion thrown for dim=3 and for +//FESystem(FE_Q (3), 2) and for FESystem(FE_Q (1), 2, +//FE_Q (3), 1) and for FESystem(FE_Q (4), 2)) build_interface_constraints (); }; @@ -1196,6 +1265,7 @@ FESystem::compute_restriction_is_additive_flags (const FiniteElement & }; + template std::vector FESystem::compute_restriction_is_additive_flags (const FiniteElement &fe1, @@ -1214,6 +1284,7 @@ FESystem::compute_restriction_is_additive_flags (const FiniteElement & }; + template std::vector FESystem::compute_restriction_is_additive_flags (const FiniteElement &fe1, @@ -1237,6 +1308,7 @@ FESystem::compute_restriction_is_additive_flags (const FiniteElement & }; + template unsigned int FESystem::memory_consumption () const @@ -1248,43 +1320,13 @@ FESystem::memory_consumption () const // after all, considering the size // of the data of the subelements unsigned int mem = (FiniteElement::memory_consumption () + - sizeof (base_elements)); + sizeof (base_elements)); for (unsigned int i=0; i -FESystem::InternalData::InternalData(const unsigned int n_base_elements): - base_fe_datas(n_base_elements), - base_fe_values_datas(n_base_elements) -{} - - - -template -FESystem::InternalData::~InternalData() -{ - // delete pointers and set them to - // zero to avoid inadvertant use - for (unsigned int i=0; i