]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
multigrid with grouped components and MGTransferSelect
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 8 Jul 2003 21:40:00 +0000 (21:40 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 8 Jul 2003 21:40:00 +0000 (21:40 +0000)
git-svn-id: https://svn.dealii.org/trunk@7861 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_renumbering.h
deal.II/deal.II/include/multigrid/mg_smoother.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/multigrid/mg_dof_tools.cc
deal.II/deal.II/source/multigrid/mg_transfer_block.cc

index d3fa5741a20da24fd2e3909faa171a7f3028e5d8..dbff2d8893470266d199f16d5ef0ba0442517d1c 100644 (file)
@@ -317,7 +317,8 @@ class DoFRenumbering
     template <int dim>
     static void
     component_wise (DoFHandler<dim>                 &dof_handler,
-                   const std::vector<unsigned int> &target_component = std::vector<unsigned int>());
+                   const std::vector<unsigned int> &target_component
+                   = std::vector<unsigned int>());
 
                                     /**
                                      * Sort the degrees of freedom by
index b6899fa323c8408bd9c1581cc38fb129de66d7da..5a1583fc287bfbd37f15c7cffd23f871f4ca8bbd 100644 (file)
@@ -532,18 +532,12 @@ MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::smooth(
   if (symmetric && (steps2 % 2 == 0))
     T = false;
 //  cerr << 'S' << level;
+//  cerr << '(' << matrices[level]->m() << ',' << matrices[level]->n() << ')';
   
   for (unsigned int i=0; i<steps2; ++i)
     {
       if (T)
        {
-                                          // This is not really the
-                                          // transposed smoother, but
-                                          // just Gauss-Seidel with
-                                          // reverse numbering.  For
-                                          // a symmetric matrix, it
-                                          // is the transpose,
-                                          // though.
 //       cerr << 'T';
          matrices[level].vmult(*r,u);
          r->sadd(-1.,1.,rhs);
index 5add525f41cbc5c8bd3bc36a14b5aa52073dbbaf..3d405e64c36bc8ec6adf3c1122d8140631887597 100644 (file)
@@ -270,22 +270,27 @@ class MGTransferBlockBase
                                      * default, all matrices are
                                      * built.
                                      *
-                                     * The last argument
-                                     * @p{target_component} allows
-                                     * grouping of components the
-                                     * same way as in
-                                     * @p{DoFRenumbering::component_wise}.
+                                     * This function is only called
+                                     * by derived classes. These can
+                                     * also set the member variable
+                                     * @p{target_component} for
+                                     * re-ordering and grouping of
+                                     * components.
                                      */
     template <int dim>
     void build_matrices (const MGDoFHandler<dim>& mg_dof,
-                        const std::vector<bool>& selected,
-                        const std::vector<unsigned int>& target_component);
+                        const std::vector<bool>& selected);
 
                                   /**
                                    * Flag of selected components.
                                    */
     std::vector<bool> selected;
 
+                                    /**
+                                     * Target component if renumbering is required.
+                                     */
+    std::vector<unsigned int> target_component;
+    
                                   /**
                                    * Sizes of the multi-level vectors.
                                    */
@@ -372,13 +377,21 @@ class MGTransferBlock : public MGTransferBase<BlockVector<number> >,
                                      * default, all matrices are
                                      * built.
                                      *
+                                     * The last argument
+                                     * @p{target_component} allows
+                                     * grouping of components the
+                                     * same way as in
+                                     * @p{DoFRenumbering::component_wise}.
+                                     * 
                                      * This function is a front-end
                                      * for the same function in
                                      * @ref{MGTransferBlockBase}.
                                      */
     template <int dim>
     void build_matrices (const MGDoFHandler<dim> &mg_dof,
-                        std::vector<bool> selected = std::vector<bool>());
+                        std::vector<bool> selected = std::vector<bool>(),
+                        const std::vector<unsigned int>& target_component
+                        =std::vector<unsigned int>());
 
                                     /**
                                      * Prolongate a vector from level
@@ -541,25 +554,15 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
                                      */
     virtual ~MGTransferSelect ();
     
-                                    /**
-                                     * Actually build the prolongation
-                                     * matrices for each level.
-                                     *
-                                     * Select the component you want
-                                     * to apply multigrid to.
-                                     *
-                                     * 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);
-
                                     /**
                                      * Actually build the prolongation
                                      * matrices for grouped components.
                                      *
+                                     * @p{selected} is the number of
+                                     * the component for which the
+                                     * transfer matrices should be
+                                     * built.
+                                     *
                                      * The argument
                                      * @p{target_component}
                                      * corresponds to the grouping
@@ -577,7 +580,8 @@ class MGTransferSelect : public MGTransferBase<Vector<number> >,
     template <int dim>
     void build_matrices (const MGDoFHandler<dim> &mg_dof,
                         unsigned int selected,
-                        const std::vector<unsigned int>& target_component);
+                        const std::vector<unsigned int>& target_component
+                        = std::vector<unsigned int>());
 
                                     /**
                                      * Prolongate a vector from level
index eee2a98470e30a87f0d377ea4257e22f4b3a522c..828e0865a0f4f7f4c171eb445b120ed869f901b1 100644 (file)
@@ -365,7 +365,7 @@ copy_to_mg (const MGDoFHandler<dim>        &mg_dof_handler,
                                           // into the level-wise one
          for (unsigned int i=0; i<dofs_per_cell; ++i)
            {
-             if (fe.system_to_component_index(i).first == selected)
+             if (target_component[fe.system_to_component_index(i).first] == selected)
                dst[level](level_dof_indices[i] - level_start)
                  = src(global_dof_indices[i] - start);
            }
@@ -450,7 +450,7 @@ copy_from_mg (const MGDoFHandler<dim>              &mg_dof_handler,
                                       // copy level-wise data to
                                       // global vector
       for (unsigned int i=0; i<dofs_per_cell; ++i)
-       if (fe.system_to_component_index(i).first == selected)
+       if (target_component[fe.system_to_component_index(i).first] == selected)
          dst(global_dof_indices[i]-start)
            = src[level](level_dof_indices[i]-level_start);
     };
@@ -512,7 +512,7 @@ copy_from_mg_add (const MGDoFHandler<dim>              &mg_dof_handler,
                                       // copy level-wise data to
                                       // global vector
       for (unsigned int i=0; i<dofs_per_cell; ++i)
-       if (fe.system_to_component_index(i).first == selected)
+       if (target_component[fe.system_to_component_index(i).first] == selected)
          dst(global_dof_indices[i]-start)
            += src[level](level_dof_indices[i]-level_start);
     };
index e692d8b4ff06948373ebbe1237cc585c3f048250..778fc5f7af267b9f494d5238291fe855db306965 100644 (file)
@@ -320,7 +320,8 @@ MGTools::count_dofs_per_component (
   for (unsigned int l=0;l<nlevels;++l)
     {
       result[l].resize (n_components);
-      
+      fill (result[l].begin(),result[l].end(), 0U);
+  
                                       // special case for only one
                                       // component. treat this first
                                       // since it does not require any
@@ -356,7 +357,7 @@ MGTools::count_dofs_per_component (
                                           // next count what we got
          for (unsigned int i=0; i<n_components; ++i)
              result[l][target_component[i]]
-               = std::count(dofs_in_component[i].begin(),
+               += std::count(dofs_in_component[i].begin(),
                             dofs_in_component[i].end(),
                             true);
          
index b9a1f0b3f81c506e006f6ed4e854fc9357b7e40a..9c41d0cd1933211307ee9de09ae9b0a570cceb5a 100644 (file)
 template <int dim>
 void MGTransferBlockBase::build_matrices (
   const MGDoFHandler<dim>& mg_dof,
-  const std::vector<bool>& select,
-  const std::vector<unsigned int>& /*target_component*/)
+  const std::vector<bool>& select)
 {
+  if (target_component.size() == 0)
+    {
+      target_component.resize(mg_dof.get_fe().n_components());
+      for (unsigned int i=0;i<target_component.size();++i)
+       target_component[i] = i;
+    } else {
+      Assert (target_component.size() == mg_dof.get_fe().n_components(),
+             ExcDimensionMismatch(target_component.size(),
+                                  mg_dof.get_fe().n_components()));
+//      unsigned int previous = 0;
+      for (unsigned int i=0;i<target_component.size();++i)
+       {
+         Assert(i<target_component.size(),
+                ExcIndexRange(i,0,target_component.size()));
+//       Assert(previous<=i,ExcIndexRange(i,previous,previous+1));
+//       Assert(i<=previous+1,ExcIndexRange(i,previous,previous+1));
+       }
+    }
+  
   const FiniteElement<dim>& fe = mg_dof.get_fe();
   const unsigned int n_components  = fe.n_components();
   const unsigned int dofs_per_cell = fe.dofs_per_cell;  
@@ -44,7 +62,7 @@ void MGTransferBlockBase::build_matrices (
          ExcDimensionMismatch(selected.size(), n_components))
 
   sizes.resize(n_levels);
-  MGTools::count_dofs_per_component(mg_dof, sizes);
+  MGTools::count_dofs_per_component(mg_dof, sizes, target_component);
   
                                   // reset the size of the array of
                                   // matrices. call resize(0) first,
@@ -158,9 +176,11 @@ void MGTransferBlockBase::build_matrices (
                  for (unsigned int j=0; j<dofs_per_cell; ++j)
                    if (prolongation(i,j) != 0)
                      {
-                       const unsigned int icomp = fe.system_to_component_index(i).first;
-                       const unsigned int jcomp = fe.system_to_component_index(j).first;
-                       if ((icomp==jcomp) && selected[icomp])
+                       const unsigned int icomp
+                         = fe.system_to_component_index(i).first;
+                       const unsigned int jcomp
+                         = fe.system_to_component_index(j).first;
+                       if ((icomp==jcomp) && selected[target_component[icomp]])
                          prolongation_sparsities[level]->add(dof_indices_child[i],
                                                              dof_indices_parent[j]);
                      };
@@ -196,7 +216,7 @@ void MGTransferBlockBase::build_matrices (
                      {
                        const unsigned int icomp = fe.system_to_component_index(i).first;
                        const unsigned int jcomp = fe.system_to_component_index(j).first;
-                       if ((icomp==jcomp) && selected[icomp])
+                       if ((icomp==jcomp) && selected[target_component[icomp]])
                          prolongation_matrices[level]->set(dof_indices_child[i],
                                                            dof_indices_parent[j],
                                                            prolongation(i,j));
@@ -249,15 +269,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,
+  const std::vector<unsigned int>& t_component)
 {
   if (select.size() == 0)
     select = std::vector<bool> (mg_dof.get_fe().n_components(), 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, select, target_component);
+
+  target_component = t_component;
+
+  MGTransferBlockBase::build_matrices (mg_dof, select);
 }
 
 
@@ -265,16 +285,24 @@ template <typename number>
 template <int dim>
 void MGTransferSelect<number>::build_matrices (
   const MGDoFHandler<dim> &mg_dof,
-  unsigned int select)
+  unsigned int select,
+  const std::vector<unsigned int>& t_component)
 {
+  unsigned int ncomp = mg_dof.get_fe().n_components();
+  
   selected = select;
-  std::vector<bool> s(mg_dof.get_fe().n_components(), false);
+  std::vector<bool> s(ncomp, 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, target_component);
+  target_component = t_component;
+  if (target_component.size() == 0)
+    {
+      target_component.resize(ncomp);
+      for (unsigned int i=0;i<ncomp;++i)
+       target_component[i] = i;
+    }
+  
+  MGTransferBlockBase::build_matrices (mg_dof, s);
 }
 
 
@@ -284,22 +312,22 @@ void MGTransferSelect<number>::build_matrices (
 template
 void MGTransferBlock<float>::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- std::vector<bool>);
+ std::vector<bool>, const std::vector<unsigned int>&);
 
 template
 void MGTransferSelect<float>::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- unsigned int);
+ unsigned int, const std::vector<unsigned int>&);
 
 template
 void MGTransferBlock<double>::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- std::vector<bool>);
+ std::vector<bool>, const std::vector<unsigned int>&);
 
 template
 void MGTransferSelect<double>::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- unsigned int);
+ unsigned int, const std::vector<unsigned int>&);
 
 
 template void

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.