]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix memory leak in Trilinos vector 2548/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 23 Apr 2016 18:28:48 +0000 (20:28 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sun, 24 Apr 2016 12:40:32 +0000 (14:40 +0200)
doc/news/changes.h
include/deal.II/lac/trilinos_vector.h
source/lac/trilinos_vector.cc
tests/lac/vector_reinit_01.cc [new file with mode: 0644]
tests/lac/vector_reinit_01.output [new file with mode: 0644]
tests/lac/vector_reinit_02.cc [new file with mode: 0644]
tests/lac/vector_reinit_02.with_trilinos=true.output [new file with mode: 0644]
tests/lac/vector_reinit_03.cc [new file with mode: 0644]
tests/lac/vector_reinit_03.with_petsc=true.output [new file with mode: 0644]

index 4878966edcb4747665c4cce1d61063b5b2a7c4d1..7ba12b45f34065c5a57d4f4c8ed24fed39bae1cd 100644 (file)
@@ -201,6 +201,15 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+ <li> Fixed: TrilinosWrappers::MPI::Vector and TrilinosWrappers::Vector could
+ access invalid memory in the reinit() method if the MPI communicator was
+ deleted before termination of the program. This usually happened when using
+ vectors from GrowingVectorMemory where a pool keeps vector alive. This has
+ been fixed.
+ <br>
+ (Martin Kronbichler, 2016/04/23)
+ </li>
+
  <li> Fixed: The methods TrilinosWrappers::SparseMatrix::(T)mmult previously
  produced invalid matrix sizes if the final matrix was non-square. This has
  been fixed.
index 2954a48f44ea6280703cf55eced38b2dad63d365..da8b33c621d7fb22f1b6b633d8efd04f386b42f8 100644 (file)
@@ -307,8 +307,12 @@ namespace TrilinosWrappers
        * dimension and the parallel distribution of the input vector, but does
        * not copy the elements in <tt>v</tt>. If <tt>omit_zeroing_entries</tt>
        * is not <tt>true</tt>, the elements in the vector are initialized with
-       * zero, otherwise the content will be left unchanged and the user has
-       * to set all elements.
+       * zero. If it is set to <tt>true</tt>, the vector entries are in an
+       * unspecified state and the user has to set all elements. In the
+       * current implementation, this method does not touch the vector entries
+       * in case the vector layout is unchanged from before, otherwise entries
+       * are set to zero.  Note that this behavior might change between
+       * releases without notification.
        *
        * This function has a third argument, <tt>allow_different_maps</tt>,
        * that allows for an exchange of data between two equal-sized vectors
@@ -567,7 +571,11 @@ namespace TrilinosWrappers
        * Reinit functionality. This function destroys the old vector content
        * and generates a new one based on the input partitioning.  The flag
        * <tt>omit_zeroing_entries</tt> determines whether the vector should be
-       * filled with zero (false) or left untouched (true).
+       * filled with zero (false). If the flag is set to <tt>true</tt>, the
+       * vector entries are in an unspecified state and the user has to set
+       * all elements. In the current implementation, this method still sets
+       * the entries to zero, but this might change between releases without
+       * notification.
        *
        *
        * Depending on whether the @p parallel_partitioning argument uniquely
@@ -796,7 +804,12 @@ namespace TrilinosWrappers
 
     /**
      * Reinit function that resizes the vector to the size specified by
-     * <tt>n</tt>.
+     * <tt>n</tt>. The flag <tt>omit_zeroing_entries</tt> specifies whether
+     * the vector entries should be initialized to zero (false). If it is set
+     * to <tt>true</tt>, the vector entries are in an unspecified state and
+     * the user has to set all elements. In the current implementation, this
+     * method still sets the entries to zero, but this might change between
+     * releases without notification.
      */
     void reinit (const size_type n,
                  const bool      omit_zeroing_entries = false);
@@ -808,7 +821,11 @@ namespace TrilinosWrappers
      * in the localized vector should be imported from a distributed vector
      * that has been initialized with the same communicator. The variable
      * <tt>omit_zeroing_entries</tt> determines whether the vector should be
-     * filled with zero or left untouched.
+     * filled with zero (false). If it is set to <tt>true</tt>, the vector
+     * entries are in an unspecified state and the user has to set all
+     * elements. In the current implementation, this method still sets the
+     * entries to zero, but this might change between releases without
+     * notification.
      *
      * Which element of the @p input_map argument are set is in fact ignored,
      * the only thing that matters is the size of the index space described by
@@ -824,7 +841,11 @@ namespace TrilinosWrappers
      * in the localized vector should be imported from a distributed vector
      * that has been initialized with the same communicator. The variable
      * <tt>omit_zeroing_entries</tt> determines whether the vector should be
-     * filled with zero (false) or left untouched (true).
+     * filled with zero (false). If it is set to <tt>true</tt>, the vector
+     * entries are in an unspecified state and the user has to set all
+     * elements. In the current implementation, this method still sets the
+     * entries to zero, but this might change between releases without
+     * notification.
      *
      * Which element of the @p input_map argument are set is in fact ignored,
      * the only thing that matters is the size of the index space described by
@@ -837,6 +858,15 @@ namespace TrilinosWrappers
     /**
      * Reinit function. Takes the information of a Vector and copies
      * everything to the calling vector, now also allowing different maps.
+     *
+     * If the optional flag <tt>omit_zeroing_entries</tt> is not
+     * <tt>true</tt>, the elements in the vector are initialized with zero. If
+     * it is set to <tt>true</tt>, the vector entries are in an unspecified
+     * state and the user has to set all elements. In the current
+     * implementation, this method does not touch the vector entries in case
+     * the vector layout is unchanged from before, otherwise entries are set
+     * to zero.  Note that this behavior might change between releases without
+     * notification.
      */
     void reinit (const VectorBase &V,
                  const bool        omit_zeroing_entries = false,
index b0cbba0a4042d33058b778e2a9b2ec685ca317f6..52479dfe80538c5bfb7cfbd2e93d1f40c683002f 100644 (file)
@@ -199,18 +199,11 @@ namespace TrilinosWrappers
 
     void
     Vector::reinit (const Epetra_Map &input_map,
-                    const bool        omit_zeroing_entries)
+                    const bool      /*omit_zeroing_entries*/)
     {
       nonlocal_vector.reset();
 
-      if (vector->Map().SameAs(input_map)==false)
-        vector.reset (new Epetra_FEVector(input_map));
-      else if (omit_zeroing_entries == false)
-        {
-          const int ierr = vector->PutScalar(0.);
-          (void)ierr;
-          Assert (ierr == 0, ExcTrilinosError(ierr));
-        }
+      vector.reset (new Epetra_FEVector(input_map));
 
       has_ghosts = vector->Map().UniqueGIDs()==false;
       last_action = Zero;
@@ -244,7 +237,19 @@ namespace TrilinosWrappers
       // with the map in v, and generate the vector.
       if (allow_different_maps == false)
         {
-          if (vector->Map().SameAs(v.vector->Map()) == false)
+          // check equality for MPI communicators: We can only choose the fast
+          // version in case the underlying Epetra_MpiComm object is the same,
+          // otherwise we might access an MPI_Comm object that has been
+          // deleted
+#ifdef DEAL_II_WITH_MPI
+          const Epetra_MpiComm *my_comm = dynamic_cast<const Epetra_MpiComm *>(&vector->Comm());
+          const Epetra_MpiComm *v_comm = dynamic_cast<const Epetra_MpiComm *>(&v.vector->Comm());
+          const bool same_communicators = my_comm != NULL && v_comm != NULL &&
+                                          my_comm->DataPtr() == v_comm->DataPtr();
+#else
+          const bool same_communicators = true;
+#endif
+          if (!same_communicators || vector->Map().SameAs(v.vector->Map()) == false)
             {
               vector.reset (new Epetra_FEVector(v.vector->Map()));
               has_ghosts = v.has_ghosts;
@@ -252,11 +257,8 @@ namespace TrilinosWrappers
             }
           else if (omit_zeroing_entries == false)
             {
-              // old and new vectors
-              // have exactly the
-              // same map, i.e. size
-              // and parallel
-              // distribution
+              // old and new vectors have exactly the same map, i.e. size and
+              // parallel distribution
               int ierr;
               ierr = vector->GlobalAssemble (last_action);
               (void)ierr;
@@ -397,11 +399,22 @@ namespace TrilinosWrappers
     Vector &
     Vector::operator = (const Vector &v)
     {
+      // check equality for MPI communicators to avoid accessing a possibly
+      // invalid MPI_Comm object
+#ifdef DEAL_II_WITH_MPI
+      const Epetra_MpiComm *my_comm = dynamic_cast<const Epetra_MpiComm *>(&vector->Comm());
+      const Epetra_MpiComm *v_comm = dynamic_cast<const Epetra_MpiComm *>(&v.vector->Comm());
+      const bool same_communicators = my_comm != NULL && v_comm != NULL &&
+                                      my_comm->DataPtr() == v_comm->DataPtr();
+#else
+      const bool same_communicators = true;
+#endif
+
       // distinguish three cases. First case: both vectors have the same
       // layout (just need to copy the local data, not reset the memory and
       // the underlying Epetra_Map). The third case means that we have to
       // rebuild the calling vector.
-      if (vector->Map().SameAs(v.vector->Map()))
+      if (same_communicators && vector->Map().SameAs(v.vector->Map()))
         {
           *vector = *v.vector;
           if (v.nonlocal_vector.get() != 0)
@@ -562,24 +575,11 @@ namespace TrilinosWrappers
 
   void
   Vector::reinit (const size_type n,
-                  const bool      omit_zeroing_entries)
+                  const bool    /*omit_zeroing_entries*/)
   {
-    if (size() != n)
-      {
-        Epetra_LocalMap map ((TrilinosWrappers::types::int_type)n, 0,
-                             Utilities::Trilinos::comm_self());
-        vector.reset (new Epetra_FEVector (map));
-      }
-    else if (omit_zeroing_entries == false)
-      {
-        int ierr;
-        ierr = vector->GlobalAssemble(last_action);
-        (void)ierr;
-        Assert (ierr == 0, ExcTrilinosError(ierr));
-
-        ierr = vector->PutScalar(0.0);
-        Assert (ierr == 0, ExcTrilinosError(ierr));
-      }
+    Epetra_LocalMap map ((TrilinosWrappers::types::int_type)n, 0,
+                         Utilities::Trilinos::comm_self());
+    vector.reset (new Epetra_FEVector (map));
 
     last_action = Zero;
   }
@@ -588,25 +588,12 @@ namespace TrilinosWrappers
 
   void
   Vector::reinit (const Epetra_Map &input_map,
-                  const bool        omit_zeroing_entries)
+                  const bool     /*omit_zeroing_entries*/)
   {
-    if (n_global_elements(vector->Map()) != n_global_elements(input_map))
-      {
-        Epetra_LocalMap map (n_global_elements(input_map),
-                             input_map.IndexBase(),
-                             input_map.Comm());
-        vector.reset (new Epetra_FEVector (map));
-      }
-    else if (omit_zeroing_entries == false)
-      {
-        int ierr;
-        ierr = vector->GlobalAssemble(last_action);
-        (void)ierr;
-        Assert (ierr == 0, ExcTrilinosError(ierr));
-
-        ierr = vector->PutScalar(0.0);
-        Assert (ierr == 0, ExcTrilinosError(ierr));
-      }
+    Epetra_LocalMap map (n_global_elements(input_map),
+                         input_map.IndexBase(),
+                         input_map.Comm());
+    vector.reset (new Epetra_FEVector (map));
 
     last_action = Zero;
   }
@@ -616,31 +603,17 @@ namespace TrilinosWrappers
   void
   Vector::reinit (const IndexSet &partitioning,
                   const MPI_Comm &communicator,
-                  const bool      omit_zeroing_entries)
+                  const bool    /*omit_zeroing_entries*/)
   {
-    if (n_global_elements(vector->Map()) !=
-        static_cast<TrilinosWrappers::types::int_type>(partitioning.size()))
-      {
-        Epetra_LocalMap map (static_cast<TrilinosWrappers::types::int_type>(partitioning.size()),
-                             0,
+    Epetra_LocalMap map (static_cast<TrilinosWrappers::types::int_type>(partitioning.size()),
+                         0,
 #ifdef DEAL_II_WITH_MPI
-                             Epetra_MpiComm(communicator));
+                         Epetra_MpiComm(communicator));
 #else
-                             Epetra_SerialComm());
-        (void)communicator;
+                         Epetra_SerialComm());
+    (void)communicator;
 #endif
-        vector.reset (new Epetra_FEVector(map));
-      }
-    else if (omit_zeroing_entries == false)
-      {
-        int ierr;
-        ierr = vector->GlobalAssemble(last_action);
-        (void)ierr;
-        Assert (ierr == 0, ExcTrilinosError(ierr));
-
-        ierr = vector->PutScalar(0.0);
-        Assert (ierr == 0, ExcTrilinosError(ierr));
-      }
+    vector.reset (new Epetra_FEVector(map));
 
     last_action = Zero;
   }
@@ -659,17 +632,25 @@ namespace TrilinosWrappers
     // the vector, initialize our
     // map with the map in v, and
     // generate the vector.
-    (void)omit_zeroing_entries;
     if (allow_different_maps == false)
       {
-        if (local_range() != v.local_range())
+        // check equality for MPI communicators
+#ifdef DEAL_II_WITH_MPI
+        const Epetra_MpiComm *my_comm = dynamic_cast<const Epetra_MpiComm *>(&vector->Comm());
+        const Epetra_MpiComm *v_comm = dynamic_cast<const Epetra_MpiComm *>(&v.vector->Comm());
+        const bool same_communicators = my_comm != NULL && v_comm != NULL &&
+                                        my_comm->DataPtr() == v_comm->DataPtr();
+#else
+        const bool same_communicators = true;
+#endif
+        if (!same_communicators || local_range() != v.local_range())
           {
             Epetra_LocalMap map (global_length(*(v.vector)),
                                  v.vector->Map().IndexBase(),
                                  v.vector->Comm());
             vector.reset (new Epetra_FEVector(map));
           }
-        else
+        else if (omit_zeroing_entries)
           {
             int ierr;
             Assert (vector->Map().SameAs(v.vector->Map()) == true,
diff --git a/tests/lac/vector_reinit_01.cc b/tests/lac/vector_reinit_01.cc
new file mode 100644 (file)
index 0000000..cd4e48a
--- /dev/null
@@ -0,0 +1,107 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// check that parallel::distributed::Vector::reinit does not carry over any
+// state that can lead to invalid memory access. In this test, the MPI
+// communicator is deleted.
+
+
+#include "../tests.h"
+#include <deal.II/base/mpi.h>
+#include <deal.II/lac/parallel_vector.h>
+#include <deal.II/lac/vector_memory.h>
+
+template <typename VectorType>
+void do_test()
+{
+  IndexSet set(5);
+  set.add_range(0,5);
+
+  VectorType v1, v2;
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v1.reinit(set, communicator);
+    deallog << "reinit: " << v1.size() << " ";
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v2.reinit(set, communicator);
+    v1.reinit(v2);
+    deallog << v1.size() << " ";
+    v1.reinit(v2);
+    deallog << v1.size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v2.reinit(set, communicator);
+    v1 = v2;
+    deallog << "assign " << v1.size() << " ";
+    v1 = v2;
+    deallog << v1.size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    GrowingVectorMemory<VectorType> memory;
+    typename VectorMemory<VectorType>::Pointer v3(memory);
+    v1.reinit(set, communicator);
+    v3->reinit(v1);
+    deallog << "reinit pool " << v1.size() << " " << v3->size() << " ";
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    GrowingVectorMemory<VectorType> memory;
+    typename VectorMemory<VectorType>::Pointer v3(memory);
+    v1.reinit(set, communicator);
+    v3->reinit(v1);
+    deallog << "reinit pool " << v3->size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+}
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  mpi_initlog();
+
+  do_test<parallel::distributed::Vector<double> >();
+}
diff --git a/tests/lac/vector_reinit_01.output b/tests/lac/vector_reinit_01.output
new file mode 100644 (file)
index 0000000..747b8af
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL::reinit: 5 5 5
+DEAL::assign 5 5
+DEAL::reinit pool 5 5 reinit pool 5
diff --git a/tests/lac/vector_reinit_02.cc b/tests/lac/vector_reinit_02.cc
new file mode 100644 (file)
index 0000000..9977c26
--- /dev/null
@@ -0,0 +1,108 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// check that TrilinosWrappers(::MPI)::Vector::reinit does not carry over any
+// state that can lead to invalid memory access. In this test, the MPI
+// communicator is deleted.
+
+
+#include "../tests.h"
+#include <deal.II/base/mpi.h>
+#include <deal.II/lac/trilinos_vector.h>
+#include <deal.II/lac/vector_memory.h>
+
+template <typename VectorType>
+void do_test()
+{
+  IndexSet set(5);
+  set.add_range(0,5);
+
+  VectorType v1, v2;
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v1.reinit(set, communicator);
+    deallog << "reinit: " << v1.size() << " ";
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v2.reinit(set, communicator);
+    v1.reinit(v2);
+    deallog << v1.size() << " ";
+    v1.reinit(v2);
+    deallog << v1.size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v2.reinit(set, communicator);
+    v1 = v2;
+    deallog << "assign " << v1.size() << " ";
+    v1 = v2;
+    deallog << v1.size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    GrowingVectorMemory<VectorType> memory;
+    typename VectorMemory<VectorType>::Pointer v3(memory);
+    v1.reinit(set, communicator);
+    v3->reinit(v1);
+    deallog << "reinit pool " << v1.size() << " " << v3->size() << " ";
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    GrowingVectorMemory<VectorType> memory;
+    typename VectorMemory<VectorType>::Pointer v3(memory);
+    v1.reinit(set, communicator);
+    v3->reinit(v1);
+    deallog << "reinit pool " << v3->size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+}
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  mpi_initlog();
+
+  do_test<TrilinosWrappers::Vector>();
+  do_test<TrilinosWrappers::MPI::Vector>();
+}
diff --git a/tests/lac/vector_reinit_02.with_trilinos=true.output b/tests/lac/vector_reinit_02.with_trilinos=true.output
new file mode 100644 (file)
index 0000000..02a81b0
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL::reinit: 5 5 5
+DEAL::assign 5 5
+DEAL::reinit pool 5 5 reinit pool 5
+DEAL::reinit: 5 5 5
+DEAL::assign 5 5
+DEAL::reinit pool 5 5 reinit pool 5
diff --git a/tests/lac/vector_reinit_03.cc b/tests/lac/vector_reinit_03.cc
new file mode 100644 (file)
index 0000000..8112357
--- /dev/null
@@ -0,0 +1,108 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// check that PETScWrappers::MPI::Vector::reinit does not carry over any
+// state that can lead to invalid memory access. In this test, the MPI
+// communicator is deleted.
+
+
+#include "../tests.h"
+#include <deal.II/base/mpi.h>
+#include <deal.II/lac/petsc_vector.h>
+#include <deal.II/lac/petsc_parallel_vector.h>
+#include <deal.II/lac/vector_memory.h>
+
+template <typename VectorType>
+void do_test()
+{
+  IndexSet set(5);
+  set.add_range(0,5);
+
+  VectorType v1, v2;
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v1.reinit(set, communicator);
+    deallog << "reinit: " << v1.size() << " ";
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v2.reinit(set, communicator);
+    v1.reinit(v2);
+    deallog << v1.size() << " ";
+    v1.reinit(v2);
+    deallog << v1.size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    v2.reinit(set, communicator);
+    v1 = v2;
+    deallog << "assign " << v1.size() << " ";
+    v1 = v2;
+    deallog << v1.size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    GrowingVectorMemory<VectorType> memory;
+    typename VectorMemory<VectorType>::Pointer v3(memory);
+    v1.reinit(set, communicator);
+    v3->reinit(v1);
+    deallog << "reinit pool " << v1.size() << " " << v3->size() << " ";
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+  {
+    MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD);
+    GrowingVectorMemory<VectorType> memory;
+    typename VectorMemory<VectorType>::Pointer v3(memory);
+    v1.reinit(set, communicator);
+    v3->reinit(v1);
+    deallog << "reinit pool " << v3->size() << std::endl;
+
+#ifdef DEAL_II_WITH_MPI
+    MPI_Comm_free (&communicator);
+#endif
+  }
+
+}
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  mpi_initlog();
+
+  do_test<PETScWrappers::MPI::Vector>();
+}
diff --git a/tests/lac/vector_reinit_03.with_petsc=true.output b/tests/lac/vector_reinit_03.with_petsc=true.output
new file mode 100644 (file)
index 0000000..747b8af
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL::reinit: 5 5 5
+DEAL::assign 5 5
+DEAL::reinit pool 5 5 reinit pool 5

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.