]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
change scope of class
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 12 Aug 2011 16:02:47 +0000 (16:02 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 12 Aug 2011 16:02:47 +0000 (16:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@24046 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/block_indices.h

index eb298bf2bea8f51c3202d10a38d33c558aae8220..0be0b045e369b5ef83ce261c677b7db9a781335a 100644 (file)
@@ -22,7 +22,19 @@ DEAL_II_NAMESPACE_OPEN
 
 
 /**
- * Class that manages the conversion of global indices into a block vector or
+ * @brief Auxiliary class aiding in the handling of block structures like in
+ * BlockVector or FESystem.
+ *
+ * The information obtained from this class falls into two
+ * groups. First, it is possible to obtain the number of blocks,
+ * namely size(), the block_size() for each block and the total_size()
+ * of the object described by the block indices, namely the length of
+ * the whole index set. These functions do not make any assumption on
+ * the ordering of the index set.
+ *
+ * If on the other hand the index set is ordered "by blocks", such
+ * that each block forms a consecutive set of indices, this
+ * class that manages the conversion of global indices into a block vector or
  * matrix to the local indices within this block. This is required, for
  * example, when you address a global element in a block vector and want to
  * know which element within which block this is. It is also useful if a
@@ -31,7 +43,7 @@ DEAL_II_NAMESPACE_OPEN
  *
  * @ingroup data
  * @see @ref GlossBlockLA "Block (linear algebra)"
- * @author Wolfgang Bangerth, Guido Kanschat, 2000, 2007
+ * @author Wolfgang Bangerth, Guido Kanschat, 2000, 2007, 2011
  */
 class BlockIndices : public Subscriptor
 {
@@ -43,7 +55,7 @@ class BlockIndices : public Subscriptor
                                      * @p n_blocks blocks and set
                                      * all block sizes to zero.
                                      */
-    BlockIndices (const unsigned int n_blocks = 0);
+    BlockIndices (/*const unsigned int n_blocks = 0*/);
 
                                     /**
                                      * Constructor. Initialize the
@@ -53,7 +65,13 @@ class BlockIndices : public Subscriptor
                                      * size of the vector
                                      */
     BlockIndices (const std::vector<unsigned int> &n);
-
+    
+                                    /**
+                                     * Specialized constructor for a
+                                     * structure with blocks of equal size.
+                                     */
+    BlockIndices(const unsigned int n_blocks, const unsigned int block_size);
+    
                                     /**
                                      * Reinitialize the number of
                                      * blocks and assign each block
@@ -74,22 +92,9 @@ class BlockIndices : public Subscriptor
     inline void reinit (const std::vector<unsigned int> &n);
     
                                     /**
-                                     * Return the block and the
-                                     * index within that block
-                                     * for the global index @p i. The
-                                     * first element of the pair is
-                                     * the block, the second the
-                                     * index within it.
-                                     */
-    std::pair<unsigned int,unsigned int>
-    global_to_local (const unsigned int i) const;
-
-                                    /**
-                                     * Return the global index of
-                                     * @p index in block @p block.
+                                     * @name Size information
                                      */
-    unsigned int local_to_global (const unsigned int block,
-                                 const unsigned int index) const;
+                                    //@{
 
                                     /**
                                      * Number of blocks in index field.
@@ -110,10 +115,46 @@ class BlockIndices : public Subscriptor
                                      */
     unsigned int block_size (const unsigned int i) const;
 
+                                    //@}
+
+                                    /**
+                                     * @name Index conversion
+                                     *
+                                     * Functions in this group
+                                     * assume an object, which
+                                     * was created after sorting by
+                                     * block, such that each block
+                                     * forms a set of consecutive
+                                     * indices in the object.
+                                     * If applied to other objects,
+                                     * the numbers obtained from
+                                     * these functions are meaningless.
+                                     */
+                                    //@{
+
+                                    /**
+                                     * Return the block and the
+                                     * index within that block
+                                     * for the global index @p i. The
+                                     * first element of the pair is
+                                     * the block, the second the
+                                     * index within it.
+                                     */
+    std::pair<unsigned int,unsigned int>
+    global_to_local (const unsigned int i) const;
+
+                                    /**
+                                     * Return the global index of
+                                     * @p index in block @p block.
+                                     */
+    unsigned int local_to_global (const unsigned int block,
+                                 const unsigned int index) const;
+
                                     /**
                                      * The start index of the ith block.
                                      */
     unsigned int block_start (const unsigned int i) const;
+                                    //@}
     
                                     /**
                                      * Copy operator.
@@ -123,8 +164,8 @@ class BlockIndices : public Subscriptor
                                     /**
                                      * Compare whether two objects
                                      * are the same, i.e. whether the
-                                     * starting indices of all blocks
-                                     * are equal.
+                                     * number of blocks and the sizes
+                                     * of all blocks are equal.
                                      */
     bool operator == (const BlockIndices &b) const;
     
@@ -162,6 +203,27 @@ class BlockIndices : public Subscriptor
 };
 
 
+/**
+ * Output operator for BlockIndices
+ *
+ * @ref BlockIndices
+ * @author Guido Kanschat
+ * @date 2011 
+ */
+template <class STREAM>
+STREAM&
+operator << (STREAM& s, const BlockIndices& bi)
+{
+  const unsigned int n = bi.size();
+  s << n << ":[";
+  if (n>0)
+    s << bi.block_size(0);
+  for (unsigned int i=1;i<n;++i)
+    s << ' ' << bi.block_size(i);
+  s << ']';
+  return s;
+}
+
 
 template <typename MatrixType> class BlockMatrixBase;
 template <typename SparsityType> class BlockSparsityPatternBase;
@@ -249,55 +311,68 @@ const bool IsBlockMatrix<MatrixType>::value;
 /* ---------------------- template and inline functions ------------------- */
 
 inline
-BlockIndices::BlockIndices (const unsigned int n_blocks)
-               :
-               n_blocks(n_blocks),
-               start_indices(n_blocks+1)
+void
+BlockIndices::reinit (const unsigned int nb,
+                     const unsigned int block_size)
 {
+  n_blocks = nb;
+  start_indices.resize(nb);
   for (unsigned int i=0; i<=n_blocks; ++i)
-    start_indices[i] = 0;
+    start_indices[i] = i * block_size;
 }
 
 
 
 inline
-BlockIndices::BlockIndices (const std::vector<unsigned int> &n)
-               :
-               n_blocks(n.size()),
-               start_indices(n.size()+1)
+void
+BlockIndices::reinit (const std::vector<unsigned int> &n)
 {
-  reinit (n);
+  if (start_indices.size() != n.size()+1)
+    {
+      n_blocks = n.size();
+      start_indices.resize(n_blocks+1);
+    }
+  start_indices[0] = 0;
+  for (unsigned int i=1; i<=n_blocks; ++i)
+    start_indices[i] = start_indices[i-1] + n[i-1];
 }
 
 
+inline
+BlockIndices::BlockIndices (/*const unsigned int n_blocks*/)
+               :
+               n_blocks(0/*n_blocks*/)//,
+//             start_indices(n_blocks+1, 0)
+{}
+
+
 
 inline
-void
-BlockIndices::reinit (const unsigned int n_blocks,
-                     const unsigned int n_elements_per_block)
+BlockIndices::BlockIndices (
+  const unsigned int n_blocks,
+  const unsigned int block_size)
+               :
+               n_blocks(n_blocks),
+               start_indices(n_blocks+1)
 {
-  const std::vector<unsigned int> v(n_blocks, n_elements_per_block);
-  reinit (v);
+  for (unsigned int i=0; i<=n_blocks; ++i)
+    start_indices[i] = i * block_size;
 }
 
 
 
 inline
-void
-BlockIndices::reinit (const std::vector<unsigned int> &n)
+BlockIndices::BlockIndices (const std::vector<unsigned int> &n)
+               :
+               n_blocks(n.size()),
+               start_indices(n.size()+1)
 {
-  if (start_indices.size() != n.size()+1)
-    {
-      n_blocks = n.size();
-      start_indices.resize(n_blocks+1);
-    }
-  start_indices[0] = 0;
-  for (unsigned int i=1; i<=n_blocks; ++i)
-    start_indices[i] = start_indices[i-1] + n[i-1];
+  reinit (n);
 }
 
 
 
+
 inline
 std::pair<unsigned int,unsigned int>
 BlockIndices::global_to_local (const unsigned int i) const 

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.