]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add constructor for FESystem with four base elements
authorThomas Wick <twick@ices.utexas.edu>
Fri, 8 Jan 2010 12:41:20 +0000 (12:41 +0000)
committerThomas Wick <twick@ices.utexas.edu>
Fri, 8 Jan 2010 12:41:20 +0000 (12:41 +0000)
git-svn-id: https://svn.dealii.org/trunk@20333 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/source/fe/fe_system.cc
deal.II/doc/news/changes.h

index 8339b49fefb9d48d35db1c292633f9fff8cca29f..4e12dc02ba74a0a31abdaa89fac6cc0cfff4f4d3 100644 (file)
@@ -129,6 +129,18 @@ class FESystem : public FiniteElement<dim,spacedim>
              const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
              const FiniteElement<dim,spacedim> &fe3, const unsigned int n3);
 
+                                    /** 
+                                     * Constructor for mixed
+                                     * discretizations with four
+                                     * base elements.
+                                     *
+                                     * See the other constructor.
+                                     */
+    FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
+             const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
+             const FiniteElement<dim,spacedim> &fe3, const unsigned int n3,
+             const FiniteElement<dim,spacedim> &fe4, const unsigned int n4);
+
                                     /**
                                      * Destructor.
                                      */
@@ -783,6 +795,19 @@ class FESystem : public FiniteElement<dim,spacedim>
                          const FiniteElementData<dim> &fe3,
                          const unsigned int            N3);
 
+                                  /**
+                                   * with 4 different sub-elements
+                                   */
+  static FiniteElementData<dim>
+  multiply_dof_numbers (const FiniteElementData<dim> &fe1,
+                       const unsigned int            N1,
+                       const FiniteElementData<dim> &fe2,
+                       const unsigned int            N2,
+                       const FiniteElementData<dim> &fe3,
+                       const unsigned int            N3,
+                       const FiniteElementData<dim> &fe4,
+                       const unsigned int            N4);
+
 
                                     /**
                                      * Helper function used in the
@@ -821,6 +846,19 @@ class FESystem : public FiniteElement<dim,spacedim>
                                           const FiniteElement<dim,spacedim> &fe3,
                                           const unsigned int        N3);
 
+                                  /**
+                                   *  with four different sub-elements
+                                   */
+    static std::vector<bool>
+    compute_restriction_is_additive_flags (const FiniteElement<dim,spacedim> &fe1,
+                                          const unsigned int        N1,
+                                          const FiniteElement<dim,spacedim> &fe2,
+                                          const unsigned int        N2,
+                                          const FiniteElement<dim,spacedim> &fe3,
+                                          const unsigned int        N3,
+                                          const FiniteElement<dim,spacedim> &fe4,
+                                          const unsigned int        N4);
+
                                     /**
                                      * Compute the named flags for a
                                      * list of finite elements with
@@ -828,7 +866,7 @@ class FESystem : public FiniteElement<dim,spacedim>
                                      * second argument. This function
                                      * is called from all the above
                                      * functions.
-                                    */
+                                     */
     static std::vector<bool>
     compute_restriction_is_additive_flags (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
                                            const std::vector<unsigned int>              &multiplicities);
@@ -867,6 +905,21 @@ class FESystem : public FiniteElement<dim,spacedim>
                                const FiniteElement<dim,spacedim> &fe3,
                                const unsigned int        N3);
 
+                                  /**
+                                   * Compute the non-zero vector
+                                   * components of a composed
+                                   * finite element.
+                                   */
+   static std::vector<std::vector<bool> >
+    compute_nonzero_components (const FiniteElement<dim,spacedim> &fe1,
+                               const unsigned int        N1,
+                               const FiniteElement<dim,spacedim> &fe2,
+                               const unsigned int        N2,
+                               const FiniteElement<dim,spacedim> &fe3,
+                               const unsigned int        N3,
+                               const FiniteElement<dim,spacedim> &fe4,
+                               const unsigned int        N4);
+
                                     /**
                                      * Compute the nonzero components
                                      * of a list of finite elements
index 71bc7606928d5994b643292c4a84caa38cc66525..98063d67219f6d7e89c21d863d855dec6c829156 100644 (file)
@@ -211,6 +211,44 @@ FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
 }
 
 
+template <int dim, int spacedim>
+FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
+                                 const unsigned int        n1,
+                                 const FiniteElement<dim,spacedim> &fe2,
+                                 const unsigned int        n2,
+                                 const FiniteElement<dim,spacedim> &fe3,
+                                 const unsigned int        n3,
+                                 const FiniteElement<dim,spacedim> &fe4,
+                                 const unsigned int        n4) :
+  FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1,
+                                                   fe2, n2,
+                                                   fe3, n3,
+                                                   fe4, n4),
+                              compute_restriction_is_additive_flags (fe1, n1,
+                                                                     fe2, n2,
+                                                                     fe3, n3,
+                                                                     fe4, n4),
+                              compute_nonzero_components(fe1, n1,
+                                                         fe2, n2,
+                                                         fe3, n3,
+                                                         fe4 ,n4)),
+  base_elements(4)
+{
+  base_elements[0] = ElementPair(fe1.clone(), n1);  
+  base_elements[0].first->subscribe (typeid(*this).name());
+  base_elements[1] = ElementPair(fe2.clone(), n2);
+  base_elements[1].first->subscribe (typeid(*this).name());
+  base_elements[2] = ElementPair(fe3.clone(), n3);
+  base_elements[2].first->subscribe (typeid(*this).name());
+  base_elements[3] = ElementPair(fe4.clone(), n4);
+  base_elements[3].first->subscribe (typeid(*this).name());
+  this->first_block_of_base_table.push_back(n1);
+  this->first_block_of_base_table.push_back(n1+n2);
+  this->first_block_of_base_table.push_back(n1+n2+n3);
+  initialize ();
+}
+
+
 
 template <int dim, int spacedim>
 FESystem<dim,spacedim>::~FESystem ()
@@ -262,23 +300,32 @@ FESystem<dim,spacedim>::clone() const
 {
   switch (n_base_elements())
     {
-      case 1:
-            return new FESystem(base_element(0),
-                                element_multiplicity(0));
-      case 2:
-            return new FESystem(base_element(0),
-                                element_multiplicity(0),
-                                base_element(1),
-                                element_multiplicity(1));
-      case 3:
-            return new FESystem(base_element(0),
-                                element_multiplicity(0),
-                                base_element(1),
-                                element_multiplicity(1),
-                                base_element(2),
-                                element_multiplicity(2));
+    case 1:
+      return new FESystem(base_element(0),
+                         element_multiplicity(0));
+    case 2:
+      return new FESystem(base_element(0),
+                         element_multiplicity(0),
+                         base_element(1),
+                         element_multiplicity(1));
+    case 3:
+      return new FESystem(base_element(0),
+                         element_multiplicity(0),
+                         base_element(1),
+                         element_multiplicity(1),
+                         base_element(2),
+                         element_multiplicity(2));
+    case 4:
+      return new FESystem(base_element(0),
+                         element_multiplicity(0),
+                         base_element(1),
+                         element_multiplicity(1),
+                         base_element(2),
+                         element_multiplicity(2),
+                         base_element(3),
+                         element_multiplicity(3));
       default:
-            Assert(false, ExcNotImplemented());
+                         Assert(false, ExcNotImplemented());
     }
   return 0;
 }
@@ -2522,6 +2569,52 @@ FESystem<dim,spacedim>::multiply_dof_numbers (const FiniteElementData<dim> &fe1,
 
 
 
+template <int dim, int spacedim>
+FiniteElementData<dim>
+FESystem<dim,spacedim>::multiply_dof_numbers (const FiniteElementData<dim> &fe1,
+                                             const unsigned int            N1,
+                                             const FiniteElementData<dim> &fe2,
+                                             const unsigned int            N2,
+                                             const FiniteElementData<dim> &fe3,
+                                             const unsigned int            N3,
+                                             const FiniteElementData<dim> &fe4,
+                                             const unsigned int            N4)
+{
+  std::vector<unsigned int> dpo;
+  dpo.push_back(fe1.dofs_per_vertex * N1 +
+               fe2.dofs_per_vertex * N2 +
+               fe3.dofs_per_vertex * N3 +
+               fe4.dofs_per_vertex * N4);
+  dpo.push_back(fe1.dofs_per_line * N1 +
+               fe2.dofs_per_line * N2 +
+               fe3.dofs_per_line * N3 +
+               fe4.dofs_per_line * N4);
+  if (dim>1) dpo.push_back(fe1.dofs_per_quad * N1 +
+                          fe2.dofs_per_quad * N2 +
+                          fe3.dofs_per_quad * N3 +
+                          fe4.dofs_per_quad * N4);
+  if (dim>2) dpo.push_back(fe1.dofs_per_hex * N1 +
+                          fe2.dofs_per_hex * N2 +
+                          fe3.dofs_per_hex * N3 +
+                          fe4.dofs_per_hex * N4);
+                                  // degree is the maximal degree of the components
+  const unsigned int
+    degree = std::max(std::max (std::max(fe1.tensor_degree(),fe2.tensor_degree()),fe3.tensor_degree()), fe4.tensor_degree());
+                      
+  return FiniteElementData<dim> (
+    dpo,
+    fe1.n_components() * N1 + fe2.n_components() * N2 + fe3.n_components() * N3 +
+    fe4.n_components() * N4,
+    degree,
+    typename FiniteElementData<dim>::Conformity(fe1.conforming_space
+                                               & fe2.conforming_space
+                                               & fe3.conforming_space
+                                               & fe4.conforming_space),
+    N1 + N2 + N3 + N4);
+}
+
+
+
 template <int dim, int spacedim>
 std::vector<bool>
 FESystem<dim,spacedim>::compute_restriction_is_additive_flags (const FiniteElement<dim,spacedim> &fe,
@@ -2584,6 +2677,36 @@ FESystem<dim,spacedim>::compute_restriction_is_additive_flags (const FiniteEleme
 }
 
 
+template <int dim, int spacedim>
+std::vector<bool>
+FESystem<dim,spacedim>::compute_restriction_is_additive_flags (const FiniteElement<dim,spacedim> &fe1,
+                                                              const unsigned int        N1,
+                                                              const FiniteElement<dim,spacedim> &fe2,
+                                                              const unsigned int        N2,
+                                                              const FiniteElement<dim,spacedim> &fe3,
+                                                              const unsigned int        N3,
+                                                              const FiniteElement<dim,spacedim> &fe4,
+                                                              const unsigned int        N4) 
+{
+  std::vector<const FiniteElement<dim,spacedim>*> fe_list;
+  std::vector<unsigned int>              multiplicities;
+
+  fe_list.push_back (&fe1);
+  multiplicities.push_back (N1);
+
+  fe_list.push_back (&fe2);
+  multiplicities.push_back (N2);
+
+  fe_list.push_back (&fe3);
+  multiplicities.push_back (N3);
+
+  fe_list.push_back (&fe4);
+  multiplicities.push_back (N4);
+  
+  return compute_restriction_is_additive_flags (fe_list, multiplicities);
+}
+
+
 
 template <int dim, int spacedim>
 std::vector<bool>
@@ -2780,6 +2903,37 @@ FESystem<dim,spacedim>::compute_nonzero_components (const FiniteElement<dim,spac
 
 
 
+template <int dim, int spacedim>
+std::vector<std::vector<bool> >
+FESystem<dim,spacedim>::compute_nonzero_components (const FiniteElement<dim,spacedim> &fe1,
+                                                   const unsigned int        N1,
+                                                   const FiniteElement<dim,spacedim> &fe2,
+                                                   const unsigned int        N2,
+                                                   const FiniteElement<dim,spacedim> &fe3,
+                                                   const unsigned int        N3,
+                                                   const FiniteElement<dim,spacedim> &fe4,
+                                                   const unsigned int        N4)
+{
+  std::vector<const FiniteElement<dim,spacedim>*> fe_list;
+  std::vector<unsigned int>              multiplicities;
+
+  fe_list.push_back (&fe1);
+  multiplicities.push_back (N1);
+
+  fe_list.push_back (&fe2);
+  multiplicities.push_back (N2);
+
+  fe_list.push_back (&fe3);
+  multiplicities.push_back (N3);
+
+  fe_list.push_back (&fe4);
+  multiplicities.push_back (N4);
+  
+  return compute_nonzero_components (fe_list, multiplicities);
+}
+
+
+
 template <int dim, int spacedim>
 std::vector<std::vector<bool> >
 FESystem<dim,spacedim>::
index 446b0ef498754ece9b4d658c6f2bb2757ae75ec3..d29ada68aa405c72e1499ed94eb3dea89dd9f2a0 100644 (file)
@@ -569,6 +569,14 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li>
+  <p>
+  New: Constructor of FESystem now also exists for four base elements.
+  <br>
+  (Thomas Wick 2010/01/08)
+  </p>
+  </li>
+
   <li>
   <p>
   New: The functions in namespace DoFRenumbering::boost are now also compiled for

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.