]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
remove competing local indices
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 7 Mar 2013 21:40:48 +0000 (21:40 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 7 Mar 2013 21:40:48 +0000 (21:40 +0000)
git-svn-id: https://svn.dealii.org/trunk@28802 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/meshworker/dof_info.h
deal.II/include/deal.II/meshworker/simple.h

index 0f5ba9d67df73fad41ac8de783ebbdb8456230dc..46116bd4e2aeb8d7102d6790716d47dd36ea49eb 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
+//    Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -100,9 +100,8 @@ namespace MeshWorker
     std::vector<unsigned int> indices;
 
     /**
-     * The DoF indices on the
-     * current cell, organized by
-     * local blocks
+     * The DoF indices on the current cell, organized by local blocks.
+     * The size of this vector is zero, unless local blocks are used.
      */
     std::vector<std::vector<unsigned int> > indices_by_block;
 
@@ -129,9 +128,7 @@ namespace MeshWorker
     void reinit(const DHCellIterator &c);
 
     /**
-     * Set the current face and
-     * fill @p indices if the #cell
-     * changed.
+     * Set the current face and fill @p indices if the #cell changed.
      */
     template <class DHCellIterator, class DHFaceIterator>
     void reinit(const DHCellIterator &c,
@@ -139,9 +136,8 @@ namespace MeshWorker
                 const unsigned int n);
 
     /**
-     * Set the current subface
-     * and fill @p indices if the
-     * #cell changed.
+     * Set the current subface and fill @p indices if the #cell
+     * changed.
      */
     template <class DHCellIterator, class DHFaceIterator>
     void reinit(const DHCellIterator &c,
@@ -150,19 +146,16 @@ namespace MeshWorker
                 const unsigned int s);
 
     /**
-     * Switch to a new face of the
-     * same cell. Does not change
-     * @p indices and does not reset
-     * data in LocalResults.
+     * Switch to a new face of the same cell. Does not change @p
+     * indices and does not reset data in LocalResults.
      */
     template <class DHFaceIterator>
     void set_face (const DHFaceIterator &f,
                    const unsigned int n);
+    
     /**
-     * Switch to a new subface of the
-     * same cell. Does not change
-     * @p indices and does not reset
-     * data in LocalResults.
+     * Switch to a new subface of the same cell. Does not change @p
+     * indices and does not reset data in LocalResults.
      */
     template <class DHFaceIterator>
     void set_subface (const DHFaceIterator &f,
@@ -178,13 +171,9 @@ namespace MeshWorker
     bool level_cell;
   private:
     /**
-     * Standard constructor, not
-     * setting any block
-     * indices. Use of this
-     * constructor is not
-     * recommended, but it is
-     * needed for the arrays in
-     * DoFInfoBox.
+     * Standard constructor, not setting any block indices. Use of
+     * this constructor is not recommended, but it is needed for the
+     * arrays in DoFInfoBox.
      */
     DoFInfo ();
 
@@ -194,9 +183,6 @@ namespace MeshWorker
     template <class DHCellIterator>
     void get_indices(const DHCellIterator &c);
 
-    /// Fill index vector with level indices
-    //void get_indices(const typename MGDoFHandler<dim, spacedim>::cell_iterator& c);
-
     /// Auxiliary vector
     std::vector<unsigned int> indices_org;
 
index 27aa1fe81abcb35f29367758a1f1f07108a42114..fb1d30fca7e26c8c57a02faeb86c23d21920309d 100644 (file)
@@ -177,6 +177,9 @@ namespace MeshWorker
       void initialize(const ConstraintMatrix &constraints);
 
       /**
+       * @deprecated This function is of no effect. Only the block info
+       * structure in DoFInfo is being used.
+       *
        * Store information on the local block structure. If the
        * assembler is inititialized with this function,
        * initialize_info() will generate one local matrix for each
@@ -196,15 +199,11 @@ namespace MeshWorker
       void initialize_local_blocks(const BlockIndices &local_indices);
 
       /**
-       * Initialize the local data
-       * in the DoFInfo object used
-       * later for assembling.
+       * Initialize the local data in the DoFInfo object used later
+       * for assembling.
        *
-       * The info object refers to
-       * a cell if
-       * <code>!face</code>, or
-       * else to an interior or
-       * boundary face.
+       * The info object refers to a cell if <code>!face</code>, or
+       * else to an interior or boundary face.
        */
       template <class DOFINFO>
       void initialize_info(DOFINFO &info, bool face) const;
@@ -242,12 +241,6 @@ namespace MeshWorker
        */
       SmartPointer<const ConstraintMatrix,MatrixSimple<MATRIX> > constraints;
 
-      /**
-       * The object containing the local block structure. Set by
-       * initialize_local_blocks() and used by assembling functions.
-       */
-      BlockIndices local_indices;
-
       /**
        * The smallest positive number that will be entered into the
        * global matrix. All smaller absolute values will be treated as
@@ -643,10 +636,8 @@ namespace MeshWorker
 
     template <class MATRIX>
     inline void
-    MatrixSimple<MATRIX>::initialize_local_blocks(const BlockIndices &b)
-    {
-      local_indices = b;
-    }
+    MatrixSimple<MATRIX>::initialize_local_blocks(const BlockIndices &)
+    {}
 
 
     template <class MATRIX >
@@ -654,7 +645,7 @@ namespace MeshWorker
     inline void
     MatrixSimple<MATRIX>::initialize_info(DOFINFO &info, bool face) const
     {
-      const unsigned int n = local_indices.size();
+      const unsigned int n = dof_info.indices_by_block.size();
 
       if (n == 0)
         info.initialize_matrices(1, face);
@@ -706,7 +697,7 @@ namespace MeshWorker
     {
       Assert(!info.level_cell, ExcMessage("Cell may not access level dofs"));
 
-      if (local_indices.size() == 0)
+      if (info.indices_by_block.size() == 0)
         assemble(info.matrix(0,false).matrix, info.indices, info.indices);
       else
         {
@@ -728,14 +719,14 @@ namespace MeshWorker
       Assert(!info1.level_cell, ExcMessage("Cell may not access level dofs"));
       Assert(!info2.level_cell, ExcMessage("Cell may not access level dofs"));
 
-      if (local_indices.size() == 0)
+      if (info1.indices_by_block.size() == 0 && info2.indices_by_block.size() == 0)
         {
           assemble(info1.matrix(0,false).matrix, info1.indices, info1.indices);
           assemble(info1.matrix(0,true).matrix, info1.indices, info2.indices);
           assemble(info2.matrix(0,false).matrix, info2.indices, info2.indices);
           assemble(info2.matrix(0,true).matrix, info2.indices, info1.indices);
         }
-      else
+      else if (info1.indices_by_block.size() != 0 && info2.indices_by_block.size() != 0)
         for (unsigned int k=0; k<info1.n_matrices(); ++k)
           {
             const unsigned int row = info1.matrix(k,false).row;
@@ -750,6 +741,10 @@ namespace MeshWorker
             assemble(info2.matrix(k,true).matrix,
                      info2.indices_by_block[row], info1.indices_by_block[column]);
           }
+      else
+       {
+         Assert(false, ExcNotImplemented());
+       }
     }
 
 

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.