]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Allow to say which components fit how into blocks.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 5 Jan 2010 17:27:07 +0000 (17:27 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 5 Jan 2010 17:27:07 +0000 (17:27 +0000)
git-svn-id: https://svn.dealii.org/trunk@20297 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_transfer.h
deal.II/deal.II/include/multigrid/mg_transfer.templates.h

index 600af6db4148c052456c573ba181e8584eb7009a..630cd295f71d6186bb031c01087577c022fa1035 100644 (file)
@@ -61,7 +61,7 @@ template <int dim, int spacedim> class MGDoFHandler;
  * @author Wolfgang Bangerth, Guido Kanschat, 1999-2004
  */
 template <class VECTOR>
-class MGTransferPrebuilt : public MGTransferBase<VECTOR> 
+class MGTransferPrebuilt : public MGTransferBase<VECTOR>
 {
   public:
                                     /**
@@ -140,12 +140,45 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                      OutVector& dst,
                      const MGLevelObject<VECTOR>& src) const;
 
+                                    /**
+                                     * If this object operates on
+                                     * BlockVector objects, we need
+                                     * to describe how the individual
+                                     * vector components are mapped
+                                     * to the blocks of a vector. For
+                                     * example, for a Stokes system,
+                                     * we have dim+1 vector
+                                     * components for velocity and
+                                     * pressure, but we may want to
+                                     * use block vectors with only
+                                     * two blocks for all velocities
+                                     * in one block, and the pressure
+                                     * variables in the other.
+                                     *
+                                     * By default, if this function
+                                     * is not called, block vectors
+                                     * have as many blocks as the
+                                     * finite element has vector
+                                     * components. However, this can
+                                     * be changed by calling this
+                                     * function with an array that
+                                     * describes how vector
+                                     * components are to be grouped
+                                     * into blocks. The meaning of
+                                     * the argument is the same as
+                                     * the one given to the
+                                     * DoFTools::count_dofs_per_component
+                                     * function.
+                                     */
+    void
+    set_component_to_block_map (const std::vector<unsigned int> &map);
+
                                     /**
                                      * Finite element does not
                                      * provide prolongation matrices.
                                      */
     DeclException0(ExcNoProlongation);
-    
+
                                     /**
                                      * Call @p build_matrices
                                      * function first.
@@ -156,7 +189,7 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      * Memory used by this object.
                                      */
     unsigned int memory_consumption () const;
-    
+
 
   private:
 
@@ -179,8 +212,8 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      * while row indices belong to the
                                      * child cell, i.e. the fine level.
                                      */
-    std::vector<std_cxx1x::shared_ptr<SparseMatrix<double> > > prolongation_matrices;    
-    
+    std::vector<std_cxx1x::shared_ptr<SparseMatrix<double> > > prolongation_matrices;
+
                                     /**
                                      * Mapping for the
                                      * <tt>copy_to/from_mg</tt>-functions.
@@ -190,6 +223,14 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
     std::vector<std::map<unsigned int, unsigned int> >
     copy_indices;
 
+                                    /**
+                                     * The vector that stores what
+                                     * has been given to the
+                                     * set_component_to_block_map()
+                                     * function.
+                                     */
+    std::vector<unsigned int> component_to_block_map;
+
                                     /**
                                      * Fill the two boolean vectors
                                      * #dofs_on_refinement_edge and
@@ -219,7 +260,7 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      *
                                      * @todo Clean up names
                                      */
-    std::vector<std::vector<bool> > dofs_on_refinement_boundary; 
+    std::vector<std::vector<bool> > dofs_on_refinement_boundary;
                                     /**
                                      * The constraints of the global
                                      * system.
index e22a2981d7b99f135265d65ea2835069140a57a2..23666ab3d614e3823feaf1bcf0a7a529c59fbaab 100644 (file)
@@ -40,7 +40,7 @@ MGTransferPrebuilt<VECTOR>::copy_to_mg (
   MGLevelObject<VECTOR>& dst,
   const InVector& src) const
 {
-  internal::mg::reinit_vector(mg_dof_handler, dst);
+  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;)
     {
@@ -112,6 +112,18 @@ MGTransferPrebuilt<VECTOR>::copy_from_mg_add (
       dst(i->first) += src[level](i->second);
 }
 
+
+
+template <class VECTOR>
+void
+MGTransferPrebuilt<VECTOR>::
+set_component_to_block_map (const std::vector<unsigned int> &map)
+{
+  component_to_block_map = map;
+}
+
+
+
 template <class VECTOR>
 template <int dim, int spacedim>
 void MGTransferPrebuilt<VECTOR>::find_dofs_on_refinement_edges (

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.