From 44bf493ea386de6e8a1f5ba7c27ec2d7918a0b2d Mon Sep 17 00:00:00 2001 From: wolf Date: Fri, 30 Mar 2001 09:20:19 +0000 Subject: [PATCH] Unify lots of code from the 6 constructors of the FE*Values classes. git-svn-id: https://svn.dealii.org/trunk@4348 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe.h | 2 + deal.II/deal.II/include/fe/fe_values.h | 44 +++++- deal.II/deal.II/source/fe/fe_values.cc | 190 +++++++++---------------- 3 files changed, 112 insertions(+), 124 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index 6e4edfc9f7..dcc57eb7a4 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -18,6 +18,7 @@ #include template class FEValuesData; +template class FEValuesBase; template class FEValues; template class FEFaceValues; template class FESubfaceValues; @@ -256,6 +257,7 @@ class FiniteElement : public FiniteElementBase * Declare some other classes as * friends of this class. */ + friend class FEValuesBase; friend class FEValues; friend class FEFaceValues; friend class FESubfaceValues; diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index 554c26ec77..cbbfdc05d8 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -615,7 +615,25 @@ class FEValuesBase : protected FEValuesData * Internal data of finite element. */ SmartPointer::InternalDataBase> fe_data; - + + /** + * Initialize some update + * flags. Called from the + * @p{initialize} functions of + * derived classes, which are in + * turn called from their + * constructors. + * + * Basically, this function finds + * out using the finite element + * and mapping object already + * stored which flags need to be + * set to compute everything the + * user wants, as expressed + * through the flags passed as + * argument. + */ + UpdateFlags compute_update_flags (const UpdateFlags update_flags) const; }; @@ -671,7 +689,13 @@ class FEValues : public FEValuesBase /** * Store the quadrature formula here. */ - const Quadrature quadrature; + const Quadrature quadrature; + + /** + * Do work common to the two + * constructors. + */ + void initialize (const UpdateFlags update_flags); }; @@ -825,6 +849,14 @@ class FEFaceValues : public FEFaceValuesBase */ void reinit (const typename DoFHandler::cell_iterator &cell, const unsigned int face_no); + + private: + + /** + * Do work common to the two + * constructors. + */ + void initialize (const UpdateFlags update_flags); }; @@ -888,6 +920,14 @@ class FESubfaceValues : public FEFaceValuesBase * Exception */ DeclException0 (ExcFaceHasNoSubfaces); + + private: + + /** + * Do work common to the two + * constructors. + */ + void initialize (const UpdateFlags update_flags); }; diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index e20c30738f..9823b5b582 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -72,7 +72,7 @@ FEValuesData::initialize (const unsigned int n_quadrature_points, if (flags & update_normal_vectors) normal_vectors.resize(n_quadrature_points); -//Todo: support points missing +//Todo:[?] support points missing } @@ -422,10 +422,33 @@ FEValuesBase::memory_consumption () const }; + +template +UpdateFlags +FEValuesBase::compute_update_flags (const UpdateFlags update_flags) const +{ + + // first find out which objects + // need to be recomputed on each + // cell we visit. this we have to + // ask the finite element and mapping. + // elements are first since they + // might require update in mapping + UpdateFlags flags = update_flags + | fe->update_once (update_flags) + | fe->update_each (update_flags); + flags |= mapping->update_once (flags) + | mapping->update_each (flags); + + return flags; +}; + + + /*------------------------------- FEValues -------------------------------*/ -// TODO:[GK,WB] unify the two constructor bodies by calling a single function that does all the work + template FEValues::FEValues (const Mapping &mapping, const FiniteElement &fe, @@ -440,34 +463,7 @@ FEValues::FEValues (const Mapping &mapping, fe), quadrature (q) { - // you can't compute normal vectors - // on cells, only on faces - Assert ((update_flags & update_normal_vectors) == false, - typename FEValuesBase::ExcInvalidUpdateFlag()); - - // first find out which objects - // need to be recomputed on each - // cell we visit. this we have to - // ask the finite element and mapping. - // elements are first since they - // might require update in mapping - UpdateFlags flags = update_flags - | fe.update_once (update_flags) - | fe.update_each (update_flags); - flags |= mapping.update_once (flags) - | mapping.update_each (flags); - - // then get objects into which the - // FE and the Mapping can store - // intermediate data used across - // calls to reinit - mapping_data = mapping.get_data(flags, quadrature); - fe_data = fe.get_data(flags, mapping, quadrature); - - // set up objects within this class - FEValuesData::initialize(n_quadrature_points, - dofs_per_cell, - flags); + initialize (update_flags); }; @@ -484,36 +480,35 @@ FEValues::FEValues (const FiniteElement &fe, mapping_q1, fe), quadrature (q) +{ + initialize (update_flags); +}; + + + +template +void +FEValues::initialize (const UpdateFlags update_flags) { // you can't compute normal vectors // on cells, only on faces Assert ((update_flags & update_normal_vectors) == false, - FEValuesBase::ExcInvalidUpdateFlag()); - - // first find out which objects - // need to be recomputed on each - // cell we visit. this we have to - // ask the finite element and mapping. - // elements are first since they - // might require update in mapping - UpdateFlags flags = update_flags - | fe.update_once (update_flags) - | fe.update_each (update_flags); - flags |= mapping->update_once (flags) - | mapping->update_each (flags); + typename FEValuesBase::ExcInvalidUpdateFlag()); - // then get objects into which the + const UpdateFlags flags = compute_update_flags (update_flags); + + // then get objects into which the // FE and the Mapping can store // intermediate data used across // calls to reinit mapping_data = mapping->get_data(flags, quadrature); - fe_data = fe.get_data(flags, *mapping, quadrature); + fe_data = fe->get_data(flags, *mapping, quadrature); // set up objects within this class FEValuesData::initialize(n_quadrature_points, dofs_per_cell, flags); -} +}; @@ -602,7 +597,6 @@ FEFaceValuesBase::get_boundary_forms () const /*------------------------------- FEFaceValues -------------------------------*/ -// TODO:[GK,WB] unify the two constructor bodies by calling a single function that does all the work template FEFaceValues::FEFaceValues (const Mapping &mapping, const FiniteElement &fe, @@ -616,29 +610,7 @@ FEFaceValues::FEFaceValues (const Mapping &mapping, mapping, fe, quadrature) { - // first find out which objects - // need to be recomputed on each - // cell we visit. this we have to - // ask the finite element and mapping. - // elements are first since they - // might require update in mapping - UpdateFlags flags = update_flags - | fe.update_once (update_flags) - | fe.update_each (update_flags); - flags |= mapping.update_once (flags) - | mapping.update_each (flags); - - // then get objects into which the - // FE and the Mapping can store - // intermediate data used across - // calls to reinit - mapping_data = mapping.get_face_data(flags, quadrature); - fe_data = fe.get_face_data(flags, mapping, quadrature); - - // set up objects within this class - FEValuesData::initialize(n_quadrature_points, - dofs_per_cell, - flags); + initialize (update_flags); }; @@ -655,25 +627,24 @@ FEFaceValues::FEFaceValues (const FiniteElement &fe, mapping_q1, fe, quadrature) { - // first find out which objects - // need to be recomputed on each - // cell we visit. this we have to - // ask the finite element and mapping. - // elements are first since they - // might require update in mapping - UpdateFlags flags = update_flags - | fe.update_once (update_flags) - | fe.update_each (update_flags); - flags |= mapping->update_once (flags) - | mapping->update_each (flags); + initialize (update_flags); +}; + + +template +void +FEFaceValues::initialize (const UpdateFlags update_flags) +{ + const UpdateFlags flags = compute_update_flags (update_flags); + // then get objects into which the // FE and the Mapping can store // intermediate data used across // calls to reinit mapping_data = mapping->get_face_data(flags, quadrature); - fe_data = fe.get_face_data(flags, *mapping, quadrature); - + fe_data = fe->get_face_data(flags, *mapping, quadrature); + // set up objects within this class FEValuesData::initialize(n_quadrature_points, dofs_per_cell, @@ -717,8 +688,6 @@ void FEFaceValues::reinit (const typename DoFHandler::cell_iterator &c /*------------------------------- FESubFaceValues -------------------------------*/ -// TODO:[GK,WB] unify the two constructor bodies by calling a single function that does all the work -// TODO:[GK,WB] FESubfaceValues constructors are exactly the same as FEFaceValues constructors with exception of calling fe_sub_face_data instead of fe_face_data. Merge! template FESubfaceValues::FESubfaceValues (const Mapping &mapping, const FiniteElement &fe, @@ -733,29 +702,7 @@ FESubfaceValues::FESubfaceValues (const Mapping &mapping, mapping, fe, quadrature) { - // first find out which objects - // need to be recomputed on each - // cell we visit. this we have to - // ask the finite element and mapping. - // elements are first since they - // might require update in mapping - UpdateFlags flags = update_flags - | fe.update_once (update_flags) - | fe.update_each (update_flags); - flags |= mapping.update_once (flags) - | mapping.update_each (flags); - - // then get objects into which the - // FE and the Mapping can store - // intermediate data used across - // calls to reinit - mapping_data = mapping.get_subface_data(flags, quadrature); - fe_data = fe.get_subface_data(flags, mapping, quadrature); - - // set up objects within this class - FEValuesData::initialize(n_quadrature_points, - dofs_per_cell, - flags); + initialize (update_flags); } @@ -773,30 +720,29 @@ FESubfaceValues::FESubfaceValues (const FiniteElement &fe, mapping_q1, fe, quadrature) { - // first find out which objects - // need to be recomputed on each - // cell we visit. this we have to - // ask the finite element and mapping. - // elements are first since they - // might require update in mapping - UpdateFlags flags = update_flags - | fe.update_once (update_flags) - | fe.update_each (update_flags); - flags |= mapping->update_once (flags) - | mapping->update_each (flags); + initialize (update_flags); +}; + + +template +void +FESubfaceValues::initialize (const UpdateFlags update_flags) +{ + const UpdateFlags flags = compute_update_flags (update_flags); + // then get objects into which the // FE and the Mapping can store // intermediate data used across // calls to reinit mapping_data = mapping->get_subface_data(flags, quadrature); - fe_data = fe.get_subface_data(flags, *mapping, quadrature); - + fe_data = fe->get_subface_data(flags, *mapping, quadrature); + // set up objects within this class FEValuesData::initialize(n_quadrature_points, dofs_per_cell, flags); -} +}; -- 2.39.5