]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clean up SUNDIALS contexts and communicators.
authorDavid Wells <drwells@email.unc.edu>
Fri, 3 Jun 2022 19:07:22 +0000 (15:07 -0400)
committerDavid Wells <drwells@email.unc.edu>
Sun, 5 Jun 2022 19:35:35 +0000 (15:35 -0400)
SUNDIALS will duplicate communicators for us, as needed. For example, in
SUNProfiler_Create, SUNDIALS implements

    #if SUNDIALS_MPI_ENABLED
    profiler->comm = NULL;
    if (comm != NULL)
    {
        profiler->comm = malloc(sizeof(MPI_Comm));
        MPI_Comm_dup(*((MPI_Comm*) comm), (MPI_Comm*) profiler->comm);
    }
    #else
    profiler->comm = comm;
    #endif

It is preferable to let SUNDIALS duplicate the communicator so that we
don't call MPI functions when we have only serial data structures (and
may not have initialized MPI).

include/deal.II/sundials/arkode.h
include/deal.II/sundials/ida.h
include/deal.II/sundials/kinsol.h
source/sundials/arkode.cc
source/sundials/ida.cc
source/sundials/kinsol.cc

index 1d8fc96e86c40853fee460654f2fcc1971eaebd6..4eaa03887ab0d117b404e09e7adc1583834f7cef 100644 (file)
@@ -473,18 +473,24 @@ namespace SUNDIALS
     };
 
     /**
-     * Constructor. It is possible to fine tune the SUNDIALS ARKode solver by
-     * passing an AdditionalData() object that sets all of the solver
-     * parameters.
+     * Constructor, with class parameters set by the AdditionalData object.
      *
-     * The MPI communicator is simply ignored in the serial case.
+     * @param data ARKode configuration data
      *
+     * @note With SUNDIALS 6 and later this constructor sets up logging
+     * objects to only work on the present processor (i.e., results are only
+     * communicated over MPI_COMM_SELF).
+     */
+    ARKode(const AdditionalData &data = AdditionalData());
+
+    /**
+     * Constructor.
      *
      * @param data ARKode configuration data
-     * @param mpi_comm MPI communicator
+     * @param mpi_comm MPI Communicator over which logging operations are
+     * computed. Only used in SUNDIALS 6 and newer.
      */
-    ARKode(const AdditionalData &data     = AdditionalData(),
-           const MPI_Comm &      mpi_comm = MPI_COMM_WORLD);
+    ARKode(const AdditionalData &data, const MPI_Comm &mpi_comm);
 
     /**
      * Destructor.
@@ -1064,9 +1070,8 @@ namespace SUNDIALS
 #  endif
 
     /**
-     * MPI communicator. SUNDIALS solver runs happily in
-     * parallel. Note that if the library is compiled without MPI
-     * support, MPI_Comm is aliased as int.
+     * MPI communicator. Only used for SUNDIALS' logging routines - the actual
+     * solve routines will use the communicator provided by the vector class.
      */
     MPI_Comm mpi_communicator;
 
index 7e374851da6348226fe3d195f4630e53274c2ea2..4ae29e54320f2b56abdb2ab18e8f56ffa6dca181 100644 (file)
@@ -590,13 +590,22 @@ namespace SUNDIALS
      * choose how these are made consistent, using the same three options that
      * you used for the initial conditions in `reset_type`.
      *
-     * The MPI communicator is simply ignored in the serial case.
-     *
      * @param data IDA configuration data
-     * @param mpi_comm MPI communicator
+     *
+     * @note With SUNDIALS 6 and later this constructor sets up logging
+     * objects to only work on the present processor (i.e., results are only
+     * communicated over MPI_COMM_SELF).
+     */
+    IDA(const AdditionalData &data = AdditionalData());
+
+    /**
+     * Constructor.
+     *
+     * @param data IDA configuration data. Same as the other constructor.
+     * @param mpi_comm MPI Communicator over which logging operations are
+     * computed. Only used in SUNDIALS 6 and newer.
      */
-    IDA(const AdditionalData &data     = AdditionalData(),
-        const MPI_Comm &      mpi_comm = MPI_COMM_WORLD);
+    IDA(const AdditionalData &data, const MPI_Comm &mpi_comm);
 
     /**
      * Destructor.
@@ -885,9 +894,8 @@ namespace SUNDIALS
 #  endif
 
     /**
-     * MPI communicator. SUNDIALS solver runs happily in
-     * parallel. Note that if the library is compiled without MPI
-     * support, MPI_Comm is aliased as int.
+     * MPI communicator. Only used for SUNDIALS' logging routines - the actual
+     * solve routines will use the communicator provided by the vector class.
      */
     MPI_Comm mpi_communicator;
 
index 9bcdfda83c1b379e6eaceea0b86d5d65b4e95555..364c82e5bf0254df36a89ad50d00fd2a9a8d9a8f 100644 (file)
@@ -373,15 +373,24 @@ namespace SUNDIALS
     };
 
     /**
-     * Constructor. It is possible to fine tune the SUNDIALS KINSOL solver by
-     * passing an AdditionalData() object that sets all of the solver
-     * parameters.
+     * Constructor, with class parameters set by the AdditionalData object.
      *
      * @param data KINSOL configuration data
-     * @param mpi_comm MPI communicator
+     *
+     * @note With SUNDIALS 6 and later this constructor sets up logging
+     * objects to only work on the present processor (i.e., results are only
+     * communicated over MPI_COMM_SELF).
+     */
+    KINSOL(const AdditionalData &data = AdditionalData());
+
+    /**
+     * Constructor.
+     *
+     * @param data KINSOL configuration data
+     * @param mpi_comm MPI Communicator over which logging operations are
+     * computed. Only used in SUNDIALS 6 and newer.
      */
-    KINSOL(const AdditionalData &data     = AdditionalData(),
-           const MPI_Comm &      mpi_comm = MPI_COMM_WORLD);
+    KINSOL(const AdditionalData &data, const MPI_Comm &mpi_comm);
 
     /**
      * Destructor.
@@ -691,6 +700,13 @@ namespace SUNDIALS
      */
     void *kinsol_mem;
 
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    /**
+     * A context object associated with the KINSOL solver.
+     */
+    SUNContext kinsol_ctx;
+#  endif
+
     /**
      * KINSOL solution vector.
      */
@@ -706,13 +722,6 @@ namespace SUNDIALS
      */
     N_Vector f_scale;
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-    /**
-     * A context object associated with the KINSOL solver.
-     */
-    SUNContext kinsol_ctx;
-#  endif
-
     /**
      * Memory pool of vectors.
      */
index 6c48c6adec08a0d148ddc7883e52b8738f7cdb21..ff4386c5bca4d2d324144cc637a82e4f77bdd30d 100644 (file)
@@ -249,17 +249,38 @@ namespace SUNDIALS
   } // namespace
 
 
+  template <typename VectorType>
+  ARKode<VectorType>::ARKode(const AdditionalData &data)
+    : ARKode(data, MPI_COMM_SELF)
+  {}
+
+
   template <typename VectorType>
   ARKode<VectorType>::ARKode(const AdditionalData &data,
                              const MPI_Comm &      mpi_comm)
     : data(data)
     , arkode_mem(nullptr)
-    , mpi_communicator(is_serial_vector<VectorType>::value ?
-                         MPI_COMM_SELF :
-                         Utilities::MPI::duplicate_communicator(mpi_comm))
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    , arkode_ctx(nullptr)
+#  endif
+    , mpi_communicator(mpi_comm)
     , last_end_time(data.initial_time)
   {
     set_functions_to_trigger_an_assert();
+
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    // SUNDIALS will always duplicate communicators if we provide them. This
+    // can cause problems if SUNDIALS is configured with MPI and we pass along
+    // MPI_COMM_SELF in a serial application as MPI won't be
+    // initialized. Hence, work around that by just not providing a
+    // communicator in that case.
+    const int status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
+                                                            &mpi_communicator,
+                        &arkode_ctx);
+    (void)status;
+    AssertARKode(status);
+#  endif
   }
 
 
@@ -267,20 +288,12 @@ namespace SUNDIALS
   template <typename VectorType>
   ARKode<VectorType>::~ARKode()
   {
-    if (arkode_mem)
-      {
-        ARKStepFree(&arkode_mem);
+    ARKStepFree(&arkode_mem);
 
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-        const int status = SUNContext_Free(&arkode_ctx);
-        (void)status;
-        AssertARKode(status);
-#  endif
-      }
-
-#  ifdef DEAL_II_WITH_MPI
-    if (is_serial_vector<VectorType>::value == false)
-      Utilities::MPI::free_communicator(mpi_communicator);
+    const int status = SUNContext_Free(&arkode_ctx);
+    (void)status;
+    AssertARKode(status);
 #  endif
   }
 
@@ -387,12 +400,13 @@ namespace SUNDIALS
     (void)status;
 
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-    if (arkode_ctx)
-      {
-        status = SUNContext_Free(&arkode_ctx);
-        AssertARKode(status);
-      }
-    status = SUNContext_Create(&mpi_communicator, &arkode_ctx);
+    status = SUNContext_Free(&arkode_ctx);
+    AssertARKode(status);
+    // Same comment applies as in class constructor:
+    status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
+                                                            &mpi_communicator,
+                        &arkode_ctx);
     AssertARKode(status);
 #  endif
 
index c07aa31977d36889f1816c04adc1bf6bcdd06ea2..0ebe08554042e894678ab4470991f14ea94a0e06 100644 (file)
@@ -119,14 +119,35 @@ namespace SUNDIALS
 
 
 
+  template <typename VectorType>
+  IDA<VectorType>::IDA(const AdditionalData &data)
+    : IDA(data, MPI_COMM_SELF)
+  {}
+
+
+
   template <typename VectorType>
   IDA<VectorType>::IDA(const AdditionalData &data, const MPI_Comm &mpi_comm)
     : data(data)
     , ida_mem(nullptr)
-    , mpi_communicator(is_serial_vector<VectorType>::value ?
-                         MPI_COMM_SELF :
-                         Utilities::MPI::duplicate_communicator(mpi_comm))
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    , ida_ctx(nullptr)
+#  endif
+    , mpi_communicator(mpi_comm)
   {
+    // SUNDIALS will always duplicate communicators if we provide them. This
+    // can cause problems if SUNDIALS is configured with MPI and we pass along
+    // MPI_COMM_SELF in a serial application as MPI won't be
+    // initialized. Hence, work around that by just not providing a
+    // communicator in that case.
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    const int status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
+                                                            &mpi_communicator,
+                        &ida_ctx);
+    (void)status;
+    AssertIDA(status);
+#  endif
     set_functions_to_trigger_an_assert();
   }
 
@@ -135,20 +156,11 @@ namespace SUNDIALS
   template <typename VectorType>
   IDA<VectorType>::~IDA()
   {
-    if (ida_mem)
-      {
-        IDAFree(&ida_mem);
-
-#  if !DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
-        const int status = SUNContext_Free(&ida_ctx);
-        (void)status;
-        AssertIDA(status);
-#  endif
-      }
-
-#  ifdef DEAL_II_WITH_MPI
-    if (is_serial_vector<VectorType>::value == false)
-      Utilities::MPI::free_communicator(mpi_communicator);
+    IDAFree(&ida_mem);
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    const int status = SUNContext_Free(&ida_ctx);
+    (void)status;
+    AssertIDA(status);
 #  endif
   }
 
@@ -220,27 +232,30 @@ namespace SUNDIALS
   {
     bool first_step = (current_time == data.initial_time);
 
+    int status;
+    (void)status;
+
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    status = SUNContext_Free(&ida_ctx);
+    AssertIDA(status);
+
+    // Same comment applies as in class constructor:
+    status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
+                                                            &mpi_communicator,
+                        &ida_ctx);
+    AssertIDA(status);
+#  endif
+
     if (ida_mem)
       {
         IDAFree(&ida_mem);
-
-#  if !DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
-        const int status = SUNContext_Free(&ida_ctx);
-        (void)status;
-        AssertIDA(status);
-#  endif
+        // Initialization is version-dependent: do that in a moment
       }
 
-    int status;
-    (void)status;
-
-
 #  if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
     ida_mem = IDACreate();
 #  else
-    status = SUNContext_Create(&mpi_communicator, &ida_ctx);
-    AssertIDA(status);
-
     ida_mem = IDACreate(ida_ctx);
 #  endif
 
index 921755853f1b98f5426c5012705cab67eea3b6bf..e5c4daa94f66fb1fb5fadd9a87b555cca4e26f09 100644 (file)
@@ -293,14 +293,22 @@ namespace SUNDIALS
 
 
 
+  template <typename VectorType>
+  KINSOL<VectorType>::KINSOL(const AdditionalData &data)
+    : KINSOL(data, MPI_COMM_SELF)
+  {}
+
+
+
   template <typename VectorType>
   KINSOL<VectorType>::KINSOL(const AdditionalData &data,
                              const MPI_Comm &      mpi_comm)
     : data(data)
-    , mpi_communicator(is_serial_vector<VectorType>::value ?
-                         MPI_COMM_SELF :
-                         Utilities::MPI::duplicate_communicator(mpi_comm))
+    , mpi_communicator(mpi_comm)
     , kinsol_mem(nullptr)
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    , kinsol_ctx(nullptr)
+#  endif
     , solution(nullptr)
     , u_scale(nullptr)
     , f_scale(nullptr)
@@ -308,7 +316,15 @@ namespace SUNDIALS
     set_functions_to_trigger_an_assert();
 
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-    const int status = SUNContext_Create(&mpi_communicator, &kinsol_ctx);
+    // SUNDIALS will always duplicate communicators if we provide them. This
+    // can cause problems if SUNDIALS is configured with MPI and we pass along
+    // MPI_COMM_SELF in a serial application as MPI won't be
+    // initialized. Hence, work around that by just not providing a
+    // communicator in that case.
+    const int status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
+                                                            &mpi_communicator,
+                        &kinsol_ctx);
     (void)status;
     AssertKINSOL(status);
 #  endif
@@ -319,22 +335,12 @@ namespace SUNDIALS
   template <typename VectorType>
   KINSOL<VectorType>::~KINSOL()
   {
-    int status = 0;
-    (void)status;
-
-    if (kinsol_mem)
-      {
-        KINFree(&kinsol_mem);
-      }
+    KINFree(&kinsol_mem);
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-    status = SUNContext_Free(&kinsol_ctx);
+    const int status = SUNContext_Free(&kinsol_ctx);
+    (void)status;
     AssertKINSOL(status);
 #  endif
-
-#  ifdef DEAL_II_WITH_MPI
-    if (is_serial_vector<VectorType>::value == false)
-      Utilities::MPI::free_communicator(mpi_communicator);
-#  endif
   }
 
 
@@ -361,20 +367,20 @@ namespace SUNDIALS
     int status = 0;
     (void)status;
 
-    if (kinsol_mem)
-      {
-        KINFree(&kinsol_mem);
-
+    KINFree(&kinsol_mem);
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-        status = SUNContext_Free(&kinsol_ctx);
-        AssertKINSOL(status);
+    status = SUNContext_Free(&kinsol_ctx);
+    AssertKINSOL(status);
 #  endif
-      }
 
 #  if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
     kinsol_mem = KINCreate();
 #  else
-    status = SUNContext_Create(&mpi_communicator, &kinsol_ctx);
+    // Same comment applies as in class constructor:
+    status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
+                                                            &mpi_communicator,
+                        &kinsol_ctx);
     AssertKINSOL(status);
 
     kinsol_mem = KINCreate(kinsol_ctx);

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.