]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Unify lots of code from the 6 constructors of the FE*Values classes.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 30 Mar 2001 09:20:19 +0000 (09:20 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 30 Mar 2001 09:20:19 +0000 (09:20 +0000)
git-svn-id: https://svn.dealii.org/trunk@4348 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/source/fe/fe_values.cc

index 6e4edfc9f72be16019844536a9a7399512b72366..dcc57eb7a444ced530037d6aafad41a0225aad00 100644 (file)
@@ -18,6 +18,7 @@
 #include <dofs/dof_handler.h>
 
 template <int dim> class FEValuesData;
+template <int dim> class FEValuesBase;
 template <int dim> class FEValues;
 template <int dim> class FEFaceValues;
 template <int dim> class FESubfaceValues;
@@ -256,6 +257,7 @@ class FiniteElement : public FiniteElementBase<dim>
                                      * Declare some other classes as
                                      * friends of this class.
                                      */
+    friend class FEValuesBase<dim>;
     friend class FEValues<dim>;
     friend class FEFaceValues<dim>;
     friend class FESubfaceValues<dim>;
index 554c26ec77fea5018fc474c9760332999baf3e41..cbbfdc05d8d8ae8fe0bc49d84994fca8494bed1d 100644 (file)
@@ -615,7 +615,25 @@ class FEValuesBase : protected FEValuesData<dim>
                                      * Internal data of finite element.
                                      */
     SmartPointer<typename Mapping<dim>::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<dim>
                                     /**
                                      * Store the quadrature formula here.
                                      */
-    const Quadrature<dim> quadrature;    
+    const Quadrature<dim> quadrature;
+
+                                    /**
+                                     * Do work common to the two
+                                     * constructors.
+                                     */
+    void initialize (const UpdateFlags update_flags);
 };
 
 
@@ -825,6 +849,14 @@ class FEFaceValues : public FEFaceValuesBase<dim>
                                      */
     void reinit (const typename DoFHandler<dim>::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<dim>
                                      * Exception
                                      */
     DeclException0 (ExcFaceHasNoSubfaces);
+
+  private:
+
+                                    /**
+                                     * Do work common to the two
+                                     * constructors.
+                                     */
+    void initialize (const UpdateFlags update_flags);
 };
 
 
index e20c30738f51b35887c712e674378df5c21d3e52..9823b5b582dc358f1bcd8a9c01a1ec1aa6f55b96 100644 (file)
@@ -72,7 +72,7 @@ FEValuesData<dim>::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<dim>::memory_consumption () const
 };
 
 
+
+template <int dim>
+UpdateFlags
+FEValuesBase<dim>::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 <int dim>
 FEValues<dim>::FEValues (const Mapping<dim>       &mapping,
                         const FiniteElement<dim> &fe,
@@ -440,34 +463,7 @@ FEValues<dim>::FEValues (const Mapping<dim>       &mapping,
                                   fe),
                quadrature (q)
 {
-                                  // you can't compute normal vectors
-                                  // on cells, only on faces
-  Assert ((update_flags & update_normal_vectors) == false,
-         typename FEValuesBase<dim>::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<dim>::initialize(n_quadrature_points,
-                               dofs_per_cell,
-                               flags);
+  initialize (update_flags);
 };
 
 
@@ -484,36 +480,35 @@ FEValues<dim>::FEValues (const FiniteElement<dim> &fe,
                                   mapping_q1,
                                   fe),
                 quadrature (q)
+{
+  initialize (update_flags);
+};
+
+
+
+template <int dim>
+void
+FEValues<dim>::initialize (const UpdateFlags update_flags)
 {
                                   // you can't compute normal vectors
                                   // on cells, only on faces
   Assert ((update_flags & update_normal_vectors) == false,
-         FEValuesBase<dim>::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<dim>::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<dim>::initialize(n_quadrature_points,
                                dofs_per_cell,
                                flags);
-}
+};
 
 
 
@@ -602,7 +597,6 @@ FEFaceValuesBase<dim>::get_boundary_forms () const
 /*------------------------------- FEFaceValues -------------------------------*/
 
 
-// TODO:[GK,WB] unify the two constructor bodies by calling a single function that does all the work
 template <int dim>
 FEFaceValues<dim>::FEFaceValues (const Mapping<dim>       &mapping,
                                 const FiniteElement<dim> &fe,
@@ -616,29 +610,7 @@ FEFaceValues<dim>::FEFaceValues (const Mapping<dim>       &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<dim>::initialize(n_quadrature_points,
-                               dofs_per_cell,
-                               flags);
+  initialize (update_flags);
 };
 
 
@@ -655,25 +627,24 @@ FEFaceValues<dim>::FEFaceValues (const FiniteElement<dim> &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 <int dim>
+void
+FEFaceValues<dim>::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<dim>::initialize(n_quadrature_points,
                                dofs_per_cell,
@@ -717,8 +688,6 @@ void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::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 <int dim>
 FESubfaceValues<dim>::FESubfaceValues (const Mapping<dim>       &mapping,
                                       const FiniteElement<dim> &fe,
@@ -733,29 +702,7 @@ FESubfaceValues<dim>::FESubfaceValues (const Mapping<dim>       &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<dim>::initialize(n_quadrature_points,
-                               dofs_per_cell,
-                               flags);
+  initialize (update_flags);
 }
 
 
@@ -773,30 +720,29 @@ FESubfaceValues<dim>::FESubfaceValues (const FiniteElement<dim> &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 <int dim>
+void
+FESubfaceValues<dim>::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<dim>::initialize(n_quadrature_points,
                                dofs_per_cell,
                                flags);
-}
+};
 
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.