]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Make DoFTools::extract_dofs and DoFTools::count_dofs_per_component work with non...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 18 Oct 2002 14:55:42 +0000 (14:55 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 18 Oct 2002 14:55:42 +0000 (14:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@6690 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc
tests/fe/Makefile
tests/fe/anna_3.cc [new file with mode: 0644]
tests/fe/anna_3.checked [new file with mode: 0644]

index e8409aaab0ad8d182c6049ef587abf9074ab6c92..336899f67a4d2ea83323a1e8e6bbe11d592867e6 100644 (file)
@@ -559,25 +559,53 @@ class DoFTools
                                   const unsigned int     component = 0);
 
                                     /**
-                                     * Extract the indices of the degrees
-                                     * of freedom belonging to certain
-                                     * components. The bit vector @p{select}
-                                     * defines, which components of an
-                                     * @ref{FESystem} are to be extracted
-                                     * from the @ref{DoFHandler} @p{dof}. The
-                                     * respective entries in @p{selected_dofs}
-                                     * are then flagged @p{true}, while all
-                                     * others are set to @p{false}.
+                                     * Extract the indices of the
+                                     * degrees of freedom belonging
+                                     * to certain vector components
+                                     * of a vector-valued finite
+                                     * element. The bit vector
+                                     * @p{select} defines, which
+                                     * components of an
+                                     * @ref{FESystem} are to be
+                                     * extracted from the
+                                     * @ref{DoFHandler} @p{dof}. The
+                                     * entries in @p{selected_dofs}
+                                     * corresponding to degrees of
+                                     * freedom belonging to these
+                                     * components are then flagged
+                                     * @p{true}, while all others are
+                                     * set to @p{false}.
                                      *
-                                     * The size of @p{component_select}
-                                     * shall equal the number of
-                                     * components in the finite
-                                     * element used by @p{dof}. The
-                                     * size of @p{selected_dofs} shall
-                                     * equal
+                                     * The size of
+                                     * @p{component_select} shall
+                                     * equal the number of components
+                                     * in the finite element used by
+                                     * @p{dof}. The size of
+                                     * @p{selected_dofs} shall equal
                                      * @p{dof_handler.n_dofs()}. Previous
                                      * contents of this array or
                                      * overwritten.
+                                     *
+                                     * If the finite element under
+                                     * consideration is not
+                                     * primitive, that is some or all
+                                     * of its shape functions are
+                                     * non-zero in more than one
+                                     * vector component (which holds,
+                                     * for example, for Nedelec or
+                                     * Raviart-Thomas elements), then
+                                     * shape functions cannot be
+                                     * associated with a single
+                                     * vector component. In this
+                                     * case, if @em{one} shape vector
+                                     * component of this element is
+                                     * flagged in
+                                     * @p{component_select}, then
+                                     * this is equivalent to
+                                     * selecting @em{all} vector
+                                     * components corresponding to
+                                     * this non-primitive base
+                                     * element.
                                      */
     template <int dim>
     static void
@@ -705,8 +733,8 @@ class DoFTools
                                      * number belong to each
                                      * components. If the number of
                                      * components the finite element
-                                     * has (i.e. you only have one
-                                     * scalar variable), then the
+                                     * has is one (i.e. you only have
+                                     * one scalar variable), then the
                                      * number in this component
                                      * obviously equals the total
                                      * number of degrees of
@@ -715,6 +743,23 @@ class DoFTools
                                      * needs to equal the total
                                      * number.
                                      *
+                                     * However, the last statement
+                                     * does not hold true if the
+                                     * finite element is not
+                                     * primitive, i.e. some or all of
+                                     * its shape functions are
+                                     * non-zero in more than one
+                                     * vector component. This
+                                     * applies, for example, to the
+                                     * Nedelec or Raviart-Thomas
+                                     * elements. In this case, a
+                                     * degree of freedom is counted
+                                     * in each component in which it
+                                     * is non-zero, so that the sum
+                                     * mentioned above is greater
+                                     * than the total number of
+                                     * degrees of freedom.
+                                     *
                                      * The result is returned in the
                                      * last argument.
                                      */
index 03594a90803bd92589a569cccb388dfbde576316..080219ddaf3abcf490636f3cb83ef6a2f26aa5dc 100644 (file)
@@ -1050,22 +1050,87 @@ DoFTools::extract_dofs (const DoFHandler<dim>      &dof,
   Assert(selected_dofs.size() == dof.n_dofs(),
         ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs()));
 
+                                   // two special cases: no component
+                                   // is selected, and all components
+                                   // are selected, both rather
+                                   // stupid, but easy to catch
+  if (std::count (component_select.begin(), component_select.end(), true)
+      == 0)
+    {
+      std::fill_n (selected_dofs.begin(), dof.n_dofs(), false);
+      return;
+    };
+  if (std::count (component_select.begin(), component_select.end(), true)
+      == static_cast<signed int>(component_select.size()))
+    {
+      std::fill_n (selected_dofs.begin(), dof.n_dofs(), true);
+      return;
+    };
+  
+  
                                   // preset all values by false
-  fill_n (selected_dofs.begin(), dof.n_dofs(), false);
+  std::fill_n (selected_dofs.begin(), dof.n_dofs(), false);
+
+                                   // next set up a table for the
+                                   // degrees of freedom on each of
+                                   // the cells whether it is
+                                   // something interesting or not
+  std::vector<bool> local_selected_dofs (fe.dofs_per_cell, false);
+  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+    if (fe.is_primitive(i))
+      local_selected_dofs[i]
+        = component_select[fe.system_to_component_index(i).first];
+    else
+                                       // if this shape function is
+                                       // not primitive, then we have
+                                       // to work harder. we have to
+                                       // find out whether _any_ of
+                                       // the vector components of
+                                       // this element is selected or
+                                       // not
+                                       //
+                                       // to do so, get the first and
+                                       // last vector components of
+                                       // the base element to which
+                                       // the local dof with index i
+                                       // belongs
+      {
+        unsigned int first_comp = 0;
+        const unsigned int this_base = fe.system_to_base_index(i).first.first;
+        const unsigned int this_multiplicity
+          = fe.system_to_base_index(i).first.second;
+        
+        for (unsigned int b=0; b<this_base; ++b)
+          first_comp += fe.base_element(b).n_components() *
+                        fe.element_multiplicity(b);
+        for (unsigned int m=0; m<this_multiplicity; ++m)
+          first_comp += fe.base_element(this_base).n_components();
+        const unsigned int end_comp = first_comp +
+                                      fe.base_element(this_base).n_components();
+
+        Assert (first_comp < fe.n_components(), ExcInternalError());
+        Assert (end_comp <= fe.n_components(),  ExcInternalError());
+
+                                         // now check whether any of
+                                         // the components in between
+                                         // is set
+        for (unsigned int c=first_comp; c<end_comp; ++c)
+          if (component_select[c] == true)
+            {
+              local_selected_dofs[i] = true;
+              break;
+            };
+      };  
   
+                                   // then loop over all cells and do
+                                   // the work
   std::vector<unsigned int> indices(fe.dofs_per_cell);
-  
   typename DoFHandler<dim>::active_cell_iterator c;
-  for (c = dof.begin_active() ; c != dof.end() ; ++ c)
+  for (c=dof.begin_active(); c!=dof.end(); ++ c)
     {
       c->get_dof_indices(indices);
-      for (unsigned int i=0;i<fe.dofs_per_cell;++i)
-       {
-         const unsigned int component = fe.system_to_component_index(i).first;
-
-         if (component_select[component] == true)
-           selected_dofs[indices[i]] = true;
-       }
+      for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+        selected_dofs[indices[i]] = local_selected_dofs[i];
     }
 }
 
@@ -1084,21 +1149,86 @@ DoFTools::extract_level_dofs(const unsigned int       level,
   Assert(selected_dofs.size() == dof.n_dofs(level),
         ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs(level)));
 
-                                  // preset all values by false
-  fill_n (selected_dofs.begin(), dof.n_dofs(level), false);
+                                   // two special cases: no component
+                                   // is selected, and all components
+                                   // are selected, both rather
+                                   // stupid, but easy to catch
+  if (std::count (component_select.begin(), component_select.end(), true)
+      == 0)
+    {
+      std::fill_n (selected_dofs.begin(), dof.n_dofs(level), false);
+      return;
+    };
+  if (std::count (component_select.begin(), component_select.end(), true)
+      == static_cast<signed int>(component_select.size()))
+    {
+      std::fill_n (selected_dofs.begin(), dof.n_dofs(level), true);
+      return;
+    };
 
-  std::vector<unsigned int> indices(fe.dofs_per_cell);
-  
+                                  // preset all values by false
+  std::fill_n (selected_dofs.begin(), dof.n_dofs(level), false);
+
+                                   // next set up a table for the
+                                   // degrees of freedom on each of
+                                   // the cells whether it is
+                                   // something interesting or not
+  std::vector<bool> local_selected_dofs (fe.dofs_per_cell, false);
+  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+    if (fe.is_primitive(i))
+      local_selected_dofs[i]
+        = component_select[fe.system_to_component_index(i).first];
+    else
+                                       // if this shape function is
+                                       // not primitive, then we have
+                                       // to work harder. we have to
+                                       // find out whether _any_ of
+                                       // the vector components of
+                                       // this element is selected or
+                                       // not
+                                       //
+                                       // to do so, get the first and
+                                       // last vector components of
+                                       // the base element to which
+                                       // the local dof with index i
+                                       // belongs
+      {
+        unsigned int first_comp = 0;
+        const unsigned int this_base = fe.system_to_base_index(i).first.first;
+        const unsigned int this_multiplicity
+          = fe.system_to_base_index(i).first.second;
+        
+        for (unsigned int b=0; b<this_base; ++b)
+          first_comp += fe.base_element(b).n_components() *
+                        fe.element_multiplicity(b);
+        for (unsigned int m=0; m<this_multiplicity; ++m)
+          first_comp += fe.base_element(this_base).n_components();
+        const unsigned int end_comp = first_comp +
+                                      fe.base_element(this_base).n_components();
+
+        Assert (first_comp < fe.n_components(), ExcInternalError());
+        Assert (end_comp <= fe.n_components(),  ExcInternalError());
+
+                                         // now check whether any of
+                                         // the components in between
+                                         // is set
+        for (unsigned int c=first_comp; c<end_comp; ++c)
+          if (component_select[c] == true)
+            {
+              local_selected_dofs[i] = true;
+              break;
+            };
+      };
+
+                                   // then loop over all cells and do
+                                   // work
+  std::vector<unsigned int> indices(fe.dofs_per_cell);  
   typename MGDoFHandler<dim>::cell_iterator c;
   for (c = dof.begin(level) ; c != dof.end(level) ; ++ c)
     {
       c->get_mg_dof_indices(indices);
-      for (unsigned int i=0;i<fe.dofs_per_cell;++i)
-       {
-         const unsigned int component = fe.system_to_component_index(i).first;
-         if (component_select[component]  == true)
-           selected_dofs[indices[i]] = true;
-       }
+      for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+        selected_dofs[indices[i]] = local_selected_dofs[i];
     }
 }
 #endif
@@ -1399,8 +1529,9 @@ DoFTools::extract_subdomain_dofs (const DoFHandler<dim> &dof_handler,
 
 template <int dim>
 void
-DoFTools::count_dofs_per_component (const DoFHandler<dim>     &dof_handler,
-                                   std::vector<unsigned int> &dofs_per_component)
+DoFTools::
+count_dofs_per_component (const DoFHandler<dim>     &dof_handler,
+                          std::vector<unsigned int> &dofs_per_component)
 {
   const unsigned int n_components = dof_handler.get_fe().n_components();
   dofs_per_component.resize (n_components);
@@ -1418,17 +1549,19 @@ DoFTools::count_dofs_per_component (const DoFHandler<dim>     &dof_handler,
                                   // otherwise determine the number
                                   // of dofs in each component
                                   // separately. do so in parallel
-  std::vector<std::vector<bool> > dofs_in_component (n_components,
-                                                    std::vector<bool>(dof_handler.n_dofs(),
-                                                                      false));
-  std::vector<std::vector<bool> > component_select (n_components,
-                                                   std::vector<bool>(n_components, false));
+  std::vector<std::vector<bool> >
+    dofs_in_component (n_components,
+                       std::vector<bool>(dof_handler.n_dofs(), false));
+  std::vector<std::vector<bool> >
+    component_select (n_components,
+                      std::vector<bool>(n_components, false));
   Threads::ThreadManager thread_manager;
   for (unsigned int i=0; i<n_components; ++i)
     {
       void (*fun_ptr) (const DoFHandler<dim>      &,
                       const std::vector<bool>    &,
-                      std::vector<bool>          &) = &DoFTools::template extract_dofs<dim>;
+                      std::vector<bool>          &)
+        = &DoFTools::template extract_dofs<dim>;
       component_select[i][i] = true;
       Threads::spawn (thread_manager,
                      Threads::encapsulate (fun_ptr)
@@ -1443,10 +1576,15 @@ DoFTools::count_dofs_per_component (const DoFHandler<dim>     &dof_handler,
                                       dofs_in_component[i].end(),
                                       true);
 
-                                  // finally sanity check
-  Assert (std::accumulate (dofs_per_component.begin(), dofs_per_component.end(), 0U)
-         ==
-         dof_handler.n_dofs(),
+                                  // finally sanity check. this is
+                                  // only valid if the finite element
+                                  // is actually primitive, so
+                                  // exclude other elements from this
+  Assert (!dof_handler.get_fe().is_primitive()
+          ||
+          (std::accumulate (dofs_per_component.begin(),
+                            dofs_per_component.end(), 0U)
+           == dof_handler.n_dofs()),
          ExcInternalError());
          
 };
index d6fcfbdb0bef7ce040efa4285842d559d073ca45..0e798005c1c6d1afa3e2a7703956c249ce26832e 100644 (file)
@@ -41,12 +41,13 @@ transfer.exe                : transfer.go            $(libraries)
 up_and_down.exe         : up_and_down.go        $(libraries)
 anna_1.exe              : anna_1.go              $(libraries)
 anna_2.exe              : anna_2.go              $(libraries)
+anna_3.exe              : anna_3.go              $(libraries)
 
 tests = derivatives fe_data_test fe_traits fe_tools_test mapping \
         mapping_c1 shapes numbering mapping_q1_eulerian \
         transfer internals \
         non_primitive_1 non_primitive_2 \
-        nedelec nedelec_2 nedelec_3 up_and_down anna_1 anna_2
+        nedelec nedelec_2 nedelec_3 up_and_down anna_1 anna_2 anna_3
 
 ############################################################
 
diff --git a/tests/fe/anna_3.cc b/tests/fe/anna_3.cc
new file mode 100644 (file)
index 0000000..dd1c25d
--- /dev/null
@@ -0,0 +1,135 @@
+//----------------------------  anna_3.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2002 by the deal.II authors and Anna Schneebeli
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  anna_3.cc  ---------------------------
+
+
+// check some things about Nedelec elements, here that
+// DoFTools::component_select and DoFTools::count_dofs_per_component
+// works
+//
+// this program is a modified version of one by Anna Schneebeli,
+// University of Basel
+
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <dofs/dof_handler.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_system.h>              
+#include <fe/fe_q.h>
+#include <fe/fe_nedelec.h>
+#include <fe/fe_base.h>
+#include <dofs/dof_tools.h>
+#include <numerics/dof_renumbering.h>
+#include <iostream>
+#include <fstream>
+
+
+template <int dim>
+class SystemTest 
+{
+  public:
+    SystemTest ();
+    void run ();    
+                                   
+  private:
+    void make_grid_and_dofs ();
+    void check ();
+
+                                   
+    Triangulation<dim>     triangulation;
+    FESystem<dim>          fe;
+    DoFHandler<dim>        dof_handler;
+
+                                  
+};
+
+template <int dim>
+SystemTest<dim>::SystemTest () :
+                fe (FE_Nedelec<dim>(1), 2,
+                    FE_Q<dim>(1), 1),
+               dof_handler (triangulation)
+{};
+
+
+template <int dim>
+void SystemTest<dim>::make_grid_and_dofs ()
+{
+                                 
+  GridGenerator::hyper_cube (triangulation, -1, 1);
+  triangulation.refine_global (0);
+  deallog << "Number of active cells: " << triangulation.n_active_cells()
+          << std::endl;
+  deallog << "Total number of cells: " << triangulation.n_cells()
+          << std::endl;
+                                 
+  dof_handler.distribute_dofs (fe);
+  deallog << "Number of degrees of freedom: " << dof_handler.n_dofs()
+          << std::endl;
+                                 
+};
+
+
+template <int dim>
+void SystemTest<dim>::check () 
+{
+  for (unsigned int c=0; c<fe.n_components(); ++c)
+    {
+      deallog << "Checking for component " << c << std::endl;
+      std::vector<bool> x(fe.n_components(), false);
+      x[c] = true;
+      std::vector<bool> sel(dof_handler.n_dofs());
+      DoFTools::extract_dofs (dof_handler, x, sel);
+
+      for (unsigned int i=0; i<sel.size(); ++i)
+        if (sel[i])
+          deallog << "  DoF " << i << std::endl;
+    };
+
+  std::vector<unsigned int> dofs_per_component (fe.n_components(), 0U);
+  DoFTools::count_dofs_per_component (dof_handler, dofs_per_component);
+  deallog << "DoFs per component: ";
+  for (unsigned int i=0; i<fe.n_components(); ++i)
+    deallog << dofs_per_component[i] << ' ';
+  deallog << std::endl;
+};
+
+
+template <int dim>
+void SystemTest<dim>::run () 
+{
+  deallog << "************* " << dim << "D *************" << std::endl;
+  make_grid_and_dofs ();
+  check();
+
+                                   // renumber degrees of freedom and
+                                   // try again
+  deallog << std::endl << "*** Renumbering ***" << std::endl;
+  DoFRenumbering::component_wise (dof_handler);
+  check ();
+};
+
+    
+
+int main () 
+{
+  std::ofstream logfile("anna_3.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  SystemTest<2>().run();
+  SystemTest<3>().run();  
+  return 0;
+};
diff --git a/tests/fe/anna_3.checked b/tests/fe/anna_3.checked
new file mode 100644 (file)
index 0000000..6cea739
--- /dev/null
@@ -0,0 +1,241 @@
+
+DEAL::************* 2D *************
+DEAL::Number of active cells: 1
+DEAL::Total number of cells: 1
+DEAL::Number of degrees of freedom: 12
+DEAL::Checking for component 0
+DEAL::  DoF 4
+DEAL::  DoF 6
+DEAL::  DoF 8
+DEAL::  DoF 10
+DEAL::Checking for component 1
+DEAL::  DoF 4
+DEAL::  DoF 6
+DEAL::  DoF 8
+DEAL::  DoF 10
+DEAL::Checking for component 2
+DEAL::  DoF 5
+DEAL::  DoF 7
+DEAL::  DoF 9
+DEAL::  DoF 11
+DEAL::Checking for component 3
+DEAL::  DoF 5
+DEAL::  DoF 7
+DEAL::  DoF 9
+DEAL::  DoF 11
+DEAL::Checking for component 4
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::DoFs per component: 4 4 4 4 4 
+DEAL::
+DEAL::*** Renumbering ***
+DEAL::Checking for component 0
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::Checking for component 1
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::Checking for component 2
+DEAL::  DoF 4
+DEAL::  DoF 5
+DEAL::  DoF 6
+DEAL::  DoF 7
+DEAL::Checking for component 3
+DEAL::  DoF 4
+DEAL::  DoF 5
+DEAL::  DoF 6
+DEAL::  DoF 7
+DEAL::Checking for component 4
+DEAL::  DoF 8
+DEAL::  DoF 9
+DEAL::  DoF 10
+DEAL::  DoF 11
+DEAL::DoFs per component: 4 4 4 4 4 
+DEAL::************* 3D *************
+DEAL::Number of active cells: 1
+DEAL::Total number of cells: 1
+DEAL::Number of degrees of freedom: 32
+DEAL::Checking for component 0
+DEAL::  DoF 8
+DEAL::  DoF 10
+DEAL::  DoF 12
+DEAL::  DoF 14
+DEAL::  DoF 16
+DEAL::  DoF 18
+DEAL::  DoF 20
+DEAL::  DoF 22
+DEAL::  DoF 24
+DEAL::  DoF 26
+DEAL::  DoF 28
+DEAL::  DoF 30
+DEAL::Checking for component 1
+DEAL::  DoF 8
+DEAL::  DoF 10
+DEAL::  DoF 12
+DEAL::  DoF 14
+DEAL::  DoF 16
+DEAL::  DoF 18
+DEAL::  DoF 20
+DEAL::  DoF 22
+DEAL::  DoF 24
+DEAL::  DoF 26
+DEAL::  DoF 28
+DEAL::  DoF 30
+DEAL::Checking for component 2
+DEAL::  DoF 8
+DEAL::  DoF 10
+DEAL::  DoF 12
+DEAL::  DoF 14
+DEAL::  DoF 16
+DEAL::  DoF 18
+DEAL::  DoF 20
+DEAL::  DoF 22
+DEAL::  DoF 24
+DEAL::  DoF 26
+DEAL::  DoF 28
+DEAL::  DoF 30
+DEAL::Checking for component 3
+DEAL::  DoF 9
+DEAL::  DoF 11
+DEAL::  DoF 13
+DEAL::  DoF 15
+DEAL::  DoF 17
+DEAL::  DoF 19
+DEAL::  DoF 21
+DEAL::  DoF 23
+DEAL::  DoF 25
+DEAL::  DoF 27
+DEAL::  DoF 29
+DEAL::  DoF 31
+DEAL::Checking for component 4
+DEAL::  DoF 9
+DEAL::  DoF 11
+DEAL::  DoF 13
+DEAL::  DoF 15
+DEAL::  DoF 17
+DEAL::  DoF 19
+DEAL::  DoF 21
+DEAL::  DoF 23
+DEAL::  DoF 25
+DEAL::  DoF 27
+DEAL::  DoF 29
+DEAL::  DoF 31
+DEAL::Checking for component 5
+DEAL::  DoF 9
+DEAL::  DoF 11
+DEAL::  DoF 13
+DEAL::  DoF 15
+DEAL::  DoF 17
+DEAL::  DoF 19
+DEAL::  DoF 21
+DEAL::  DoF 23
+DEAL::  DoF 25
+DEAL::  DoF 27
+DEAL::  DoF 29
+DEAL::  DoF 31
+DEAL::Checking for component 6
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::  DoF 4
+DEAL::  DoF 5
+DEAL::  DoF 6
+DEAL::  DoF 7
+DEAL::DoFs per component: 12 12 12 12 12 12 8 
+DEAL::
+DEAL::*** Renumbering ***
+DEAL::Checking for component 0
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::  DoF 4
+DEAL::  DoF 5
+DEAL::  DoF 6
+DEAL::  DoF 7
+DEAL::  DoF 8
+DEAL::  DoF 9
+DEAL::  DoF 10
+DEAL::  DoF 11
+DEAL::Checking for component 1
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::  DoF 4
+DEAL::  DoF 5
+DEAL::  DoF 6
+DEAL::  DoF 7
+DEAL::  DoF 8
+DEAL::  DoF 9
+DEAL::  DoF 10
+DEAL::  DoF 11
+DEAL::Checking for component 2
+DEAL::  DoF 0
+DEAL::  DoF 1
+DEAL::  DoF 2
+DEAL::  DoF 3
+DEAL::  DoF 4
+DEAL::  DoF 5
+DEAL::  DoF 6
+DEAL::  DoF 7
+DEAL::  DoF 8
+DEAL::  DoF 9
+DEAL::  DoF 10
+DEAL::  DoF 11
+DEAL::Checking for component 3
+DEAL::  DoF 12
+DEAL::  DoF 13
+DEAL::  DoF 14
+DEAL::  DoF 15
+DEAL::  DoF 16
+DEAL::  DoF 17
+DEAL::  DoF 18
+DEAL::  DoF 19
+DEAL::  DoF 20
+DEAL::  DoF 21
+DEAL::  DoF 22
+DEAL::  DoF 23
+DEAL::Checking for component 4
+DEAL::  DoF 12
+DEAL::  DoF 13
+DEAL::  DoF 14
+DEAL::  DoF 15
+DEAL::  DoF 16
+DEAL::  DoF 17
+DEAL::  DoF 18
+DEAL::  DoF 19
+DEAL::  DoF 20
+DEAL::  DoF 21
+DEAL::  DoF 22
+DEAL::  DoF 23
+DEAL::Checking for component 5
+DEAL::  DoF 12
+DEAL::  DoF 13
+DEAL::  DoF 14
+DEAL::  DoF 15
+DEAL::  DoF 16
+DEAL::  DoF 17
+DEAL::  DoF 18
+DEAL::  DoF 19
+DEAL::  DoF 20
+DEAL::  DoF 21
+DEAL::  DoF 22
+DEAL::  DoF 23
+DEAL::Checking for component 6
+DEAL::  DoF 24
+DEAL::  DoF 25
+DEAL::  DoF 26
+DEAL::  DoF 27
+DEAL::  DoF 28
+DEAL::  DoF 29
+DEAL::  DoF 30
+DEAL::  DoF 31
+DEAL::DoFs per component: 12 12 12 12 12 12 8 

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.