]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Move functions into .cc files that are only used in these very files.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 6 Jan 2010 09:32:48 +0000 (09:32 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 6 Jan 2010 09:32:48 +0000 (09:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@20302 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_base.h
deal.II/deal.II/include/multigrid/mg_transfer.templates.h
deal.II/deal.II/include/multigrid/mg_transfer_block.templates.h
deal.II/deal.II/include/multigrid/mg_transfer_component.templates.h
deal.II/deal.II/source/multigrid/mg_base.cc
deal.II/deal.II/source/multigrid/mg_transfer_block.cc
deal.II/deal.II/source/multigrid/mg_transfer_component.cc
deal.II/deal.II/source/multigrid/mg_transfer_prebuilt.cc

index 8bea78bf4ec5543b1e754b0cb4d3a631c9923174..de6ea5961202397d619a902cfd8aa55cead59a77 100644 (file)
@@ -31,154 +31,6 @@ template <int, int> class MGDoFHandler;
 template <typename> class MGLevelObject;
 
 
-namespace internal
-{
-  namespace mg
-  {
-                                    /**
-                                     * Ajust vectors on all levels to
-                                     * correct size.  Here, we just
-                                     * count the numbers of degrees
-                                     * 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);
-
-                                    /**
-                                     * Ajust vectors on all levels to
-                                     * correct size.  Here, we just
-                                     * count the numbers of degrees
-                                     * of freedom on each level and
-                                     * @p reinit each level vector
-                                     * 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);
-
-
-                                    /**
-                                     * Adjust vectors on all levels
-                                     * to correct size. The degrees
-                                     * of freedom on each level are
-                                     * counted by block and only the
-                                     * block selected is used.
-                                     */
-    template <int dim, typename number, int spacedim>
-    static void
-    reinit_vector_by_blocks (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<dealii::Vector<number> > &v,
-      const unsigned int selected,
-      std::vector<std::vector<unsigned int> >& cached_sizes);
-
-                                    /**
-                                     * Adjust block vectors on all
-                                     * levels to correct size. The
-                                     * degrees of freedom on each
-                                     * level are counted by block.
-                                     */
-    template <int dim, typename number, int spacedim>
-    static void
-    reinit_vector_by_blocks (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<BlockVector<number> > &v,
-      const std::vector<bool> &selected,
-      std::vector<std::vector<unsigned int> >& cached_sizes);
-
-                                    /**
-                                     * Adjust block-vectors on all
-                                     * levels to correct size.  Count
-                                     * the numbers of degrees of
-                                     * freedom on each level
-                                     * component-wise. Then, assign
-                                     * each block of @p vector the
-                                     * corresponding size.
-                                     *
-                                     * The boolean field @p selected
-                                     * allows restricting this
-                                     * operation to certain
-                                     * components. In this case, @p
-                                     * vector will only have as many
-                                     * blocks as there are true
-                                     * values in @p selected (no
-                                     * blocks of length zero are
-                                     * padded in). If this argument
-                                     * is omitted, all blocks will be
-                                     * considered.
-                                     *
-                                     * Degrees of freedom must be
-                                     * sorted by component in order
-                                     * to obtain reasonable results
-                                     * from this function.
-                                     *
-                                     * The argument
-                                     * @p target_component allows to
-                                     * re-sort and group components
-                                     * as in
-                                     * DoFRenumbering::component_wise.
-                                     *
-                                     *
-                                     */
-    template <int dim, typename number, int spacedim>
-    static void
-    reinit_vector_by_components (const dealii::MGDoFHandler<dim,spacedim>& mg_dof,
-                                MGLevelObject<BlockVector<number> >& v,
-                                const std::vector<bool>& selected,
-                                const std::vector<unsigned int>& target_component,
-                                std::vector<std::vector<unsigned int> >& cached_sizes);
-                                    /**
-                                     * Adjust vectors on all levels
-                                     * to correct size.  Count the
-                                     * numbers of degrees of freedom
-                                     * on each level component-wise
-                                     * in a single component. Then,
-                                     * assign @p vector the
-                                     * corresponding size.
-                                     *
-                                     * The boolean field @p selected
-                                     * may be nonzero in a single
-                                     * component, indicating the
-                                     * block of a block vector the
-                                     * argument @p v corresponds to.
-                                     *
-                                     * Degrees of freedom must be
-                                     * sorted by component in order
-                                     * to obtain reasonable results
-                                     * from this function.
-                                     *
-                                     * The argument
-                                     * @p target_component allows to
-                                     * re-sort and group components
-                                     * as in
-                                     * DoFRenumbering::component_wise.
-                                     */
-    template <int dim, typename number, int spacedim>
-    static void
-    reinit_vector_by_components (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<dealii::Vector<number> > &v,
-      const std::vector<bool> &selected,
-      const std::vector<unsigned int> &target_component,
-      std::vector<std::vector<unsigned int> >& cached_sizes);
-
-  }
-}
-
-
 
 /*!@addtogroup mg */
 /*@{*/
index a756eff5fb6596df5306f1ea79e8d4f1f7cde75a..6ddff575484a1f732d3641e4ea26d0c3f62b20f8 100644 (file)
@@ -29,37 +29,6 @@ DEAL_II_NAMESPACE_OPEN
 /* --------------------- MGTransferPrebuilt -------------- */
 
 
-template <class VECTOR>
-template <int dim, class InVector, int spacedim>
-void
-MGTransferPrebuilt<VECTOR>::copy_to_mg (
-  const MGDoFHandler<dim,spacedim>& mg_dof_handler,
-  MGLevelObject<VECTOR>& dst,
-  const InVector& src) const
-{
-  internal::mg::reinit_vector(mg_dof_handler, component_to_block_map, dst);
-  bool first = true;
-  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
-    {
-      --level;
-
-      typedef std::map<unsigned int, unsigned int>::const_iterator IT;
-      for (IT i= copy_indices[level].begin();
-          i != copy_indices[level].end();++i)
-       dst[level](i->second) = src(i->first);
-                                      // For non-DG: degrees of
-                                      // freedom in the refinement
-                                      // face may need special
-                                      // attention, since they belong
-                                      // to the coarse level, but
-                                      // have fine level basis
-                                      // functions
-      if (!first)
-       restrict_and_add (level+1, dst[level], dst[level+1]);
-      first = false;
-    }
-}
-
 
 
 template <class VECTOR>
index f9a903efac55c9eb7862d5dbc81e1a4d843d38c6..b80d4d652803dcf986c3641180a844c229311df5 100644 (file)
@@ -34,63 +34,6 @@ typedef std::map<unsigned int, unsigned int>::const_iterator IT;
 
 
 
-template <typename number>
-template <int dim, typename number2, int spacedim>
-void
-MGTransferBlockSelect<number>::copy_to_mg (
-  const MGDoFHandler<dim,spacedim>        &mg_dof_handler,
-  MGLevelObject<Vector<number> > &dst,
-  const BlockVector<number2>     &src) const
-{
-  internal::mg::reinit_vector_by_blocks(mg_dof_handler, dst,
-                                       selected_block, sizes);
-                                  // For MGTransferBlockSelect, the
-                                  // multilevel block is always the
-                                  // first, since only one block is
-                                  // selected.
-  bool first = true;
-  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
-    {
-      --level;
-      for (IT i= copy_indices[selected_block][level].begin();
-          i != copy_indices[selected_block][level].end();++i)
-       dst[level](i->second) = src.block(selected_block)(i->first);
-      if (!first)
-       restrict_and_add (level+1, dst[level], dst[level+1]);
-      first = false;
-    }
-}
-
-
-
-template <typename number>
-template <int dim, typename number2, int spacedim>
-void
-MGTransferBlockSelect<number>::copy_to_mg (
-  const MGDoFHandler<dim,spacedim>        &mg_dof_handler,
-  MGLevelObject<Vector<number> > &dst,
-  const Vector<number2>          &src) const
-{
-  internal::mg::reinit_vector_by_blocks(mg_dof_handler, dst,
-                                       selected_block, sizes);
-                                  // For MGTransferBlockSelect, the
-                                  // multilevel block is always the
-                                  // first, since only one block is selected.
-  bool first = true;
-  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
-    {
-      --level;
-      for (IT i= copy_indices[selected_block][level].begin();
-          i != copy_indices[selected_block][level].end();++i)
-       dst[level](i->second) = src(i->first);
-      if (!first)
-       restrict_and_add (level+1, dst[level], dst[level+1]);
-      first = false;
-    }
-}
-
-
-
 template <typename number>
 template <int dim, typename number2, int spacedim>
 void
@@ -167,34 +110,6 @@ MGTransferBlockSelect<number>::memory_consumption () const
 
 
 
-template <typename number>
-template <int dim, typename number2, int spacedim>
-void
-MGTransferBlock<number>::copy_to_mg (
-  const MGDoFHandler<dim,spacedim>& mg_dof_handler,
-  MGLevelObject<BlockVector<number> >& dst,
-  const BlockVector<number2>& src) const
-{
-  internal::mg::reinit_vector_by_blocks(mg_dof_handler, dst,
-                                       selected, sizes);
-  bool first = true;
-  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
-    {
-      --level;
-      for (unsigned int block=0;block<selected.size();++block)
-       if (selected[block])
-         for (IT i= copy_indices[block][level].begin();
-              i != copy_indices[block][level].end();++i)
-           dst[level].block(mg_block[block])(i->second) = src.block(block)(i->first);
-      if (!first)
-       restrict_and_add (level+1, dst[level], dst[level+1]);
-      first = false;
-    }
-}
-
-
-
-
 template <typename number>
 template <int dim, typename number2, int spacedim>
 void
index 65a2b9c0f21fca1f9502cf922e53c287de99197a..4765f323c9df1d977d9bc0a963ef363d502a4c49 100644 (file)
@@ -57,106 +57,6 @@ MGTransferSelect<number>::copy_to_mg (
 
 
 
-template <typename number>
-template <int dim, class InVector, int spacedim>
-void
-MGTransferSelect<number>::do_copy_to_mg (
-  const MGDoFHandler<dim,spacedim>&        mg_dof_handler,
-  MGLevelObject<Vector<number> >& dst,
-  const InVector&                 src,
-  const unsigned int              offset) const
-{
-  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
-
-  const unsigned int dofs_per_cell = fe.dofs_per_cell;
-
-                                  // set the elements of the vectors
-                                  // on all levels to zero
-  unsigned int minlevel = dst.get_minlevel();
-  unsigned int maxlevel = dst.get_maxlevel();
-
-  dst=0;
-
-  Assert(sizes.size()==mg_dof_handler.get_tria().n_levels(),
-        ExcMatricesNotBuilt());
-
-  internal::mg::reinit_vector_by_components(
-    mg_dof_handler, dst, mg_selected, mg_target_component, sizes);
-
-  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
-  std::vector<unsigned int> level_dof_indices  (dofs_per_cell);
-
-                                  // Build a vector of the selected
-                                  // indices, since traversing all
-                                  // indices on each cell is too
-                                  // slow.
-  std::vector<unsigned int> selected_indices;
-  selected_indices.reserve(dofs_per_cell);
-  for (unsigned int i=0; i<dofs_per_cell; ++i)
-    if (mg_selected[mg_target_component[fe.system_to_component_index(i).first]])
-      selected_indices.push_back(i);
-  unsigned int selected_component
-    = mg_target_component[fe.system_to_component_index(selected_indices[0]).first];
-
-                                  // traverse the grid top-down
-                                  // (i.e. starting with the most
-                                  // refined grid). this way, we can
-                                  // always get that part of one
-                                  // level of the output vector which
-                                  // corresponds to a region which is
-                                  // more refined, by restriction of
-                                  // the respective vector on the
-                                  // next finer level, which we then
-                                  // already have built.
-  for (int level=maxlevel; level>=static_cast<signed int>(minlevel); --level)
-    {
-
-  typename MGDoFHandler<dim,spacedim>::active_cell_iterator
-       level_cell = mg_dof_handler.begin_active(level);
-      const typename MGDoFHandler<dim,spacedim>::active_cell_iterator
-       level_end  = mg_dof_handler.end_active(level);
-                                      // Compute coarse level right hand side
-                                      // by restricting from fine level.
-      for (; level_cell!=level_end; ++level_cell)
-       {
-                                          // get the dof numbers of
-                                          // this cell for the global
-                                          // and the level-wise
-                                          // numbering
-         level_cell->get_dof_indices(global_dof_indices);
-         level_cell->get_mg_dof_indices (level_dof_indices);
-
-                                          // transfer the global
-                                          // defect in the vector
-                                          // into the level-wise one
-         const unsigned int level_start
-           = mg_component_start[level][selected_component];
-         const typename std::vector<unsigned int>::const_iterator
-           end = selected_indices.end();
-
-         for (typename std::vector<unsigned int>::const_iterator
-                i=selected_indices.begin();
-              i != end; ++i)
-           {
-             dst[level](level_dof_indices[*i] - level_start)
-               = src(global_dof_indices[*i] - offset);
-           }
-       }
-
-                                      // for that part of the level
-                                      // which is further refined:
-                                      // get the defect by
-                                      // restriction of the defect on
-                                      // one level higher
-      if (static_cast<unsigned int>(level) < maxlevel)
-       {
-         restrict_and_add (level+1, dst[level], dst[level+1]);
-       }
-    };
-}
-
-
-
 template <typename number>
 template <int dim, typename number2, int spacedim>
 void
index 77e6ea3a0e2070314329320346754672b2a10b54..90476e759b2a0630c26fcc916c8411dcf805eba7 100644 (file)
 #include <multigrid/mg_level_object.h>
 #include <multigrid/mg_tools.h>
 
-#include <numeric>
 
 DEAL_II_NAMESPACE_OPEN
 
 
-namespace internal
-{
-  namespace mg
-  {
-    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();
-          level<=v.get_maxlevel();++level)
-       {
-         unsigned int n = mg_dof.n_dofs (level);
-         v[level].reinit(n);
-       }
-
-    }
-
-
-    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, target_component);
-
-      for (unsigned int level=v.get_minlevel();
-          level<=v.get_maxlevel();++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();
-       }
-    }
-
-
-    template <int dim, typename number, int spacedim>
-    void
-    reinit_vector_by_components (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<BlockVector<number> > &v,
-      const std::vector<bool> &sel,
-      const std::vector<unsigned int> &target_comp,
-      std::vector<std::vector<unsigned int> >& ndofs)
-    {
-      std::vector<bool> selected=sel;
-      std::vector<unsigned int> target_component=target_comp;
-      const unsigned int ncomp = mg_dof.get_fe().n_components();
-
-                                      // If the selected and
-                                      // target_component have size 0,
-                                      // they must be replaced by default
-                                      // values.
-                                      //
-                                      // Since we already made copies
-                                      // directly after this function was
-                                      // called, we use the arguments
-                                      // directly.
-      if (target_component.size() == 0)
-       {
-         target_component.resize(ncomp);
-         for (unsigned int i=0;i<ncomp;++i)
-           target_component[i] = i;
-       }
-
-                                      // If selected is an empty vector,
-                                      // all components are selected.
-      if (selected.size() == 0)
-       {
-         selected.resize(target_component.size());
-         std::fill_n (selected.begin(), ncomp, false);
-         for (unsigned int i=0;i<target_component.size();++i)
-           selected[target_component[i]] = true;
-       }
-
-      Assert (selected.size() == target_component.size(),
-             ExcDimensionMismatch(selected.size(), target_component.size()));
-
-                                      // Compute the number of blocks needed
-      const unsigned int n_selected
-       = std::accumulate(selected.begin(),
-                         selected.end(),
-                         0U);
-
-      if (ndofs.size() == 0)
-       {
-         std::vector<std::vector<unsigned int> >
-           new_dofs(mg_dof.get_tria().n_levels(),
-                    std::vector<unsigned int>(target_component.size()));
-         std::swap(ndofs, new_dofs);
-         MGTools::count_dofs_per_component (mg_dof, ndofs,
-                                            true, target_component);
-       }
-
-      for (unsigned int level=v.get_minlevel();
-          level<=v.get_maxlevel();++level)
-       {
-         v[level].reinit(n_selected, 0);
-         unsigned int k=0;
-         for (unsigned int i=0;i<selected.size() && (k<v[level].n_blocks());++i)
-           {
-             if (selected[i])
-               {
-                 v[level].block(k++).reinit(ndofs[level][i]);
-               }
-             v[level].collect_sizes();
-           }
-       }
-    }
-
-
-    template <int dim, typename number, int spacedim>
-    void
-    reinit_vector_by_components (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<dealii::Vector<number> > &v,
-      const std::vector<bool> &selected,
-      const std::vector<unsigned int> &target_component,
-      std::vector<std::vector<unsigned int> >& ndofs)
-    {
-      Assert (selected.size() == target_component.size(),
-             ExcDimensionMismatch(selected.size(), target_component.size()));
-
-                                      // Compute the number of blocks needed
-#ifdef DEBUG
-      const unsigned int n_selected
-       = std::accumulate(selected.begin(),
-                         selected.end(),
-                         0U);
-      Assert(n_selected == 1, ExcDimensionMismatch(n_selected, 1));
-#endif
-
-      unsigned int selected_block = 0;
-      while (!selected[selected_block])
-       ++selected_block;
-
-      if (ndofs.size() == 0)
-       {
-         std::vector<std::vector<unsigned int> >
-           new_dofs(mg_dof.get_tria().n_levels(),
-                    std::vector<unsigned int>(target_component.size()));
-         std::swap(ndofs, new_dofs);
-         MGTools::count_dofs_per_component (mg_dof, ndofs,
-                                            true, target_component);
-       }
-
-      for (unsigned int level=v.get_minlevel();
-          level<=v.get_maxlevel();++level)
-       {
-         v[level].reinit(ndofs[level][selected_block]);
-       }
-    }
-
-
-    template <int dim, typename number, int spacedim>
-    void
-    reinit_vector_by_blocks (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<BlockVector<number> > &v,
-      const std::vector<bool> &sel,
-      std::vector<std::vector<unsigned int> >& ndofs)
-    {
-      std::vector<bool> selected=sel;
-                                      // Compute the number of blocks needed
-      const unsigned int n_selected
-       = std::accumulate(selected.begin(),
-                         selected.end(),
-                         0U);
-
-      if (ndofs.size() == 0)
-       {
-         std::vector<std::vector<unsigned int> >
-           new_dofs(mg_dof.get_tria().n_levels(),
-                    std::vector<unsigned int>(selected.size()));
-         std::swap(ndofs, new_dofs);
-         MGTools::count_dofs_per_block (mg_dof, ndofs);
-       }
-
-      for (unsigned int level=v.get_minlevel();
-          level<=v.get_maxlevel();++level)
-       {
-         v[level].reinit(n_selected, 0);
-         unsigned int k=0;
-         for (unsigned int i=0;i<selected.size() && (k<v[level].n_blocks());++i)
-           {
-             if (selected[i])
-               {
-                 v[level].block(k++).reinit(ndofs[level][i]);
-               }
-             v[level].collect_sizes();
-           }
-       }
-    }
-
-
-    template <int dim, typename number, int spacedim>
-    void
-    reinit_vector_by_blocks (
-      const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
-      MGLevelObject<dealii::Vector<number> > &v,
-      const unsigned int selected_block,
-      std::vector<std::vector<unsigned int> >& ndofs)
-    {
-      const unsigned int n_blocks = mg_dof.get_fe().n_blocks();
-      Assert(selected_block < n_blocks, ExcIndexRange(selected_block, 0, n_blocks));
-
-      std::vector<bool> selected(n_blocks, false);
-      selected[selected_block] = true;
-
-      if (ndofs.size() == 0)
-       {
-         std::vector<std::vector<unsigned int> >
-           new_dofs(mg_dof.get_tria().n_levels(),
-                    std::vector<unsigned int>(selected.size()));
-         std::swap(ndofs, new_dofs);
-         MGTools::count_dofs_per_block (mg_dof, ndofs);
-       }
-
-      for (unsigned int level=v.get_minlevel();
-          level<=v.get_maxlevel();++level)
-       {
-         v[level].reinit(ndofs[level][selected_block]);
-       }
-    }
-  }
-}
-
-
-
 template <class VECTOR>
 MGTransferBase<VECTOR>::~MGTransferBase()
 {}
@@ -297,67 +47,6 @@ MGCoarseGridBase<VECTOR>::~MGCoarseGridBase()
 // Explicit instantiations
 
 //TODO: Use the template expander script for this
-namespace internal
-{
-  namespace mg
-  {
-    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> (
-      const dealii::MGDoFHandler<deal_II_dimension>&,
-      MGLevelObject<BlockVector<double> >&,
-      const std::vector<bool> &,
-      const std::vector<unsigned int> &,
-      std::vector<std::vector<unsigned int> >&);
-    template void reinit_vector_by_components<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&,
-      MGLevelObject<BlockVector<float> >&,
-      const std::vector<bool> &,
-      const std::vector<unsigned int> &,
-      std::vector<std::vector<unsigned int> >&);
-
-    template void reinit_vector_by_components<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&, MGLevelObject<dealii::Vector<double> >&,
-      const std::vector<bool>&, const std::vector<unsigned int>&,
-      std::vector<std::vector<unsigned int> >&);
-    template void reinit_vector_by_components<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&, MGLevelObject<dealii::Vector<float> >&,
-      const std::vector<bool>&, const std::vector<unsigned int>&,
-      std::vector<std::vector<unsigned int> >&);
-
-    template void reinit_vector_by_blocks<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&, MGLevelObject<BlockVector<double> >&,
-      const std::vector<bool> &, std::vector<std::vector<unsigned int> >&);
-    template void reinit_vector_by_blocks<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&, MGLevelObject<BlockVector<float> >&,
-      const std::vector<bool> &, std::vector<std::vector<unsigned int> >&);
-
-    template void reinit_vector_by_blocks<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&, MGLevelObject<dealii::Vector<double> >&,
-      unsigned int, std::vector<std::vector<unsigned int> >&);
-    template void reinit_vector_by_blocks<deal_II_dimension> (
-      const dealii::MGDoFHandler<deal_II_dimension>&, MGLevelObject<dealii::Vector<float> >&,
-      unsigned int, std::vector<std::vector<unsigned int> >&);
-  }
-}
-
-
 template class MGTransferBase<dealii::Vector<double> >;
 template class MGTransferBase<dealii::Vector<float> >;
 template class MGTransferBase<BlockVector<double> >;
index 4ac5aaa0d79351eaca371f9997e91ecb004bcc96..9684b3668c99a90366ab51835b476c4f4f5f06f4 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 #include <multigrid/mg_tools.h>
 
 #include <algorithm>
+#include <numeric>
 #include <iostream>
 
 DEAL_II_NAMESPACE_OPEN
 
+namespace
+{
+                                    /**
+                                     * Adjust vectors on all levels
+                                     * to correct size. The degrees
+                                     * of freedom on each level are
+                                     * counted by block and only the
+                                     * block selected is used.
+                                     */
+  template <int dim, typename number, int spacedim>
+  void
+  reinit_vector_by_blocks (
+    const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+    MGLevelObject<BlockVector<number> > &v,
+    const std::vector<bool> &sel,
+    std::vector<std::vector<unsigned int> >& ndofs)
+  {
+    std::vector<bool> selected=sel;
+                                    // Compute the number of blocks needed
+    const unsigned int n_selected
+      = std::accumulate(selected.begin(),
+                       selected.end(),
+                       0U);
+
+    if (ndofs.size() == 0)
+      {
+       std::vector<std::vector<unsigned int> >
+         new_dofs(mg_dof.get_tria().n_levels(),
+                  std::vector<unsigned int>(selected.size()));
+       std::swap(ndofs, new_dofs);
+       MGTools::count_dofs_per_block (mg_dof, ndofs);
+      }
+
+    for (unsigned int level=v.get_minlevel();
+        level<=v.get_maxlevel();++level)
+      {
+       v[level].reinit(n_selected, 0);
+       unsigned int k=0;
+       for (unsigned int i=0;i<selected.size() && (k<v[level].n_blocks());++i)
+         {
+           if (selected[i])
+             {
+               v[level].block(k++).reinit(ndofs[level][i]);
+             }
+           v[level].collect_sizes();
+         }
+      }
+  }
+
+
+                                    /**
+                                     * Adjust block vectors on all
+                                     * levels to correct size. The
+                                     * degrees of freedom on each
+                                     * level are counted by block.
+                                     */
+  template <int dim, typename number, int spacedim>
+  void
+  reinit_vector_by_blocks (
+    const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+    MGLevelObject<dealii::Vector<number> > &v,
+    const unsigned int selected_block,
+    std::vector<std::vector<unsigned int> >& ndofs)
+  {
+    const unsigned int n_blocks = mg_dof.get_fe().n_blocks();
+    Assert(selected_block < n_blocks, ExcIndexRange(selected_block, 0, n_blocks));
+
+    std::vector<bool> selected(n_blocks, false);
+    selected[selected_block] = true;
+
+    if (ndofs.size() == 0)
+      {
+       std::vector<std::vector<unsigned int> >
+         new_dofs(mg_dof.get_tria().n_levels(),
+                  std::vector<unsigned int>(selected.size()));
+       std::swap(ndofs, new_dofs);
+       MGTools::count_dofs_per_block (mg_dof, ndofs);
+      }
+
+    for (unsigned int level=v.get_minlevel();
+        level<=v.get_maxlevel();++level)
+      {
+       v[level].reinit(ndofs[level][selected_block]);
+      }
+  }
+}
+
+
+template <typename number>
+template <int dim, typename number2, int spacedim>
+void
+MGTransferBlockSelect<number>::copy_to_mg (
+  const MGDoFHandler<dim,spacedim>        &mg_dof_handler,
+  MGLevelObject<Vector<number> > &dst,
+  const BlockVector<number2>     &src) const
+{
+  reinit_vector_by_blocks(mg_dof_handler, dst, selected_block, sizes);
+                                  // For MGTransferBlockSelect, the
+                                  // multilevel block is always the
+                                  // first, since only one block is
+                                  // selected.
+  bool first = true;
+  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
+    {
+      --level;
+      for (IT i= copy_indices[selected_block][level].begin();
+          i != copy_indices[selected_block][level].end();++i)
+       dst[level](i->second) = src.block(selected_block)(i->first);
+      if (!first)
+       restrict_and_add (level+1, dst[level], dst[level+1]);
+      first = false;
+    }
+}
+
+
+
+template <typename number>
+template <int dim, typename number2, int spacedim>
+void
+MGTransferBlockSelect<number>::copy_to_mg (
+  const MGDoFHandler<dim,spacedim>        &mg_dof_handler,
+  MGLevelObject<Vector<number> > &dst,
+  const Vector<number2>          &src) const
+{
+  reinit_vector_by_blocks(mg_dof_handler, dst, selected_block, sizes);
+                                  // For MGTransferBlockSelect, the
+                                  // multilevel block is always the
+                                  // first, since only one block is selected.
+  bool first = true;
+  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
+    {
+      --level;
+      for (IT i= copy_indices[selected_block][level].begin();
+          i != copy_indices[selected_block][level].end();++i)
+       dst[level](i->second) = src(i->first);
+      if (!first)
+       restrict_and_add (level+1, dst[level], dst[level+1]);
+      first = false;
+    }
+}
+
+
+
+template <typename number>
+template <int dim, typename number2, int spacedim>
+void
+MGTransferBlock<number>::copy_to_mg (
+  const MGDoFHandler<dim,spacedim>& mg_dof_handler,
+  MGLevelObject<BlockVector<number> >& dst,
+  const BlockVector<number2>& src) const
+{
+  reinit_vector_by_blocks(mg_dof_handler, dst, selected, sizes);
+  bool first = true;
+  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
+    {
+      --level;
+      for (unsigned int block=0;block<selected.size();++block)
+       if (selected[block])
+         for (IT i= copy_indices[block][level].begin();
+              i != copy_indices[block][level].end();++i)
+           dst[level].block(mg_block[block])(i->second) = src.block(block)(i->first);
+      if (!first)
+       restrict_and_add (level+1, dst[level], dst[level+1]);
+      first = false;
+    }
+}
+
 
 
 template <int dim, int spacedim>
@@ -43,10 +211,10 @@ void MGTransferBlockBase::build_matrices (
   const unsigned int n_blocks  = fe.n_blocks();
   const unsigned int dofs_per_cell = fe.dofs_per_cell;
   const unsigned int n_levels      = mg_dof.get_tria().n_levels();
-  
+
   Assert (selected.size() == n_blocks,
          ExcDimensionMismatch(selected.size(), n_blocks));
-  
+
                                   // Compute the mapping between real
                                   // blocks and blocks used for
                                   // multigrid computations.
@@ -57,11 +225,11 @@ void MGTransferBlockBase::build_matrices (
       mg_block[i] = n_mg_blocks++;
     else
       mg_block[i] = numbers::invalid_unsigned_int;
-  
+
                                   // Compute the lengths of all blocks
   sizes.resize(n_levels);
   MGTools::count_dofs_per_block(mg_dof, sizes);
-  
+
                                   // Fill some index vectors
                                   // for later use.
   mg_block_start = sizes;
@@ -76,7 +244,7 @@ void MGTransferBlockBase::build_matrices (
          k += t;
        }
     }
-  
+
   block_start.resize(n_blocks);
   DoFTools::count_dofs_per_block (static_cast<const DoFHandler<dim,spacedim>&>(mg_dof), block_start);
 
@@ -97,9 +265,9 @@ void MGTransferBlockBase::build_matrices (
   for (unsigned int block=0;block<n_blocks;++block)
     if (selected[block])
       copy_indices[block].resize(n_levels);
-  
+
 // Building the prolongation matrices starts here!
-  
+
                                   // reset the size of the array of
                                   // matrices. call resize(0) first,
                                   // in order to delete all elements
@@ -126,7 +294,7 @@ void MGTransferBlockBase::build_matrices (
                                   // for a cell and one of its children
   std::vector<unsigned int> dof_indices_parent (dofs_per_cell);
   std::vector<unsigned int> dof_indices_child (dofs_per_cell);
-  
+
                                   // for each level: first build the
                                   // sparsity pattern of the matrices
                                   // and then build the matrices
@@ -159,7 +327,7 @@ void MGTransferBlockBase::build_matrices (
                      0, false);
 
       prolongation_sparsities[level]->collect_sizes();
-      
+
       for (typename MGDoFHandler<dim,spacedim>::cell_iterator cell=mg_dof.begin(level);
           cell != mg_dof.end(level); ++cell)
        if (cell->has_children())
@@ -175,7 +343,7 @@ void MGTransferBlockBase::build_matrices (
                                                 // this child
                const FullMatrix<double> &prolongation
                  = mg_dof.get_fe().get_prolongation_matrix (child, cell->refinement_case());
-           
+
                cell->child(child)->get_mg_dof_indices (dof_indices_child);
 
                                                 // now tag the entries in the
@@ -214,7 +382,7 @@ void MGTransferBlockBase::build_matrices (
                                                 // this child
                const FullMatrix<double> &prolongation
                  = mg_dof.get_fe().get_prolongation_matrix (child, cell->refinement_case());
-           
+
                cell->child(child)->get_mg_dof_indices (dof_indices_child);
 
                                                 // now set the entries in the
@@ -246,13 +414,13 @@ void MGTransferBlockSelect<number>::build_matrices (
 {
   const FiniteElement<dim>& fe = mg_dof.get_fe();
   unsigned int n_blocks = mg_dof.get_fe().n_blocks();
-  
+
   selected_block = select;
   selected.resize(n_blocks, false);
   selected[select] = true;
-  
+
   MGTransferBlockBase::build_matrices (dof, mg_dof);
-  
+
   std::vector<unsigned int> global_dof_indices (fe.dofs_per_cell);
   std::vector<unsigned int> level_dof_indices  (fe.dofs_per_cell);
   for (int level=dof.get_tria().n_levels()-1; level>=0; --level)
@@ -261,7 +429,7 @@ void MGTransferBlockSelect<number>::build_matrices (
        level_cell = mg_dof.begin_active(level);
       const typename MGDoFHandler<dim,spacedim>::active_cell_iterator
        level_end  = mg_dof.end_active(level);
-      
+
                                       // Compute coarse level right hand side
                                       // by restricting from fine level.
       for (; level_cell!=level_end; ++level_cell)
@@ -273,7 +441,7 @@ void MGTransferBlockSelect<number>::build_matrices (
                                           // numbering
          global_cell.get_dof_indices(global_dof_indices);
          level_cell->get_mg_dof_indices (level_dof_indices);
-         
+
          for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
            {
              const unsigned int block = fe.system_to_block_index(i).first;
@@ -298,18 +466,18 @@ void MGTransferBlock<number>::build_matrices (
 {
   const FiniteElement<dim>& fe = mg_dof.get_fe();
   unsigned int n_blocks = mg_dof.get_fe().n_blocks();
-  
+
   if (sel.size() != 0)
     {
       Assert(sel.size() == n_blocks,
             ExcDimensionMismatch(sel.size(), n_blocks));
       selected = sel;
-    }  
+    }
   if (selected.size() == 0)
     selected = std::vector<bool> (n_blocks, true);
 
   MGTransferBlockBase::build_matrices (dof, mg_dof);
-  
+
   std::vector<unsigned int> global_dof_indices (fe.dofs_per_cell);
   std::vector<unsigned int> level_dof_indices  (fe.dofs_per_cell);
   for (int level=dof.get_tria().n_levels()-1; level>=0; --level)
@@ -318,7 +486,7 @@ void MGTransferBlock<number>::build_matrices (
        level_cell = mg_dof.begin_active(level);
       const typename MGDoFHandler<dim,spacedim>::active_cell_iterator
        level_end  = mg_dof.end_active(level);
-      
+
                                       // Compute coarse level right hand side
                                       // by restricting from fine level.
       for (; level_cell!=level_end; ++level_cell)
@@ -330,7 +498,7 @@ void MGTransferBlock<number>::build_matrices (
                                           // numbering
          global_cell.get_dof_indices(global_dof_indices);
          level_cell->get_mg_dof_indices (level_dof_indices);
-         
+
          for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
            {
              const unsigned int block = fe.system_to_block_index(i).first;
index 7eadf5861b3cab71bb59d67586eebffd3d591afc..7d2ae176eb4ada772521cc55ef3b8cf4977473f7 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 #include <multigrid/mg_tools.h>
 
 #include <algorithm>
+#include <numeric>
 
 DEAL_II_NAMESPACE_OPEN
 
 
+namespace
+{
+                                    /**
+                                     * Adjust block-vectors on all
+                                     * levels to correct size.  Count
+                                     * the numbers of degrees of
+                                     * freedom on each level
+                                     * component-wise. Then, assign
+                                     * each block of @p vector the
+                                     * corresponding size.
+                                     *
+                                     * The boolean field @p selected
+                                     * allows restricting this
+                                     * operation to certain
+                                     * components. In this case, @p
+                                     * vector will only have as many
+                                     * blocks as there are true
+                                     * values in @p selected (no
+                                     * blocks of length zero are
+                                     * padded in). If this argument
+                                     * is omitted, all blocks will be
+                                     * considered.
+                                     *
+                                     * Degrees of freedom must be
+                                     * sorted by component in order
+                                     * to obtain reasonable results
+                                     * from this function.
+                                     *
+                                     * The argument
+                                     * @p target_component allows to
+                                     * re-sort and group components
+                                     * as in
+                                     * DoFRenumbering::component_wise.
+                                     *
+                                     *
+                                     */
+  template <int dim, typename number, int spacedim>
+  void
+  reinit_vector_by_components (
+    const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+    MGLevelObject<BlockVector<number> > &v,
+    const std::vector<bool> &sel,
+    const std::vector<unsigned int> &target_comp,
+    std::vector<std::vector<unsigned int> >& ndofs)
+  {
+    std::vector<bool> selected=sel;
+    std::vector<unsigned int> target_component=target_comp;
+    const unsigned int ncomp = mg_dof.get_fe().n_components();
+
+                                    // If the selected and
+                                    // target_component have size 0,
+                                    // they must be replaced by default
+                                    // values.
+                                    //
+                                    // Since we already made copies
+                                    // directly after this function was
+                                    // called, we use the arguments
+                                    // directly.
+    if (target_component.size() == 0)
+      {
+       target_component.resize(ncomp);
+       for (unsigned int i=0;i<ncomp;++i)
+         target_component[i] = i;
+      }
+
+                                    // If selected is an empty vector,
+                                    // all components are selected.
+    if (selected.size() == 0)
+      {
+       selected.resize(target_component.size());
+       std::fill_n (selected.begin(), ncomp, false);
+       for (unsigned int i=0;i<target_component.size();++i)
+         selected[target_component[i]] = true;
+      }
+
+    Assert (selected.size() == target_component.size(),
+           ExcDimensionMismatch(selected.size(), target_component.size()));
+
+                                    // Compute the number of blocks needed
+    const unsigned int n_selected
+      = std::accumulate(selected.begin(),
+                       selected.end(),
+                       0U);
+
+    if (ndofs.size() == 0)
+      {
+       std::vector<std::vector<unsigned int> >
+         new_dofs(mg_dof.get_tria().n_levels(),
+                  std::vector<unsigned int>(target_component.size()));
+       std::swap(ndofs, new_dofs);
+       MGTools::count_dofs_per_component (mg_dof, ndofs,
+                                          true, target_component);
+      }
+
+    for (unsigned int level=v.get_minlevel();
+        level<=v.get_maxlevel();++level)
+      {
+       v[level].reinit(n_selected, 0);
+       unsigned int k=0;
+       for (unsigned int i=0;i<selected.size() && (k<v[level].n_blocks());++i)
+         {
+           if (selected[i])
+             {
+               v[level].block(k++).reinit(ndofs[level][i]);
+             }
+           v[level].collect_sizes();
+         }
+      }
+  }
+
+
+                                    /**
+                                     * Adjust vectors on all levels
+                                     * to correct size.  Count the
+                                     * numbers of degrees of freedom
+                                     * on each level component-wise
+                                     * in a single component. Then,
+                                     * assign @p vector the
+                                     * corresponding size.
+                                     *
+                                     * The boolean field @p selected
+                                     * may be nonzero in a single
+                                     * component, indicating the
+                                     * block of a block vector the
+                                     * argument @p v corresponds to.
+                                     *
+                                     * Degrees of freedom must be
+                                     * sorted by component in order
+                                     * to obtain reasonable results
+                                     * from this function.
+                                     *
+                                     * The argument
+                                     * @p target_component allows to
+                                     * re-sort and group components
+                                     * as in
+                                     * DoFRenumbering::component_wise.
+                                     */
+  template <int dim, typename number, int spacedim>
+  void
+  reinit_vector_by_components (
+    const dealii::MGDoFHandler<dim,spacedim> &mg_dof,
+    MGLevelObject<dealii::Vector<number> > &v,
+    const std::vector<bool> &selected,
+    const std::vector<unsigned int> &target_component,
+    std::vector<std::vector<unsigned int> >& ndofs)
+  {
+    Assert (selected.size() == target_component.size(),
+           ExcDimensionMismatch(selected.size(), target_component.size()));
+
+                                    // Compute the number of blocks needed
+#ifdef DEBUG
+    const unsigned int n_selected
+      = std::accumulate(selected.begin(),
+                       selected.end(),
+                       0U);
+    Assert(n_selected == 1, ExcDimensionMismatch(n_selected, 1));
+#endif
+
+    unsigned int selected_block = 0;
+    while (!selected[selected_block])
+      ++selected_block;
+
+    if (ndofs.size() == 0)
+      {
+       std::vector<std::vector<unsigned int> >
+         new_dofs(mg_dof.get_tria().n_levels(),
+                  std::vector<unsigned int>(target_component.size()));
+       std::swap(ndofs, new_dofs);
+       MGTools::count_dofs_per_component (mg_dof, ndofs,
+                                          true, target_component);
+      }
+
+    for (unsigned int level=v.get_minlevel();
+        level<=v.get_maxlevel();++level)
+      {
+       v[level].reinit(ndofs[level][selected_block]);
+      }
+  }
+}
+
+
+template <typename number>
+template <int dim, class InVector, int spacedim>
+void
+MGTransferSelect<number>::do_copy_to_mg (
+  const MGDoFHandler<dim,spacedim>&        mg_dof_handler,
+  MGLevelObject<Vector<number> >& dst,
+  const InVector&                 src,
+  const unsigned int              offset) const
+{
+  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
+
+  const unsigned int dofs_per_cell = fe.dofs_per_cell;
+
+                                  // set the elements of the vectors
+                                  // on all levels to zero
+  unsigned int minlevel = dst.get_minlevel();
+  unsigned int maxlevel = dst.get_maxlevel();
+
+  dst=0;
+
+  Assert(sizes.size()==mg_dof_handler.get_tria().n_levels(),
+        ExcMatricesNotBuilt());
+
+  reinit_vector_by_components(mg_dof_handler, dst, mg_selected,
+                             mg_target_component, sizes);
+
+  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
+  std::vector<unsigned int> level_dof_indices  (dofs_per_cell);
+
+                                  // Build a vector of the selected
+                                  // indices, since traversing all
+                                  // indices on each cell is too
+                                  // slow.
+  std::vector<unsigned int> selected_indices;
+  selected_indices.reserve(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    if (mg_selected[mg_target_component[fe.system_to_component_index(i).first]])
+      selected_indices.push_back(i);
+  unsigned int selected_component
+    = mg_target_component[fe.system_to_component_index(selected_indices[0]).first];
+
+                                  // traverse the grid top-down
+                                  // (i.e. starting with the most
+                                  // refined grid). this way, we can
+                                  // always get that part of one
+                                  // level of the output vector which
+                                  // corresponds to a region which is
+                                  // more refined, by restriction of
+                                  // the respective vector on the
+                                  // next finer level, which we then
+                                  // already have built.
+  for (int level=maxlevel; level>=static_cast<signed int>(minlevel); --level)
+    {
+
+      typename MGDoFHandler<dim,spacedim>::active_cell_iterator
+       level_cell = mg_dof_handler.begin_active(level);
+      const typename MGDoFHandler<dim,spacedim>::active_cell_iterator
+       level_end  = mg_dof_handler.end_active(level);
+                                      // Compute coarse level right hand side
+                                      // by restricting from fine level.
+      for (; level_cell!=level_end; ++level_cell)
+       {
+                                          // get the dof numbers of
+                                          // this cell for the global
+                                          // and the level-wise
+                                          // numbering
+         level_cell->get_dof_indices(global_dof_indices);
+         level_cell->get_mg_dof_indices (level_dof_indices);
+
+                                          // transfer the global
+                                          // defect in the vector
+                                          // into the level-wise one
+         const unsigned int level_start
+           = mg_component_start[level][selected_component];
+         const typename std::vector<unsigned int>::const_iterator
+           end = selected_indices.end();
+
+         for (typename std::vector<unsigned int>::const_iterator
+                i=selected_indices.begin();
+              i != end; ++i)
+           {
+             dst[level](level_dof_indices[*i] - level_start)
+               = src(global_dof_indices[*i] - offset);
+           }
+       }
+
+                                      // for that part of the level
+                                      // which is further refined:
+                                      // get the defect by
+                                      // restriction of the defect on
+                                      // one level higher
+      if (static_cast<unsigned int>(level) < maxlevel)
+       {
+         restrict_and_add (level+1, dst[level], dst[level+1]);
+       }
+    }
+}
+
 
 template <int dim, int spacedim>
 void MGTransferComponentBase::build_matrices (
@@ -51,7 +331,7 @@ void MGTransferComponentBase::build_matrices (
       Assert (target_component.size() == mg_dof.get_fe().n_components(),
              ExcDimensionMismatch(target_component.size(),
                                   mg_dof.get_fe().n_components()));
-      
+
       for (unsigned int i=0;i<target_component.size();++i)
        {
          Assert(i<target_component.size(),
@@ -70,14 +350,14 @@ void MGTransferComponentBase::build_matrices (
       Assert (mg_target_component.size() == mg_dof.get_fe().n_components(),
              ExcDimensionMismatch(mg_target_component.size(),
                                   mg_dof.get_fe().n_components()));
-      
+
       for (unsigned int i=0;i<mg_target_component.size();++i)
        {
          Assert(i<mg_target_component.size(),
                 ExcIndexRange(i,0,mg_target_component.size()));
        }
     }
-  
+
   const FiniteElement<dim>& fe = mg_dof.get_fe();
 
                                   // Effective number of components
@@ -89,10 +369,10 @@ void MGTransferComponentBase::build_matrices (
     *std::max_element(mg_target_component.begin(), mg_target_component.end()) + 1;
   const unsigned int dofs_per_cell = fe.dofs_per_cell;
   const unsigned int n_levels      = mg_dof.get_tria().n_levels();
-  
+
   Assert (mg_selected.size() == fe.n_components(),
          ExcDimensionMismatch(mg_selected.size(), fe.n_components()));
-  
+
                                   // Compute the lengths of all blocks
   sizes.resize(n_levels);
   MGTools::count_dofs_per_component(mg_dof, sizes, true, mg_target_component);
@@ -134,9 +414,9 @@ void MGTransferComponentBase::build_matrices (
                                   // too slow
 
   copy_to_and_from_indices.resize(n_levels);
-  
+
 // Building the prolongation matrices starts here!
-  
+
                                   // reset the size of the array of
                                   // matrices. call resize(0) first,
                                   // in order to delete all elements
@@ -163,7 +443,7 @@ void MGTransferComponentBase::build_matrices (
                                   // for a cell and one of its children
   std::vector<unsigned int> dof_indices_parent (dofs_per_cell);
   std::vector<unsigned int> dof_indices_child (dofs_per_cell);
-  
+
                                   // for each level: first build the
                                   // sparsity pattern of the matrices
                                   // and then build the matrices
@@ -196,7 +476,7 @@ void MGTransferComponentBase::build_matrices (
                      0, false);
 
       prolongation_sparsities[level]->collect_sizes();
-      
+
       for (typename MGDoFHandler<dim,spacedim>::cell_iterator cell=mg_dof.begin(level);
           cell != mg_dof.end(level); ++cell)
        if (cell->has_children())
@@ -210,7 +490,7 @@ void MGTransferComponentBase::build_matrices (
                                                 // this child
                const FullMatrix<double> &prolongation
                  = mg_dof.get_fe().get_prolongation_matrix (child, cell->refinement_case());
-           
+
                cell->child(child)->get_mg_dof_indices (dof_indices_child);
 
                                                 // now tag the entries in the
@@ -247,7 +527,7 @@ void MGTransferComponentBase::build_matrices (
                                                 // this child
                const FullMatrix<double> &prolongation
                  = mg_dof.get_fe().get_prolongation_matrix (child, cell->refinement_case());
-           
+
                cell->child(child)->get_mg_dof_indices (dof_indices_child);
 
                                                 // now set the entries in the
@@ -281,10 +561,10 @@ void MGTransferSelect<number>::build_matrices (
 {
   const FiniteElement<dim>& fe = mg_dof.get_fe();
   unsigned int ncomp = mg_dof.get_fe().n_components();
-  
+
   target_component = t_component;
   mg_target_component = mg_t_component;
-  
+
   selected_component = select;
   mg_selected_component = mg_select;
   selected.resize(ncomp, false);
@@ -303,7 +583,7 @@ void MGTransferSelect<number>::build_matrices (
          break;
        }
     }
-    
+
   for (unsigned int i=0;i<mg_target_component.size();++i)
     {
       if (mg_target_component[i] == mg_select)
@@ -311,9 +591,9 @@ void MGTransferSelect<number>::build_matrices (
          mg_selected_component = i;
          break;
        }
-    }    
+    }
   MGTransferComponentBase::build_matrices (dof, mg_dof);
-  
+
   std::vector<unsigned int> global_dof_indices (fe.dofs_per_cell);
   std::vector<unsigned int> level_dof_indices  (fe.dofs_per_cell);
   for (int level=dof.get_tria().n_levels()-1; level>=0; --level)
@@ -322,7 +602,7 @@ void MGTransferSelect<number>::build_matrices (
        level_cell = mg_dof.begin_active(level);
       const typename MGDoFHandler<dim,spacedim>::active_cell_iterator
        level_end  = mg_dof.end_active(level);
-      
+
                                       // Compute coarse level right hand side
                                       // by restricting from fine level.
       for (; level_cell!=level_end; ++level_cell)
@@ -334,7 +614,7 @@ void MGTransferSelect<number>::build_matrices (
                                           // numbering
          global_cell.get_dof_indices(global_dof_indices);
          level_cell->get_mg_dof_indices (level_dof_indices);
-         
+
          for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
            {
              const unsigned int component
index 99866e5c78f3403bde8bb40321c1643b51654088..b2ea34d8da2ed4607211054a5092cfcb309c561f 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+namespace
+{
+                                    /**
+                                     * Adjust vectors on all levels to
+                                     * correct size.  Here, we just
+                                     * count the numbers of degrees
+                                     * 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>
+  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();
+        level<=v.get_maxlevel();++level)
+      {
+       unsigned int n = mg_dof.n_dofs (level);
+       v[level].reinit(n);
+      }
+
+  }
+
+
+                                    /**
+                                     * Adjust vectors on all levels to
+                                     * correct size.  Here, we just
+                                     * count the numbers of degrees
+                                     * of freedom on each level and
+                                     * @p reinit each level vector
+                                     * 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>
+  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, target_component);
+
+    for (unsigned int level=v.get_minlevel();
+        level<=v.get_maxlevel();++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();
+      }
+  }
+}
+
+
+
+template <class VECTOR>
+template <int dim, class InVector, int spacedim>
+void
+MGTransferPrebuilt<VECTOR>::copy_to_mg (
+  const MGDoFHandler<dim,spacedim>& mg_dof_handler,
+  MGLevelObject<VECTOR>& dst,
+  const InVector& src) const
+{
+  reinit_vector(mg_dof_handler, component_to_block_map, dst);
+  bool first = true;
+  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
+    {
+      --level;
+
+      typedef std::map<unsigned int, unsigned int>::const_iterator IT;
+      for (IT i= copy_indices[level].begin();
+          i != copy_indices[level].end();++i)
+       dst[level](i->second) = src(i->first);
+                                      // For non-DG: degrees of
+                                      // freedom in the refinement
+                                      // face may need special
+                                      // attention, since they belong
+                                      // to the coarse level, but
+                                      // have fine level basis
+                                      // functions
+      if (!first)
+       restrict_and_add (level+1, dst[level], dst[level+1]);
+      first = false;
+    }
+}
+
 
 
 template <typename VECTOR>
@@ -233,6 +342,7 @@ void MGTransferPrebuilt<VECTOR>::build_matrices (
 }
 
 
+//TODO: Use template expander script
 
 template
 void MGTransferPrebuilt<Vector<float> >::build_matrices<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.