]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add code to serialize PETSc vectors. 15945/head
authorDavid Wells <drwells@email.unc.edu>
Wed, 4 Oct 2023 15:28:55 +0000 (11:28 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 4 Oct 2023 15:28:55 +0000 (11:28 -0400)
include/deal.II/lac/petsc_vector_base.h
tests/petsc/serialize_vector_01.cc [new file with mode: 0644]
tests/petsc/serialize_vector_01.mpirun=2.output [new file with mode: 0644]

index a2dec0f63e2034600f0f7c9ee8b74ce3578d1164..d936d92aa5dcbb213170955386082c6b3a03d70e 100644 (file)
@@ -28,6 +28,8 @@
 #  include <deal.II/lac/vector.h>
 #  include <deal.II/lac/vector_operation.h>
 
+#  include <boost/serialization/split_member.hpp>
+
 #  include <petscvec.h>
 
 #  include <utility>
@@ -724,6 +726,41 @@ namespace PETScWrappers
           const bool         scientific = true,
           const bool         across     = true) const;
 
+    /**
+     * Write the data of this object to a stream for the purpose of
+     * serialization using the [BOOST serialization
+     * library](https://www.boost.org/doc/libs/1_74_0/libs/serialization/doc/index.html).
+     *
+     * @note Each processor only serializes its own locally owned values.
+     */
+    template <class Archive>
+    void
+    save(Archive &ar, const unsigned int version) const;
+
+    /**
+     * Read the data of this object from a stream for the purpose of
+     * serialization using the [BOOST serialization
+     * library](https://www.boost.org/doc/libs/1_74_0/libs/serialization/doc/index.html).
+     */
+    template <class Archive>
+    void
+    load(Archive &ar, const unsigned int version);
+
+#  ifdef DOXYGEN
+    /**
+     * Write and read the data of this object from a stream for the purpose
+     * of serialization using the [BOOST serialization
+     * library](https://www.boost.org/doc/libs/1_74_0/libs/serialization/doc/index.html).
+     */
+    template <class Archive>
+    void
+    serialize(Archive &archive, const unsigned int version);
+#  else
+    // This macro defines the serialize() method that is compatible with
+    // the templated save() and load() method that have been implemented.
+    BOOST_SERIALIZATION_SPLIT_MEMBER()
+#  endif
+
     /**
      * Swap the contents of this vector and the other vector @p v. One could
      * do this operation with a temporary variable and copying over the data
@@ -1271,6 +1308,65 @@ namespace PETScWrappers
       }
   }
 
+  template <class Archive>
+  inline void
+  VectorBase::save(Archive &ar, const unsigned int) const
+  {
+    // forward to serialization function in the base class.
+    ar &static_cast<const Subscriptor &>(*this);
+    ar &size();
+    ar &local_range();
+
+    const PetscScalar *array = nullptr;
+    int                ierr  = VecGetArrayRead(*this, &array);
+    AssertThrow(ierr == 0, ExcPETScError(ierr));
+
+    boost::serialization::array_wrapper<const double> wrapper(
+      array, locally_owned_size());
+    ar &wrapper;
+
+    ierr = VecRestoreArrayRead(*this, &array);
+    AssertThrow(ierr == 0, ExcPETScError(ierr));
+  }
+
+
+
+  template <class Archive>
+  inline void
+  VectorBase::load(Archive &ar, const unsigned int)
+  {
+    ar &static_cast<Subscriptor &>(*this);
+
+    size_type                       size = 0;
+    std::pair<size_type, size_type> local_range;
+
+    ar &size;
+    Assert(size == this->size(),
+           ExcMessage("The serialized value of size (" + std::to_string(size) +
+                      ") does not match the current size (" +
+                      std::to_string(this->size()) + ")"));
+    (void)size;
+    ar &local_range;
+    Assert(local_range == this->local_range(),
+           ExcMessage("The serialized value of local_range (" +
+                      std::to_string(local_range.first) + ", " +
+                      std::to_string(local_range.second) +
+                      ") does not match the current local_range (" +
+                      std::to_string(this->local_range().first) + ", " +
+                      std::to_string(this->local_range().second) + ")"));
+    (void)local_range;
+
+    PetscScalar *array = nullptr;
+    int          ierr  = VecGetArray(petsc_vector(), &array);
+    AssertThrow(ierr == 0, ExcPETScError(ierr));
+
+    boost::serialization::array_wrapper<double> velocity_wrapper(
+      array, locally_owned_size());
+    ar &velocity_wrapper;
+
+    ierr = VecRestoreArray(petsc_vector(), &array);
+    AssertThrow(ierr == 0, ExcPETScError(ierr));
+  }
 #  endif // DOXYGEN
 } // namespace PETScWrappers
 
diff --git a/tests/petsc/serialize_vector_01.cc b/tests/petsc/serialize_vector_01.cc
new file mode 100644 (file)
index 0000000..5a51af4
--- /dev/null
@@ -0,0 +1,63 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Check that we can correctly serialize a vector via boost.
+
+#include <deal.II/base/index_set.h>
+#include <deal.II/base/mpi.h>
+
+#include <deal.II/lac/petsc_vector.h>
+
+#include "../tests.h"
+
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_context(argc, argv, 1);
+  MPILogInitAll                    mpi_log;
+
+  const auto rank            = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  const auto n_mpi_processes = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+
+  constexpr types::global_dof_index dofs_per_process = 5;
+  IndexSet owned_indices(dofs_per_process * n_mpi_processes);
+  owned_indices.add_range(dofs_per_process * rank,
+                          dofs_per_process * (rank + 1));
+  owned_indices.compress();
+
+  owned_indices.print(deallog);
+
+  PETScWrappers::MPI::Vector petsc_vector(owned_indices, MPI_COMM_WORLD);
+  for (const auto dof : owned_indices)
+    petsc_vector(dof) = double(dof);
+
+  petsc_vector.compress(VectorOperation::insert);
+
+  std::stringstream               stream;
+  boost::archive::binary_oarchive out_archive(stream);
+  out_archive << petsc_vector;
+
+  PETScWrappers::MPI::Vector petsc_copy(owned_indices, MPI_COMM_WORLD);
+
+  boost::archive::binary_iarchive in_archive(stream);
+  in_archive >> petsc_copy;
+
+  petsc_copy.print(deallog.get_file_stream());
+
+  AssertThrow(petsc_copy == petsc_vector, ExcInternalError());
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/petsc/serialize_vector_01.mpirun=2.output b/tests/petsc/serialize_vector_01.mpirun=2.output
new file mode 100644 (file)
index 0000000..3f56274
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL:0::{[0,4]}
+[Proc0 0-4] 0.000e+00 1.000e+00 2.000e+00 3.000e+00 4.000e+00 
+DEAL:0::OK
+
+DEAL:1::{[5,9]}
+[Proc1 5-9] 5.000e+00 6.000e+00 7.000e+00 8.000e+00 9.000e+00 
+DEAL:1::OK
+

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.