]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new local renumbering
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 25 May 2004 12:19:30 +0000 (12:19 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 25 May 2004 12:19:30 +0000 (12:19 +0000)
git-svn-id: https://svn.dealii.org/trunk@9317 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 84c9a09058e380c18cd12aeaafb2b0ef3f1cfc04..f08a7056840fffa6fd852d4ba91afc101becd7a7 100644 (file)
@@ -41,11 +41,38 @@ class ConstraintMatrix;
  * $id-I_h$ that is needed for evaluating $(id-I_h)z$ for e.g. the
  * dual solution $z$.
  *
- * @author Ralf Hartmann, 2000
+ * @author Ralf Hartmann, Guido Kanschat, 2000, 2004
  */
 class FETools
 {
   public:
+                                    /**
+                                     * Compute the vector required to
+                                     * renumber the dofs of a cell by
+                                     * component. Furthermore,
+                                     * compute the vector storing the
+                                     * start indices of each
+                                     * component in the local block
+                                     * vector.
+                                     *
+                                     * The second vector is organized
+                                     * such that there is a vector
+                                     * for each base element
+                                     * containing the start index for
+                                     * each component served by this
+                                     * base element.
+                                     *
+                                     * While the first vector is
+                                     * checked to have the correct
+                                     * size, the second one is
+                                     * reinitialized for convenience.
+                                     */
+    template<int dim>
+    static void compute_component_wise(
+      const FiniteElement<dim>&                fe,
+      std::vector<unsigned int>&               renumbering,
+      std::vector<std::vector<unsigned int> >& start_indices);
+    
                                     /**
                                      * Gives the interpolation matrix
                                      * that interpolates a @p fe1-
index ff26769d583f3e05244b744b4ee34a77a503f163..28eeb32df805b3cf88a67cfc49332efb632f045d 100644 (file)
@@ -196,6 +196,47 @@ namespace
   }
 }
 
+template<int dim>
+void FETools::compute_component_wise(
+  const FiniteElement<dim>& element,
+  std::vector<unsigned int>& renumbering,
+  std::vector<std::vector<unsigned int> >& comp_start)
+{
+  Assert(renumbering.size() == element.dofs_per_cell,
+        ExcDimensionMismatch(renumbering.size(),
+                             element.dofs_per_cell));
+  
+  comp_start.resize(element.n_base_elements());
+  
+  unsigned int k=0;
+  for (unsigned int i=0;i<comp_start.size();++i)
+    {
+      comp_start[i].resize(element.element_multiplicity(i));
+      const unsigned int increment
+       = element.base_element(i).dofs_per_cell;
+      
+      for (unsigned int j=0;j<comp_start[i].size();++j)
+       {
+         comp_start[i][j] = k;
+         k += increment;
+       }
+    }
+  
+                                  // For each index i of the
+                                  // unstructured cellwise
+                                  // numbering, renumbering
+                                  // contains the index of the
+                                  // cell-block numbering
+  for (unsigned int i=0;i<element.dofs_per_cell;++i)
+    {
+      std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
+       indices = element.system_to_base_index(i);
+      renumbering[i] = comp_start[indices.first.first][indices.first.second]
+                            +indices.second;
+    }
+}
+
+
 
 template <int dim, typename number>
 void FETools::get_interpolation_matrix (const FiniteElement<dim> &fe1,
@@ -1539,6 +1580,11 @@ FETools::get_fe_from_name_aux (const std::string &name)
 
 /*-------------- Explicit Instantiations -------------------------------*/
 
+
+template
+void FETools::compute_component_wise(
+  const FiniteElement<deal_II_dimension>& element,
+  std::vector<unsigned int>&, std::vector<std::vector<unsigned int> >&);
 template
 void FETools::get_interpolation_matrix<deal_II_dimension>
 (const FiniteElement<deal_II_dimension> &,

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.