]> https://gitweb.dealii.org/ - dealii.git/commitdiff
renumbering by block
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 10 Mar 2006 14:52:37 +0000 (14:52 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 10 Mar 2006 14:52:37 +0000 (14:52 +0000)
git-svn-id: https://svn.dealii.org/trunk@12576 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 39967576bcb674af0b90aa3ea3c719ffec45c473..a97f56d041783b4f34c4cca0bab4ddf586bceb6f 100644 (file)
@@ -58,6 +58,10 @@ class FETools
 {
   public:
                                     /**
+                                     * @warning In most cases, you
+                                     * will probably want to use
+                                     * compute_base_renumbering().
+                                     *
                                      * Compute the vector required to
                                      * renumber the dofs of a cell by
                                      * component. Furthermore,
@@ -84,6 +88,23 @@ class FETools
       std::vector<unsigned int>&               renumbering,
       std::vector<std::vector<unsigned int> >& start_indices);
 
+                                    /**
+                                     * Compute the vector required to
+                                     * renumber the dofs of a cell by
+                                     * block. Furthermore,
+                                     * compute the vector storing the
+                                     * start indices of each local block
+                                     * vector.
+                                     *
+                                     * @todo Which way does this
+                                     * vector map the numbers?
+                                     */
+    template<int dim>
+    static void compute_block_renumbering (
+      const FiniteElement<dim>&  fe,
+      std::vector<unsigned int>& renumbering,
+      std::vector<unsigned int>& start_indices);
+    
                                     /**
                                      * @name Generation of local matrices
                                      * @{
index b1d88164f00b30ba7b8b955bd631c25f41448e02..ab5a8f8086d44d126e90779ce750cbb171ff3d0f 100644 (file)
@@ -198,6 +198,40 @@ void FETools::compute_component_wise(
 
 
 
+template<int dim>
+void FETools::compute_block_renumbering (
+  const FiniteElement<dim>& element,
+  std::vector<unsigned int>& renumbering,
+  std::vector<unsigned int>& start_indices)
+{
+  Assert(renumbering.size() == element.dofs_per_cell,
+        ExcDimensionMismatch(renumbering.size(),
+                             element.dofs_per_cell));
+  Assert(start_indices.size() == element.n_blocks(),
+        ExcDimensionMismatch(start_indices.size(),
+                             element.n_blocks()));
+  
+  unsigned int k=0;
+  unsigned int i=0;
+  for (unsigned int b=0;b<element.n_base_elements();++b)
+    for (unsigned int m=0;m<element.element_multiplicity(b);++m)
+      {
+       start_indices[i++] = k;
+       k += element.base_element(k).n_dofs_per_cell();
+      }
+  Assert (i == element.n_blocks(), ExcInternalError());
+  
+  for (unsigned int i=0;i<element.dofs_per_cell;++i)
+    {
+      std::pair<unsigned int, unsigned int>
+       indices = element.system_to_block_index(i);
+      renumbering[i] = start_indices[indices.first]
+                      +indices.second;
+    }
+}
+
+
+
 template <int dim, typename number>
 void FETools::get_interpolation_matrix (const FiniteElement<dim> &fe1,
                                         const FiniteElement<dim> &fe2,

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.