]> https://gitweb.dealii.org/ - dealii.git/commitdiff
changed the reinit_vector function
authorBaerbel Jannsen <baerbel.janssen@gmail.com>
Wed, 6 Jan 2010 08:13:21 +0000 (08:13 +0000)
committerBaerbel Jannsen <baerbel.janssen@gmail.com>
Wed, 6 Jan 2010 08:13:21 +0000 (08:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@20299 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_base.h
deal.II/deal.II/include/multigrid/mg_smoother.h
deal.II/deal.II/source/multigrid/mg_base.cc

index 0c8a42d5e9ff410b0109e5bf953e028334a625ac..8bea78bf4ec5543b1e754b0cb4d3a631c9923174 100644 (file)
@@ -42,10 +42,15 @@ namespace internal
                                      * of freedom on each level and
                                      * @p reinit each level vector
                                      * to this length.
+                                      * For compatibility reasons with 
+                                      * the next function 
+                                      * the target_component is added 
+                                      * here but is not used.
                                      */
     template <int dim, typename number, int spacedim>
     static void
     reinit_vector (const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+                   std::vector<unsigned int> target_component,
                   MGLevelObject<dealii::Vector<number> > &vector);
 
                                     /**
@@ -54,11 +59,14 @@ namespace internal
                                      * count the numbers of degrees
                                      * of freedom on each level and
                                      * @p reinit each level vector
-                                     * to this length.
+                                     * to this length. The target_component 
+                                      * is handed to MGTools::count_dofs_per_block.
+                                      * See for documentation there.
                                      */
     template <int dim, typename number, int spacedim>
     static void
     reinit_vector (const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+                   std::vector<unsigned int> target_component,
                   MGLevelObject<BlockVector<number> > &vector);
 
 
index fdcb4b729f3ce1c44cafe68cc729dd7aa0ee25ca..ac6b7078c58bf25c9550f9af62ad45fcf5993d7f 100644 (file)
@@ -641,6 +641,11 @@ MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::initialize (
 {
   const unsigned int min = m.get_minlevel();
   const unsigned int max = m.get_maxlevel();
+
+  Assert (data.get_minlevel() == min, 
+      ExcDimensionMismatch(data.get_minlevel(), min));
+  Assert (data.get_maxlevel() == max, 
+      ExcDimensionMismatch(data.get_maxlevel(), max));
   
   matrices.resize(min, max);
   smoothers.resize(min, max);
@@ -686,6 +691,11 @@ MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::initialize (
   const unsigned int min = m.get_minlevel();
   const unsigned int max = m.get_maxlevel();
   
+  Assert (data.get_minlevel() == min, 
+      ExcDimensionMismatch(data.get_minlevel(), min));
+  Assert (data.get_maxlevel() == max, 
+      ExcDimensionMismatch(data.get_maxlevel(), max));
+
   matrices.resize(min, max);
   smoothers.resize(min, max);
 
@@ -862,6 +872,11 @@ MGSmootherPrecondition<MATRIX, RELAX, VECTOR>::initialize (
   const unsigned int min = m.get_minlevel();
   const unsigned int max = m.get_maxlevel();
   
+  Assert (data.get_minlevel() == min, 
+      ExcDimensionMismatch(data.get_minlevel(), min));
+  Assert (data.get_maxlevel() == max, 
+      ExcDimensionMismatch(data.get_maxlevel(), max));
+
   matrices.resize(min, max);
   smoothers.resize(min, max);
 
@@ -906,6 +921,11 @@ MGSmootherPrecondition<MATRIX, RELAX, VECTOR>::initialize (
   const unsigned int min = m.get_minlevel();
   const unsigned int max = m.get_maxlevel();
   
+  Assert (data.get_minlevel() == min, 
+      ExcDimensionMismatch(data.get_minlevel(), min));
+  Assert (data.get_maxlevel() == max, 
+      ExcDimensionMismatch(data.get_maxlevel(), max));
+
   matrices.resize(min, max);
   smoothers.resize(min, max);
 
index ab2cc9bec19ffd30a8cb1b0477313c04d9ffbb42..77e6ea3a0e2070314329320346754672b2a10b54 100644 (file)
@@ -32,6 +32,7 @@ namespace internal
     template <int dim, typename number, int spacedim>
     void
     reinit_vector (const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+                   std::vector<unsigned int> ,
                   MGLevelObject<dealii::Vector<number> > &v)
     {
       for (unsigned int level=v.get_minlevel();
@@ -47,18 +48,35 @@ namespace internal
     template <int dim, typename number, int spacedim>
     void
     reinit_vector (const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+                   std::vector<unsigned int> target_component,
                   MGLevelObject<BlockVector<number> > &v)
     {
       const unsigned int n_blocks = mg_dof.get_fe().n_blocks();
+      if (target_component.size()==0)
+      {
+        target_component.resize(n_blocks);
+        for (unsigned int i=0;i<n_blocks;++i)
+          target_component[i] = i;
+      }
+      Assert(target_component.size()==n_blocks,
+          ExcDimensionMismatch(target_component.size(),n_blocks));
+      const unsigned int max_block
+        = *std::max_element (target_component.begin(),
+            target_component.end());
+      const unsigned int n_target_blocks = max_block + 1;
+
       std::vector<std::vector<unsigned int> >
        ndofs(mg_dof.get_tria().n_levels(),
              std::vector<unsigned int>(n_blocks));
-      MGTools::count_dofs_per_block (mg_dof, ndofs);
+      MGTools::count_dofs_per_block (mg_dof, ndofs, target_component);
 
       for (unsigned int level=v.get_minlevel();
           level<=v.get_maxlevel();++level)
        {
-         v[level].reinit(ndofs[level]);
+         v[level].reinit(n_target_blocks);
+          for (unsigned int b=0; b<n_target_blocks; ++b)
+            v[level].block(b).reinit(ndofs[level][b]);
+          v[level].collect_sizes();
        }
     }
 
@@ -285,16 +303,20 @@ namespace internal
   {
     template void reinit_vector<deal_II_dimension> (
       const dealii::MGDoFHandler<deal_II_dimension>&,
+      std::vector<unsigned int>,
       MGLevelObject<dealii::Vector<double> >&);
     template void reinit_vector<deal_II_dimension> (
       const dealii::MGDoFHandler<deal_II_dimension>&,
+      std::vector<unsigned int>,
       MGLevelObject<dealii::Vector<float> >&);
 
     template void reinit_vector<deal_II_dimension> (
       const dealii::MGDoFHandler<deal_II_dimension>&,
+      std::vector<unsigned int>,
       MGLevelObject<BlockVector<double> >&);
     template void reinit_vector<deal_II_dimension> (
       const dealii::MGDoFHandler<deal_II_dimension>&,
+      std::vector<unsigned int>,
       MGLevelObject<BlockVector<float> >&);
 
     template void reinit_vector_by_components<deal_II_dimension> (

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.