]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add locally_owned_elements() to all vector classes.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 24 May 2013 19:02:27 +0000 (19:02 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 24 May 2013 19:02:27 +0000 (19:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@29578 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/lac/block_vector_base.h
deal.II/include/deal.II/lac/parallel_vector.h
deal.II/include/deal.II/lac/petsc_vector_base.h
deal.II/include/deal.II/lac/trilinos_vector_base.h
deal.II/include/deal.II/lac/vector.h
deal.II/include/deal.II/lac/vector.templates.h

index 366a4d36187928c4542bd91ff48cf9266ba4b0cb..fcb38a37b012611862a2b254e622c7eea6edd1df 100644 (file)
@@ -65,6 +65,15 @@ this function.
 
 <ol>
 
+  <li> New: All vector classes now have a member function
+  <code>locally_owned_elements</code> that returns an index
+  set indicating which elements of this vector the current
+  processor owns.
+  <br>
+  (Wolfgang Bangerth, 2013/05/24)
+  </li>
+
+
   <li> New: A new element FE_Q_iso_Q1 has been implemented that is defined by
   a subdivision of the element into smaller Q1 elements. An element of order
   @p p is similar to FE_Q of degree @p p with the same numbering of degrees of
index 68b469c76f0ecdcda04edba7ec19b179adf4ea18..a512e1aa4424ac67cf4c2dbf816a9c31f37d1867 100644 (file)
@@ -846,6 +846,26 @@ public:
    */
   unsigned int size () const;
 
+  /**
+   * Return an index set that describes which elements of this vector
+   * are owned by the current processor. Note that this index set does
+   * not include elements this vector may store locally as ghost
+   * elements but that are in fact owned by another processor.
+   * As a consequence, the index sets returned on different
+   * processors if this is a distributed vector will form disjoint
+   * sets that add up to the complete index set.
+   * Obviously, if a vector is created on only one processor, then
+   * the result would satisfy
+   * @code
+   *   vec.locally_owned_elements() == complete_index_set (vec.size())
+   * @endcode
+   *
+   * For block vectors, this function returns the union of the
+   * locally owned elements of the individual blocks, shifted by
+   * their respective index offsets.
+   */
+  IndexSet locally_owned_elements () const;
+
   /**
    * Return an iterator pointing to
    * the first element.
@@ -1783,6 +1803,32 @@ BlockVectorBase<VectorType>::size () const
 
 
 
+template <class VectorType>
+inline
+IndexSet
+BlockVectorBase<VectorType>::locally_owned_elements () const
+{
+  IndexSet is (size());
+
+  // copy index sets from blocks into the global one
+  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.compress();
+
+  return is;
+}
+
+
+
 template <class VectorType>
 inline
 unsigned int
index aeef6b51f36cae9776209c3a3c813565b4a11cc5..29d0b864b6b1273b459fd0e9fc5fdce252bd117b 100644 (file)
@@ -458,6 +458,22 @@ namespace parallel
        */
       bool in_local_range (const types::global_dof_index global_index) const;
 
+      /**
+       * Return an index set that describes which elements of this vector
+       * are owned by the current processor. Note that this index set does
+       * not include elements this vector may store locally as ghost
+       * elements but that are in fact owned by another processor.
+       * As a consequence, the index sets returned on different
+       * processors if this is a distributed vector will form disjoint
+       * sets that add up to the complete index set.
+       * Obviously, if a vector is created on only one processor, then
+       * the result would satisfy
+       * @code
+       *   vec.locally_owned_elements() == complete_index_set (vec.size())
+       * @endcode
+       */
+      IndexSet locally_owned_elements () const;
+
       /**
        * Returns the number of ghost elements present on the vector.
        */
@@ -1429,6 +1445,20 @@ namespace parallel
     }
 
 
+    template <typename Number>
+    inline
+    IndexSet
+    Vector<Number>::locally_owned_elements() const
+    {
+      IndexSet is (size());
+
+      const std::pair<types::global_dof_index,types::global_dof_index> x = local_range();
+      is.add_range (x.first, x.second);
+
+      return is;
+    }
+
+
 
     template <typename Number>
     inline
index 49ad8674add393d83cd07552b9180890ba3325da..4222be10e12efe1abbe92b6a50bf4fec5aac7286 100644 (file)
@@ -392,6 +392,22 @@ namespace PETScWrappers
      */
     bool in_local_range (const unsigned int index) const;
 
+    /**
+     * Return an index set that describes which elements of this vector
+     * are owned by the current processor. Note that this index set does
+     * not include elements this vector may store locally as ghost
+     * elements but that are in fact owned by another processor.
+     * As a consequence, the index sets returned on different
+     * processors if this is a distributed vector will form disjoint
+     * sets that add up to the complete index set.
+     * Obviously, if a vector is created on only one processor, then
+     * the result would satisfy
+     * @code
+     *   vec.locally_owned_elements() == complete_index_set (vec.size())
+     * @endcode
+     */
+    IndexSet locally_owned_elements () const;
+
     /**
      * Return if the vector contains ghost
      * elements.
@@ -1127,6 +1143,21 @@ namespace PETScWrappers
             (index < static_cast<unsigned int>(end)));
   }
 
+
+  inline
+  IndexSet
+  VectorBase::locally_owned_elements() const
+  {
+    IndexSet is (size());
+
+    // PETSc only allows for contiguous local ranges, so this is simple
+    const std::pair<unsigned int, unsigned int> x = local_range();
+    is.add_range (x.first, x.second);
+    return is;
+  }
+
+
+
   inline
   bool
   VectorBase::has_ghost_elements() const
index bd180641b9b28151eba034eb3e19edc222131d0a..05d5f2411716271c932b32204cc68206ba133916 100644 (file)
@@ -521,6 +521,22 @@ namespace TrilinosWrappers
      */
     bool in_local_range (const unsigned int index) const;
 
+    /**
+     * Return an index set that describes which elements of this vector
+     * are owned by the current processor. Note that this index set does
+     * not include elements this vector may store locally as ghost
+     * elements but that are in fact owned by another processor.
+     * As a consequence, the index sets returned on different
+     * processors if this is a distributed vector will form disjoint
+     * sets that add up to the complete index set.
+     * Obviously, if a vector is created on only one processor, then
+     * the result would satisfy
+     * @code
+     *   vec.locally_owned_elements() == complete_index_set (vec.size())
+     * @endcode
+     */
+    IndexSet locally_owned_elements () const;
+
     /**
      * Return if the vector contains ghost
      * elements. This answer is true if there
@@ -1226,6 +1242,23 @@ namespace TrilinosWrappers
 
 
 
+  inline
+  IndexSet
+  VectorBase::locally_owned_elements() const
+  {
+    IndexSet is (size());
+
+    // TODO: Trilinos does allow non-contiguous ranges for locally
+    // owned elements. if that is the case for the current
+    // vector, local_range() will throw an exception. Fix this.
+    const std::pair<unsigned int, unsigned int> x = local_range();
+    is.add_range (x.first, x.second);
+
+    return is;
+  }
+
+
+
   inline
   bool
   VectorBase::has_ghost_elements() const
index 33f3a5ef47cc32b1cc0c16ca6a494143145d29ed..6c7fe3e9ced2288f81db3ba729e7fbd44b88ec1d 100644 (file)
@@ -17,6 +17,7 @@
 #include <deal.II/base/logstream.h>
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/subscriptor.h>
+#include <deal.II/base/index_set.h>
 #include <boost/serialization/array.hpp>
 #include <boost/serialization/split_member.hpp>
 
@@ -600,6 +601,26 @@ public:
    */
   bool in_local_range (const types::global_dof_index global_index) const;
 
+  /**
+   * Return an index set that describes which elements of this vector
+   * are owned by the current processor. Note that this index set does
+   * not include elements this vector may store locally as ghost
+   * elements but that are in fact owned by another processor.
+   * As a consequence, the index sets returned on different
+   * processors if this is a distributed vector will form disjoint
+   * sets that add up to the complete index set.
+   * Obviously, if a vector is created on only one processor, then
+   * the result would satisfy
+   * @code
+   *   vec.locally_owned_elements() == complete_index_set (vec.size())
+   * @endcode
+   *
+   * Since the current data type does not support parallel data storage
+   * across different processors, the returned index set is the
+   * complete index set.
+   */
+  IndexSet locally_owned_elements () const;
+
   /**
    * Return dimension of the vector.
    */
index 82437f140cbe0a2c9fd4dcea9c5b3843b5b943e3..5cb3ba459aed03a972071a0a252bebf18e168248 100644 (file)
@@ -1434,6 +1434,15 @@ void Vector<Number>::block_read (std::istream &in)
 
 
 
+template <typename Number>
+IndexSet
+Vector<Number>::locally_owned_elements() const
+{
+  return complete_index_set(size());
+}
+
+
+
 template <typename Number>
 std::size_t
 Vector<Number>::memory_consumption () 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.