]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move several functions from the .h to the .cc file as they are no
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 7 Mar 2001 17:26:07 +0000 (17:26 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 7 Mar 2001 17:26:07 +0000 (17:26 +0000)
longer templated on other types. Minor clean-ups. Fix a bug where the
wrong array was resized. Fix a bug where clone would not have worked
because of missing 'break's in a switch statement.

git-svn-id: https://svn.dealii.org/trunk@4155 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 295ebc49622151ae04c82336cea7c778d08dd8eb..54bae219a13bc475d10f139d915929f1f9267596 100644 (file)
@@ -14,7 +14,7 @@
 #define __deal2__fe_system_h
 
 
-/*----------------------------   fe_lib.system.h     ---------------------------*/
+/*----------------------------   fe_system.h     ---------------------------*/
 
 
 #include <fe/fe.h>
@@ -194,12 +194,13 @@ class FESystem : public FiniteElement<dim>
                                      * 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 after
-                                     * calling this function.
+                                     * 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_support_points (typename std::vector<Point<dim> > &) const;    
 
@@ -572,151 +573,41 @@ class FESystem : public FiniteElement<dim>
 };
 
 
-/* ------------------------- template functions ------------------------- */
+/* ------------------------- inline functions ------------------------- */
 
 template<int dim>
 inline unsigned int
 FESystem<dim>::n_base_elements() const
 {
   return base_elements.size();
-}
-
-
-template <int dim>
-FESystem<dim>::FESystem (const FiniteElement<dim> &fe, const unsigned int n_elements) :
-               FiniteElement<dim> (multiply_dof_numbers(fe, n_elements),
-                                   compute_restriction_is_additive_flags (fe, n_elements)),
-  base_elements(1)
-{
-  base_elements[0] = ElementPair(fe.clone(), n_elements);
-  base_elements[0].first -> subscribe ();
-  initialize ();
 };
 
 
-template <int dim>
-FESystem<dim>::FESystem (const FiniteElement<dim> &fe1, const unsigned int n1,
-                        const FiniteElement<dim> &fe2, const unsigned int n2)
-               :
-               FiniteElement<dim> (multiply_dof_numbers(fe1, n1, fe2, n2),
-                                   compute_restriction_is_additive_flags (fe1, n1,
-                                                                          fe2, n2)),
-  base_elements(2)
-{
-  base_elements[0] = ElementPair(fe1.clone(), n1);
-  base_elements[0].first -> subscribe ();
-  base_elements[1] = ElementPair(fe2.clone(), n2);
-  base_elements[1].first -> subscribe ();
-  initialize ();
-};
-
-
-template <int dim>
-FESystem<dim>::FESystem (const FiniteElement<dim> &fe1, const unsigned int n1,
-                        const FiniteElement<dim> &fe2, const unsigned int n2,
-                        const FiniteElement<dim> &fe3, const unsigned int n3)
-               :
-               FiniteElement<dim> (multiply_dof_numbers(fe1, n1,
-                                                        fe2, n2,
-                                                        fe3, n3),
-                                   compute_restriction_is_additive_flags (fe1, n1,
-                                                                          fe2, n2,
-                                                                          fe3, n3)),
-  base_elements(3)
-{
-  base_elements[0] = ElementPair(fe1.clone(), n1);  
-  base_elements[0].first -> subscribe ();
-  base_elements[1] = ElementPair(fe2.clone(), n2);
-  base_elements[1].first -> subscribe ();
-  base_elements[2] = ElementPair(fe3.clone(), n3);
-  base_elements[2].first -> subscribe ();
-  initialize ();
-};
-
 
 template<int dim>
 inline unsigned int
-FESystem<dim>::element_multiplicity(unsigned int index) const
+FESystem<dim>::element_multiplicity (const unsigned int index) const
 {
   Assert (index < base_elements.size(), 
          ExcIndexRange(index, 0, base_elements.size()));
   return base_elements[index].second;
-}
+};
+
 
 
 template <int dim>
-inline const FiniteElement<dim>&
-FESystem<dim>::base_element(unsigned int index) const
+inline const FiniteElement<dim> &
+FESystem<dim>::base_element (const unsigned int index) const
 {
   Assert (index < base_elements.size(), 
          ExcIndexRange(index, 0, base_elements.size()));
   return *base_elements[index].first;
-}
-
-
-template <int dim>
-inline FiniteElementBase<dim>::InternalDataBase &
-FESystem<dim>::
-InternalData::get_fe_data(unsigned int base_no) const
-{
-  Assert(base_no<base_fe_datas.size(),
-        ExcIndexRange(base_no,0,base_fe_datas.size()));
-  return *base_fe_datas[base_no];
-}
-
-
-
-template <int dim>
-inline void
-FESystem<dim>::
-InternalData::set_fe_data(unsigned int base_no,
-                         FiniteElementBase<dim>::InternalDataBase *ptr)
-{
-  Assert(base_no<base_fe_datas.size(),
-        ExcIndexRange(base_no,0,base_fe_datas.size()));
-  base_fe_datas[base_no]=ptr;
-}
-
-
-template <int dim>
-inline FEValuesData<dim> &
-FESystem<dim>::
-InternalData::get_fe_values_data(unsigned int base_no) const
-{
-  Assert(base_no<base_fe_values_datas.size(),
-        ExcIndexRange(base_no,0,base_fe_values_datas.size()));
-  Assert(base_fe_values_datas[base_no]!=0, ExcInternalError());
-  return *base_fe_values_datas[base_no];
-}
-
-
-
-template <int dim>
-inline void
-FESystem<dim>::
-InternalData::set_fe_values_data(unsigned int base_no,
-                                FEValuesData<dim> *ptr)
-{
-  Assert(base_no<base_fe_values_datas.size(),
-        ExcIndexRange(base_no,0,base_fe_values_datas.size()));
-  base_fe_values_datas[base_no]=ptr;
-}
+};
 
 
-template <int dim>
-inline void
-FESystem<dim>::
-InternalData::delete_fe_values_data(unsigned int base_no)
-{
-  Assert(base_no<base_fe_values_datas.size(),
-        ExcIndexRange(base_no,0,base_fe_values_datas.size()));
-  Assert(base_fe_values_datas[base_no]!=0, ExcInternalError());
-  delete base_fe_values_datas[base_no];
-  base_fe_values_datas[base_no]=0;
-}
 
 
-/*----------------------------  fe_lib.system.h  ---------------------------*/
 
+/*----------------------------  fe_system.h  ---------------------------*/
 #endif
-/*----------------------------  fe_lib.system.h  ---------------------------*/
+/*----------------------------  fe_system.h  ---------------------------*/
index 18e0c872960e9a8e5242b6bba3ed3b3888c0d72a..8dec1d71685d87575193d1a7e4923ed6390fdb18 100644 (file)
@@ -28,70 +28,212 @@ using namespace std;
 #endif
 
 
+/* ----------------------- FESystem::InternalData ------------------- */
+
+
+template <int dim>
+inline FiniteElementBase<dim>::InternalDataBase &
+FESystem<dim>::
+InternalData::get_fe_data (const unsigned int base_no) const
+{
+  Assert(base_no<base_fe_datas.size(),
+        ExcIndexRange(base_no,0,base_fe_datas.size()));
+  return *base_fe_datas[base_no];
+};
+
+
+
+template <int dim>
+inline void
+FESystem<dim>::
+InternalData::set_fe_data (const unsigned int base_no,
+                          FiniteElementBase<dim>::InternalDataBase *ptr)
+{
+  Assert(base_no<base_fe_datas.size(),
+        ExcIndexRange(base_no,0,base_fe_datas.size()));
+  base_fe_datas[base_no]=ptr;
+};
+
+
+
+template <int dim>
+inline FEValuesData<dim> &
+FESystem<dim>::
+InternalData::get_fe_values_data (const unsigned int base_no) const
+{
+  Assert(base_no<base_fe_values_datas.size(),
+        ExcIndexRange(base_no,0,base_fe_values_datas.size()));
+  Assert(base_fe_values_datas[base_no]!=0, ExcInternalError());
+  return *base_fe_values_datas[base_no];
+};
+
+
+
+template <int dim>
+inline void
+FESystem<dim>::
+InternalData::set_fe_values_data (const unsigned int base_no,
+                                 FEValuesData<dim> *ptr)
+{
+  Assert(base_no<base_fe_values_datas.size(),
+        ExcIndexRange(base_no,0,base_fe_values_datas.size()));
+  base_fe_values_datas[base_no]=ptr;
+};
+
+
+
+template <int dim>
+inline void
+FESystem<dim>::
+InternalData::delete_fe_values_data (const unsigned int base_no)
+{
+  Assert(base_no<base_fe_values_datas.size(),
+        ExcIndexRange(base_no,0,base_fe_values_datas.size()));
+  Assert(base_fe_values_datas[base_no]!=0, ExcInternalError());
+  delete base_fe_values_datas[base_no];
+  base_fe_values_datas[base_no]=0;
+}
+
+
+
+/* ---------------------------------- FESystem ------------------- */
+
+
+template <int dim>
+FESystem<dim>::FESystem (const FiniteElement<dim> &fe, const unsigned int n_elements) :
+               FiniteElement<dim> (multiply_dof_numbers(fe, n_elements),
+                                   compute_restriction_is_additive_flags (fe, n_elements)),
+                base_elements(1)
+{
+  base_elements[0] = ElementPair(fe.clone(), n_elements);
+  base_elements[0].first->subscribe ();
+  initialize ();
+};
+
+
+
+template <int dim>
+FESystem<dim>::FESystem (const FiniteElement<dim> &fe1, const unsigned int n1,
+                        const FiniteElement<dim> &fe2, const unsigned int n2) :
+               FiniteElement<dim> (multiply_dof_numbers(fe1, n1, fe2, n2),
+                                   compute_restriction_is_additive_flags (fe1, n1,
+                                                                          fe2, n2)),
+                base_elements(2)
+{
+  base_elements[0] = ElementPair(fe1.clone(), n1);
+  base_elements[0].first->subscribe ();
+  base_elements[1] = ElementPair(fe2.clone(), n2);
+  base_elements[1].first->subscribe ();
+  initialize ();
+};
+
+
+
+template <int dim>
+FESystem<dim>::FESystem (const FiniteElement<dim> &fe1, const unsigned int n1,
+                        const FiniteElement<dim> &fe2, const unsigned int n2,
+                        const FiniteElement<dim> &fe3, const unsigned int n3) :
+               FiniteElement<dim> (multiply_dof_numbers(fe1, n1,
+                                                        fe2, n2,
+                                                        fe3, n3),
+                                   compute_restriction_is_additive_flags (fe1, n1,
+                                                                          fe2, n2,
+                                                                          fe3, n3)),
+                base_elements(3)
+{
+  base_elements[0] = ElementPair(fe1.clone(), n1);  
+  base_elements[0].first->subscribe ();
+  base_elements[1] = ElementPair(fe2.clone(), n2);
+  base_elements[1].first->subscribe ();
+  base_elements[2] = ElementPair(fe3.clone(), n3);
+  base_elements[2].first->subscribe ();
+  initialize ();
+};
+
+
 
 template <int dim>
 FESystem<dim>::~FESystem ()
 {
-  for (unsigned i=0;i<base_elements.size();++i)
+                                  // delete base elements created in
+                                  // the constructor
+  for (unsigned i=0; i<base_elements.size(); ++i)
     {
-      base_elements[i].first -> unsubscribe ();
+      base_elements[i].first->unsubscribe();
       delete base_elements[i].first;
+      base_elements[i].first = 0;
     }
 };
 
 
+
 template <int dim>
 FiniteElement<dim>*
 FESystem<dim>::clone() const
 {
-  FiniteElement<dim> *fe=0;
   switch (n_base_elements())
     {
       case 1:
-           fe=new FESystem(base_element(0),
-                           element_multiplicity(0));
+           return new FESystem(base_element(0),
+                               element_multiplicity(0));
       case 2:
-           fe=new FESystem(base_element(0),
-                           element_multiplicity(0),
-                           base_element(1),
-                           element_multiplicity(1));
+           return new FESystem(base_element(0),
+                               element_multiplicity(0),
+                               base_element(1),
+                               element_multiplicity(1));
       case 3:
-           fe=new FESystem(base_element(0),
-                           element_multiplicity(0),
-                           base_element(1),
-                           element_multiplicity(1),
-                           base_element(2),
-                           element_multiplicity(2));
+           return new FESystem(base_element(0),
+                               element_multiplicity(0),
+                               base_element(1),
+                               element_multiplicity(1),
+                               base_element(2),
+                               element_multiplicity(2));
       default:
            Assert(false, ExcNotImplemented());
     }
-  return fe;
+  return 0;
 }
 
 
+
 template <int dim>
-void FESystem<dim>::get_unit_support_points (
-  typename std::vector<Point<dim> > &unit_support_points) const
+void
+FESystem<dim>::
+get_unit_support_points (typename std::vector<Point<dim> > &unit_support_points) const
 {
+                                  // generate unit support points
+                                  // from unit support points of sub
+                                  // elements
   unit_support_points.resize(dofs_per_cell);
   
   typename std::vector<Point<dim> > base_unit_support_points (base_element(0).dofs_per_cell);
   unsigned int comp = 0;
-  for (unsigned int base_el=0 ; base_el<n_base_elements(); ++base_el)
+  for (unsigned int base_el=0; base_el<n_base_elements(); ++base_el)
     {
-      const unsigned int base_element_dofs_per_cell
-       =base_element(base_el).dofs_per_cell;
+      const unsigned int
+       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
+                                      // support points for the
+                                      // composed element, so return
+                                      // an empty array to
+                                      // demonstrate that fact
       if (base_unit_support_points.size()==0)
        {
-         base_unit_support_points.resize(0);
+         unit_support_points.resize(0);
          return;
        }
-      
+
+                                      // otherwise distribute the
+                                      // support of this base element
+                                      // to all degrees of freedom
+                                      // contributed by it
       Assert(base_unit_support_points.size()==base_element_dofs_per_cell,
             ExcInternalError());
-      for (unsigned int n = 0 ; n < element_multiplicity(base_el); ++n, ++comp)
+      for (unsigned int n=0; n<element_multiplicity(base_el); ++n, ++comp)
        for (unsigned int i=0; i<base_element_dofs_per_cell; ++i)
          unit_support_points[component_to_system_index(comp,i)]
            = base_unit_support_points[i];

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.