Add an argument to IndexSet::add_indices. Use this to make a function in BlockVector...

git-svn-id: https://svn.dealii.org/trunk@29607 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index a29684ee60baabe19dbcb45118ccdc959688421e..291868facf2874b9d35048b1b6b112161cd499e7 100644 (file)
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -130,6 +130,13 @@ this function.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: The IndexSet::add_indices function that takes another IndexSet
+object now has an additional argument <code>offset</code> that can be used
+to offset the indices of first argument.
+<br>
+(Wolfgang Bangerth, 2013/05/25)
+</li>
+
 <li> New: ConstraintMatrix::distribute is now also implemented for
 arguments of type PETScWrappers::MPI::BlockVector.
 <br>

diff --git a/deal.II/include/deal.II/base/index_set.h b/deal.II/include/deal.II/base/index_set.h
index aff873cba72e080a40943198687a2193719e0c4f..52a82c6ba16b6fccf018509513c2c3e275e892ee 100644 (file)
--- a/deal.II/include/deal.II/base/index_set.h
+++ b/deal.II/include/deal.II/base/index_set.h
@@ -128,8 +128,23 @@ public:
    * Add the given IndexSet @p other to the
    * current one, constructing the union of
    * *this and @p other.
+   *
+   * If the @p offset argument is nonzero, then every
+   * index in @p other is shifted by @p offset before being
+   * added to the current index set. This allows to construct,
+   * for example, one index set from several others that are
+   * supposed to represent index sets corresponding to
+   * different ranges (e.g., when constructing the set of
+   * nonzero entries of a block vector from the sets of nonzero
+   * elements of the individual blocks of a vector).
+   *
+   * This function will generate an exception if any of the
+   * (possibly shifted) indices of the @p other index set
+   * lie outside the range <code>[0,size())</code> represented
+   * by the current object.
    */
-  void add_indices(const IndexSet &other);
+  void add_indices(const IndexSet &other,
+                   const unsigned int offset = 0);
 
   /**
    * Return whether the specified
@@ -767,16 +782,17 @@ IndexSet::add_indices (const ForwardIterator &begin,
 
 inline
 void
-IndexSet::add_indices(const IndexSet &other)
+IndexSet::add_indices(const IndexSet &other,
+                      const unsigned int offset)
 {
-  if (this == &other)
+  if ((this == &other) && (offset == 0))
     return;
 
   for (std::vector<Range>::iterator range = other.ranges.begin();
        range != other.ranges.end();
        ++range)
     {
-      add_range(range->begin, range->end);
+      add_range(range->begin+offset, range->end+offset);
     }
 
   compress();

diff --git a/deal.II/include/deal.II/lac/block_vector_base.h b/deal.II/include/deal.II/lac/block_vector_base.h
index a512e1aa4424ac67cf4c2dbf816a9c31f37d1867..9a2344c042f65f827c4053aac15989d444711c99 100644 (file)
--- a/deal.II/include/deal.II/lac/block_vector_base.h
+++ b/deal.II/include/deal.II/lac/block_vector_base.h
@@ -1810,16 +1810,12 @@ BlockVectorBase<VectorType>::locally_owned_elements () const
 {
   IndexSet is (size());
 
-  // copy index sets from blocks into the global one
+  // copy index sets from blocks into the global one, shifted
+  // by the appropriate amount for each block
   for (unsigned int b=0; b<n_blocks(); ++b)
     {
       IndexSet x = block(b).locally_owned_elements();
-
-      //TODO: This can surely be made more efficient by just shifting
-      // x
-      for (unsigned int i=0; i<block(b).size(); ++i)
-        if (x.is_element(i))
-          is.add_index(block_indices.local_to_global(b,i));
+      is.add_indices(x, block_indices.block_start(b));
     }
 
   is.compress();

diff --git a/tests/base/index_set_20_offset.cc b/tests/base/index_set_20_offset.cc
new file mode 100644 (file)
index 0000000..567f62e
--- /dev/null
+++ b/tests/base/index_set_20_offset.cc
@@ -0,0 +1,94 @@
+//-----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2010, 2013 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-----------------------------------------------------------------------------
+
+// test IndexSet::add_indices(IndexSet)
+
+#include "../tests.h"
+#include <iomanip>
+#include <fstream>
+#include <cmath>
+#include <stdlib.h>
+
+#include <deal.II/base/index_set.h>
+
+void testor(IndexSet & a,
+           IndexSet & other,
+           unsigned int offset,
+           bool verbose)
+{
+  IndexSet merged(a);
+
+  merged.add_indices(other, offset);
+
+  if (verbose)
+    {
+      deallog << "Original index set: " << std::endl;
+      a.print(deallog);
+      deallog << "other index set: " << std::endl;
+      other.print(deallog);
+      deallog << "merged index set: " << std::endl;
+      merged.print(deallog);
+    }
+
+  for (unsigned int i=0;i<merged.size();++i)
+    {
+      Assert(
+       merged.is_element(i)
+       ==
+       (a.is_element(i) || (i>=offset && other.is_element(i-offset))),
+       ExcInternalError());
+    }
+}
+
+
+
+
+void test()
+{
+  const int size = 10;
+  
+  IndexSet empty(size);
+  IndexSet id(size);
+
+  id.add_index(3);
+  id.add_index(4);
+  id.add_index(7);
+  
+  deallog << "* add empty: " << std::endl;
+  testor(id, empty, 2, true);
+
+  deallog << "* add self: " << std::endl;
+  testor(id, id, 2, true);
+
+  deallog << "* add id2: " << std::endl;
+  IndexSet id2(size);
+  id2.add_index(0);
+  id2.add_index(2);
+  id2.add_index(3);
+  testor(id, id2, 3, true);
+}
+
+
+
+
+
+
+int main()
+{
+  std::ofstream logfile("index_set_20_offset/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test ();
+}

diff --git a/tests/base/index_set_20_offset/cmp/generic b/tests/base/index_set_20_offset/cmp/generic
new file mode 100644 (file)
index 0000000..7d32353
--- /dev/null
+++ b/tests/base/index_set_20_offset/cmp/generic
@@ -0,0 +1,22 @@
+
+DEAL::* add empty: 
+DEAL::Original index set: 
+DEAL::{[3,4], 7}
+DEAL::other index set: 
+DEAL::{}
+DEAL::merged index set: 
+DEAL::{[3,4], 7}
+DEAL::* add self: 
+DEAL::Original index set: 
+DEAL::{[3,4], 7}
+DEAL::other index set: 
+DEAL::{[3,4], 7}
+DEAL::merged index set: 
+DEAL::{[3,7], 9}
+DEAL::* add id2: 
+DEAL::Original index set: 
+DEAL::{[3,4], 7}
+DEAL::other index set: 
+DEAL::{0, [2,3]}
+DEAL::merged index set: 
+DEAL::{[3,7]}