]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement renumbering by component.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 22 Apr 1999 15:10:57 +0000 (15:10 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 22 Apr 1999 15:10:57 +0000 (15:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@1201 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/dof_renumbering.h
deal.II/deal.II/source/dofs/dof_renumbering.cc

index 1b8bc3d71d60dfd940d5122dd3d39896cdd496b6..4e8beacb8759e4985d252d03428a420d2a1160b1 100644 (file)
  * may be difficult, however, and in many cases will not justify the effort.
  *
  *
- * \subsection{Multigrid DoF numbering}
+ * \section{Componentwise numbering}
+ *
+ * For finite elements composed of several base elements using the #FESystem#
+ * class, or for elements which provide several components themselves, it
+ * may be of interest to sort the DoF indices by component. This will then
+ * bring out the block matrix structure, since otherwise the degrees of freedom
+ * are numbered cell-wise without taking into account that they may belong to
+ * different components.
+ *
+ * This kind of numbering may be obtained by calling the #component_wise# function
+ * of this class. Since it does not touch the order of indices within each, it
+ * may be worthwhile to first renumber using the Cuthill-McKee or a similar
+ * algorithm and afterwards renumbering component-wise. This will bring out the
+ * matrix structure and additionally have a good numbering within each block.
+ *
+ *
+ * \section{Multigrid DoF numbering}
  *
  * Most algorithms also work on multigrid degree of freedom numberings. Refer
  * to the actual function declarations to get more information on this.
@@ -185,11 +201,52 @@ class DoFRenumbering
                               const bool         reversed_numbering = false,
                               const vector<int> &starting_indices   = vector<int> ());
 
+                                    /**
+                                     * Sort the degrees of freedom by
+                                     * component. The numbering within
+                                     * each component is not touched,
+                                     * so a degree of freedom with index
+                                     * $i$, belonging to some component,
+                                     * and another degree of freedom
+                                     * with index $j$ belonging to the same
+                                     * component will be assigned new
+                                     * indices $n(i)$ and $n(j)$ with
+                                     * $n(i)<n(j)$ if $i<j$ and
+                                     * $n(i)>n(j)$ if $i>j$.
+                                     *
+                                     * You may want to give the order in
+                                     * which the components are to be ordered
+                                     * (e.g. if the second argument contains
+                                     * the numbers #(0, 3, 2, 1)#, then all
+                                     * indices of component #0# will be
+                                     * before those of component #3#, before
+                                     * those of component #2#, ...). The
+                                     * length of this list has to be the
+                                     * same as the number of components
+                                     * in the finite element, and has to
+                                     * contain all numbers counted from
+                                     * zero onwards. If
+                                     * you ommit this argument, the same
+                                     * order as given by the finite element
+                                     * is used.
+                                     *
+                                     * For finite elements with only one
+                                     * component, this function is the
+                                     * identity operation.
+                                     */
+    template <int dim>
+    static void component_wise (DoFHandler<dim>            &dof_handler,
+                               const vector<unsigned int> &component_order = vector<unsigned int>());
+
     
                                     /**
                                      * Exception
                                      */
     DeclException0 (ExcRenumberingIncomplete);
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcInvalidComponentOrder);
 };
 
 
index 8a50499ae900d2aebc9e5063207f0ffd6b686da5..dfeaf76f9460dbbcd2336c8721f3995834136a02 100644 (file)
@@ -1,8 +1,11 @@
 /* $Id$ */
 
+#include <grid/dof_accessor.h>
+#include <grid/tria_iterator.h>
 #include <grid/dof.h>
 #include <grid/mg_dof.h>
 #include <grid/dof_constraints.h>
+#include <fe/fe.h>
 #include <numerics/dof_renumbering.h>
 #include <lac/sparsematrix.h>
 
@@ -361,6 +364,102 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler<dim>      &dof_handler,
 
 
 
+template <int dim>
+void DoFRenumbering::component_wise (DoFHandler<dim>            &dof_handler,
+                                    const vector<unsigned int> &component_order_arg)
+{
+  const unsigned int total_dofs = dof_handler.get_fe().total_dofs;
+
+                                  // do nothing if the FE has only
+                                  // one component
+  if (dof_handler.get_fe().n_components == 1)
+    return;
+  
+  vector<unsigned int> component_order (component_order_arg);
+  if (component_order.size() == 0)
+    for (unsigned int i=0; i<dof_handler.get_fe().n_components; ++i)
+      component_order.push_back (i);
+
+                                  // check whether the component list has
+                                  // the right length and contains all
+                                  // component numbers
+  Assert (component_order.size() == dof_handler.get_fe().n_components,
+         ExcInvalidComponentOrder());
+  for (unsigned int i=0; i<dof_handler.get_fe().n_components; ++i)
+    Assert (find (component_order.begin(), component_order.end(), i)
+           != component_order.end(),
+           ExcInvalidComponentOrder ());
+
+                                  // vector to hold the dof indices on
+                                  // the cell we visit at a time
+  vector<int> local_dof_indices(total_dofs);
+                                  // prebuilt list to which component
+                                  // a given dof on a cell belongs
+  vector<unsigned int> component_list (total_dofs);
+  for (unsigned int i=0; i<component_list.size(); ++i)
+    component_list[i] = dof_handler.get_fe().system_to_component_index(i).first;
+  
+                                  // set up a map where for each component
+                                  // the respective degrees of freedom are
+                                  // collected.
+                                  //
+                                  // note that this map is sorted by component
+                                  // but that within each component it is NOT
+                                  // sorted by dof index. note also that some
+                                  // dof indices are entered multiply, so we
+                                  // will have to take care of that
+  vector<vector<int> > component_to_dof_map (dof_handler.get_fe().n_components);
+  for (typename DoFHandler<dim>::active_cell_iterator cell=dof_handler.begin_active();
+       cell!=dof_handler.end(); ++cell)
+    {
+                                      // on each cell: get dof indices
+                                      // and insert them into the global
+                                      // list using their component
+      cell->get_dof_indices (local_dof_indices);
+      for (unsigned int i=0; i<total_dofs; ++i)
+       component_to_dof_map[component_list[i]].push_back (local_dof_indices[i]);
+    };
+  
+                                  // now we've got all indices sorted into
+                                  // buckets labelled with their component
+                                  // number. we've only got to traverse this
+                                  // list and assign the new indices
+                                  //
+                                  // however, we first want to sort the
+                                  // indices entered into the buckets to
+                                  // preserve the order within each component
+                                  // and during this also remove duplicate
+                                  // entries
+  for (unsigned int component=0; component<dof_handler.get_fe().n_components; ++component)
+    {
+      sort (component_to_dof_map[component].begin(),
+           component_to_dof_map[component].end());
+      component_to_dof_map[component].erase (unique (component_to_dof_map[component].begin(),
+                                                    component_to_dof_map[component].end()),
+                                            component_to_dof_map[component].end());
+    };
+  
+  int next_free_index = 0;
+  vector<int> new_indices (dof_handler.n_dofs(), -1);
+  for (unsigned int component=0; component<dof_handler.get_fe().n_components; ++component)
+    {
+      const typename vector<int>::const_iterator
+       begin_of_component = component_to_dof_map[component].begin(),
+       end_of_component   = component_to_dof_map[component].end();
+            
+      for (typename vector<int>::const_iterator dof_index = begin_of_component;
+          dof_index != end_of_component; ++dof_index)
+       new_indices[*dof_index] = next_free_index++;
+    };
+
+  Assert (next_free_index == static_cast<int>(dof_handler.n_dofs()),
+         ExcInternalError());
+
+  dof_handler.renumber_dofs (new_indices);
+};
+
+  
+  
 
 
 
@@ -377,3 +476,6 @@ void DoFRenumbering::Cuthill_McKee (MGDoFHandler<deal_II_dimension> &dof_handler
                                    const bool                       reversed_numbering,
                                    const vector<int>               &starting_indices);
 
+template
+void DoFRenumbering::component_wise (DoFHandler<deal_II_dimension> &dof_handler,
+                                    const vector<unsigned int>    &component_order_arg);

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.