]> https://gitweb.dealii.org/ - dealii.git/commitdiff
introduce component renumbering and grouping
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 8 Jul 2003 10:56:42 +0000 (10:56 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 8 Jul 2003 10:56:42 +0000 (10:56 +0000)
git-svn-id: https://svn.dealii.org/trunk@7851 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_renumbering.h
deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/include/multigrid/mg_dof_tools.h
deal.II/deal.II/include/multigrid/mg_transfer.h
deal.II/deal.II/include/multigrid/mg_transfer.templates.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/multigrid/mg_dof_tools.cc
deal.II/deal.II/source/multigrid/mg_transfer_block.cc

index da9735c5cb3f6ca9378a3c975eea80b39f15b5d4..d3fa5741a20da24fd2e3909faa171a7f3028e5d8 100644 (file)
@@ -279,7 +279,7 @@ class DoFRenumbering
                                      * different way than suggested
                                      * by the @p{FESystem} object you
                                      * use. To this end, Set up the
-                                     * vector @p{component_order}
+                                     * vector @p{target_component}
                                      * such that the entry at index
                                      * @p{i} denotes the number of
                                      * the target component for dofs
@@ -317,13 +317,27 @@ class DoFRenumbering
     template <int dim>
     static void
     component_wise (DoFHandler<dim>                 &dof_handler,
-                   const std::vector<unsigned int> &component_order = std::vector<unsigned int>());
+                   const std::vector<unsigned int> &target_component = std::vector<unsigned int>());
+
+                                    /**
+                                     * Sort the degrees of freedom by
+                                     * component. It does the same
+                                     * thing as the above function,
+                                     * only that it does this for one
+                                     * level of a multi-level
+                                     * discretization.
+                                     */
+    template <int dim>
+    static void
+    component_wise (MGDoFHandler<dim>&               dof_handler,
+                   unsigned int                     level,
+                   const std::vector<unsigned int>& target_component = std::vector<unsigned int>());
 
                                     /**
                                      * Computes the renumbering
                                      * vector needed by the
                                      * @p{component_wise}
-                                     * function. Does not perform the
+                                     * functions. Does not perform the
                                      * renumbering on the
                                      * @p{DoFHandler} dofs but
                                      * returns the renumbering
@@ -334,22 +348,8 @@ class DoFRenumbering
     compute_component_wise (std::vector<unsigned int>& new_index,
                            ITERATOR& start,
                            const ENDITERATOR& end,
-                           const std::vector<unsigned int> &component_order = std::vector<unsigned int>());
+                           const std::vector<unsigned int> &target_component);
     
-                                    /**
-                                     * Sort the degrees of freedom by
-                                     * component. It does the same
-                                     * thing as the above function,
-                                     * only that it does this for one
-                                     * level of a multi-level
-                                     * discretization.
-                                     */
-    template <int dim>
-    static void
-    component_wise (MGDoFHandler<dim>&               dof_handler,
-                   unsigned int                     level,
-                   const std::vector<unsigned int>& component_order = std::vector<unsigned int>());
-
                                     /**
                                      * Cell-wise renumbering for DG
                                      * elements.  This function takes
index fd128f439095f238cd8e76417ea318298dd73ddf..924f5ccb8a391bd508786a2fb91c03b8c4a9051e 100644 (file)
@@ -775,7 +775,7 @@ class DoFTools
                                      * Count how many degrees of
                                      * freedom out of the total
                                      * number belong to each
-                                     * components. If the number of
+                                     * component. If the number of
                                      * components the finite element
                                      * has is one (i.e. you only have
                                      * one scalar variable), then the
@@ -804,13 +804,26 @@ class DoFTools
                                      * than the total number of
                                      * degrees of freedom.
                                      *
-                                     * The result is returned in the
-                                     * last argument.
+                                     * The additional optional
+                                     * argument @p{target_component}
+                                     * allows for a re-sorting and
+                                     * grouping of components. To
+                                     * this end, it contains for each
+                                     * component the component number
+                                     * it shall be counted as. Having
+                                     * the same number entered
+                                     * several times sums up several
+                                     * components as the same.
+                                     *
+                                     * The result is returned in
+                                     * @p{dofs_per_component}.
                                      */
     template <int dim>
     static void
-    count_dofs_per_component (const DoFHandler<dim>     &dof_handler,
-                             std::vector<unsigned int> &dofs_per_component);
+    count_dofs_per_component (const DoFHandler<dim>&     dof_handler,
+                             std::vector<unsigned int>& dofs_per_component,
+                             std::vector<unsigned int>  target_component
+                             =std::vector<unsigned int>());
     
                                     /**
                                      * This function can be used when
index 57edad521e9ff47b5d182087f84e9fde0396685e..738b96204d4245a35f8bb9907d1b1985d4dc071b 100644 (file)
@@ -36,7 +36,7 @@ template <class number> class FullMatrix;
  * All member functions are static, so there is no need to create an
  * object of class @p{MGTools}.
  *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999, 2000
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999-2003
  */
 class MGTools
 {
@@ -126,7 +126,9 @@ class MGTools
                                      */
     template <int dim>
       static void count_dofs_per_component (const MGDoFHandler<dim>& mg_dof,
-                                           std::vector<std::vector<unsigned int> >& result);
+                                           std::vector<std::vector<unsigned int> >& result,
+                                           std::vector<unsigned int> target_component
+                                           = std::vector<unsigned int>());
     
     
                                     /**
index ad9c4878717358bfa7b9f22905c61f48599f8720..5add525f41cbc5c8bd3bc36a14b5aa52073dbbaf 100644 (file)
@@ -41,7 +41,7 @@ template <int dim> class MGDoFHandler;
  * once by looping over all cells and storing the result in a matrix for
  * each level, but requires additional memory.
  *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999, 2000, 2001, 2002
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999-2003
  */
 template <typename number>
 class MGTransferPrebuilt : public MGTransferBase<Vector<number> > 
@@ -255,7 +255,7 @@ class MGTransferPrebuilt : public MGTransferBase<Vector<number> >
  * Implementation of matrix generation for @ref{MGTransferBlock} and
  * @p{MGTransferSelected}.
  *
- * @author Guido Kanschat, 2001
+ * @author Guido Kanschat, 2001-2003
  */
 class MGTransferBlockBase
 {
@@ -269,10 +269,17 @@ class MGTransferBlockBase
                                      * components are to be built. By
                                      * default, all matrices are
                                      * built.
+                                     *
+                                     * The last argument
+                                     * @p{target_component} allows
+                                     * grouping of components the
+                                     * same way as in
+                                     * @p{DoFRenumbering::component_wise}.
                                      */
     template <int dim>
-    void build_matrices (const MGDoFHandler<dim> &mg_dof,
-                        std::vector<bool> selected);
+    void build_matrices (const MGDoFHandler<dim>& mg_dof,
+                        const std::vector<bool>& selected,
+                        const std::vector<unsigned int>& target_component);
 
                                   /**
                                    * Flag of selected components.
@@ -518,10 +525,11 @@ class MGTransferBlock : public MGTransferBase<BlockVector<number> >,
 /**
  * Implementation of the @p{MGTransferBase} interface for block
  * matrices and simple vectors. This class uses @ref{MGTransferBlock}
- * selecting a single component. The transfer operators themselves are
- * implemented for simple vectors again.
+ * selecting a single component or grouping several components into a
+ * single block. The transfer operators themselves are implemented for
+ * simple vectors again.
  *
- * @author Guido Kanschat, 2001, 2002
+ * @author Guido Kanschat, 2001, 2002, 2003
  */
 template <typename number>
 class MGTransferSelect : public MGTransferBase<Vector<number> >,
@@ -548,6 +556,29 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
     void build_matrices (const MGDoFHandler<dim> &mg_dof,
                         unsigned int selected);
 
+                                    /**
+                                     * Actually build the prolongation
+                                     * matrices for grouped components.
+                                     *
+                                     * The argument
+                                     * @p{target_component}
+                                     * corresponds to the grouping
+                                     * mechanism of
+                                     * @p{DoFRenumbering::component_wise(...)},
+                                     * which should be used to create
+                                     * the corresponding block
+                                     * structure in matrices and
+                                     * vectors.
+                                     *
+                                     * This function is a front-end
+                                     * for the same function in
+                                     * @ref{MGTransferBlockBase}.
+                                     */
+    template <int dim>
+    void build_matrices (const MGDoFHandler<dim> &mg_dof,
+                        unsigned int selected,
+                        const std::vector<unsigned int>& target_component);
+
                                     /**
                                      * Prolongate a vector from level
                                      * @p{to_level-1} to level
index 0f384aa95863de0c9b63468a94a8131acf18b648..eee2a98470e30a87f0d377ea4257e22f4b3a522c 100644 (file)
@@ -26,8 +26,6 @@
 #include <algorithm>
 
 
-//TODO:[GK] This file is included from nowhere. It should either go away entirely, or be included at least from some .h or .cc file
-
 /* --------------------- MGTransferPrebuilt -------------- */
 
 
index a36065742aca6e30e175ce81c157b57ea4a5524a..9ae23dd5593bd149f0a5bb0df9c833dea531645b 100644 (file)
@@ -1674,11 +1674,26 @@ 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)
+count_dofs_per_component (const DoFHandler<dim>&     dof_handler,
+                          std::vector<unsigned int>& dofs_per_component,
+                         std::vector<unsigned int>  target_component)
 {
   const unsigned int n_components = dof_handler.get_fe().n_components();
   dofs_per_component.resize (n_components);
+  fill (dofs_per_component.begin(), dofs_per_component.end(), 0U);
+  
+                                  // If the empty vector was given as
+                                  // default argument, set up this
+                                  // vector as identity.
+  if (target_component.size()==0)
+    {
+      target_component.resize(n_components);
+      for (unsigned int i=0;i<n_components;++i)
+       target_component[i] = i;
+    }
+  
+  Assert(target_component.size()==n_components,
+        ExcDimensionMismatch(target_component.size(),n_components));
 
                                   // special case for only one
                                   // component. treat this first
@@ -1714,9 +1729,10 @@ count_dofs_per_component (const DoFHandler<dim>     &dof_handler,
 
                                   // next count what we got
   for (unsigned int i=0; i<n_components; ++i)
-    dofs_per_component[i] = std::count(dofs_in_component[i].begin(),
-                                      dofs_in_component[i].end(),
-                                      true);
+    dofs_per_component[target_component[i]]
+      += std::count(dofs_in_component[i].begin(),
+                   dofs_in_component[i].end(),
+                   true);
 
                                   // finally sanity check. this is
                                   // only valid if the finite element
index 505deea3e340629c89646e3372f41a61176491a8..a5a1335eb4d0b773847a68ac6ba4ac5f291b607b 100644 (file)
 
 
 template <int dim, class SparsityPattern>
-void MGTools::make_sparsity_pattern (const MGDoFHandler<dim> &dof,
-                                       SparsityPattern         &sparsity,
-                                       const unsigned int       level)
+void MGTools::make_sparsity_pattern (
+  const MGDoFHandler<dim> &dof,
+  SparsityPattern         &sparsity,
+  const unsigned int       level)
 {
   const unsigned int n_dofs = dof.n_dofs(level);
 
@@ -65,9 +66,10 @@ void MGTools::make_sparsity_pattern (const MGDoFHandler<dim> &dof,
 
 template<int dim, class SparsityPattern>
 void
-MGTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
-                                       SparsityPattern       &sparsity,
-                                       const unsigned int level)
+MGTools::make_flux_sparsity_pattern (
+  const MGDoFHandler<dim> &dof,
+  SparsityPattern       &sparsity,
+  const unsigned int level)
 {
   const unsigned int n_dofs = dof.n_dofs(level);
   
@@ -121,9 +123,10 @@ MGTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
 
 template<int dim, class SparsityPattern>
 void
-MGTools::make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof,
-                                         SparsityPattern       &sparsity,
-                                         const unsigned int level)
+MGTools::make_flux_sparsity_pattern_edge (
+  const MGDoFHandler<dim> &dof,
+  SparsityPattern       &sparsity,
+  const unsigned int level)
 {
   Assert ((level>=1) && (level<dof.get_tria().n_levels()),
          ExcIndexRange(level, 1, dof.get_tria().n_levels()));
@@ -177,14 +180,15 @@ MGTools::make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof,
 
 
 
-//TODO[GK]: Replace FullMatrix by Table<2,bool>
+//TODO[GK]: Replace FullMatrix by Table<2,char>
 template<int dim, class SparsityPattern>
 void
-MGTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
-                                     SparsityPattern       &sparsity,
-                                     const unsigned int level,
-                                     const FullMatrix<double>& int_mask,
-                                     const FullMatrix<double>& flux_mask)
+MGTools::make_flux_sparsity_pattern (
+  const MGDoFHandler<dim> &dof,
+  SparsityPattern       &sparsity,
+  const unsigned int level,
+  const FullMatrix<double>& int_mask,
+  const FullMatrix<double>& flux_mask)
 {
   const unsigned int n_dofs = dof.n_dofs(level);
   const unsigned int n_comp = dof.get_fe().n_components();
@@ -291,8 +295,10 @@ MGTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
 
 template <int dim>
 void
-MGTools::count_dofs_per_component (const MGDoFHandler<dim>& dof_handler,
-                                     std::vector<std::vector<unsigned int> >& result)
+MGTools::count_dofs_per_component (
+  const MGDoFHandler<dim>& dof_handler,
+  std::vector<std::vector<unsigned int> >& result,
+  std::vector<unsigned int> target_component)
 {
   const unsigned int nlevels = dof_handler.get_tria().n_levels();
   
@@ -300,6 +306,17 @@ MGTools::count_dofs_per_component (const MGDoFHandler<dim>& dof_handler,
          ExcDimensionMismatch(result.size(), nlevels));
 
   const unsigned int n_components = dof_handler.get_fe().n_components();
+
+  if (target_component.size() == 0)
+    {
+      target_component.resize(n_components);
+      for (unsigned int i=0;i<n_components;++i)
+       target_component[i] = i;
+    }
+
+  Assert(target_component.size() == n_components,
+        ExcDimensionMismatch(target_component.size(), n_components));
+
   for (unsigned int l=0;l<nlevels;++l)
     {
       result[l].resize (n_components);
@@ -338,9 +355,10 @@ MGTools::count_dofs_per_component (const MGDoFHandler<dim>& dof_handler,
          
                                           // next count what we got
          for (unsigned int i=0; i<n_components; ++i)
-             result[l][i] = std::count(dofs_in_component[i].begin(),
-                                       dofs_in_component[i].end(),
-                                       true);
+             result[l][target_component[i]]
+               = std::count(dofs_in_component[i].begin(),
+                            dofs_in_component[i].end(),
+                            true);
          
                                           // finally sanity check
          Assert (std::accumulate (result[l].begin(),
index bd83460beb8f357460246b741ed6308260a81c29..47344d53b7bd21d92152e5f15fefdf4da949158e 100644 (file)
@@ -29,8 +29,9 @@
 
 template <int dim>
 void MGTransferBlockBase::build_matrices (
-  const MGDoFHandler<dim> &mg_dof,
-  std::vector<bool> select) 
+  const MGDoFHandler<dim>& mg_dof,
+  const std::vector<bool>& select,
+  const std::vector<unsigned int>& target_component)
 {
   const FiniteElement<dim>& fe = mg_dof.get_fe();
   const unsigned int n_components  = fe.n_components();
@@ -248,12 +249,15 @@ template <typename number>
 template <int dim>
 void MGTransferBlock<number>::build_matrices (
   const MGDoFHandler<dim> &mg_dof,
-  std::vector<bool> select) 
+  std::vector<bool> select)
 {
   if (select.size() == 0)
     select = std::vector<bool> (mg_dof.get_fe().n_components(), true);
-
-  MGTransferBlockBase::build_matrices (mg_dof, select);
+  std::vector<unsigned int> target_component(mg_dof.get_fe().n_components());
+  for (unsigned int i=0;i<target_component.size();++i)
+    target_component[i] = i;
+  
+  MGTransferBlockBase::build_matrices (mg_dof, select, target_component);
 }
 
 
@@ -266,8 +270,11 @@ void MGTransferSelect<number>::build_matrices (
   selected = select;
   std::vector<bool> s(mg_dof.get_fe().n_components(), false);
   s[select] = true;
+  std::vector<unsigned int> target_component(mg_dof.get_fe().n_components());
+  for (unsigned int i=0;i<target_component.size();++i)
+    target_component[i] = i;
 
-  MGTransferBlockBase::build_matrices (mg_dof, s);
+  MGTransferBlockBase::build_matrices (mg_dof, s, target_component);
 }
 
 

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.