const MPI_Comm mpi_comm) :
data(data),
arkode_mem(nullptr),
- communicator(Utilities::MPI::duplicate_communicator(mpi_comm))
+ communicator(is_serial_vector<VectorType>::value ?
+ MPI_COMM_SELF :
+ Utilities::MPI::duplicate_communicator(mpi_comm))
{
set_functions_to_trigger_an_assert();
}
if (arkode_mem)
ARKodeFree(&arkode_mem);
#ifdef DEAL_II_WITH_MPI
- MPI_Comm_free(&communicator);
+ if (is_serial_vector<VectorType>::value == false)
+ {
+ const int ierr = MPI_Comm_free(&communicator);
+ AssertThrowMPI(ierr);
+ }
#endif
}
const MPI_Comm mpi_comm) :
data(data),
ida_mem(nullptr),
- communicator(Utilities::MPI::duplicate_communicator(mpi_comm))
+ communicator(is_serial_vector<VectorType>::value ?
+ MPI_COMM_SELF :
+ Utilities::MPI::duplicate_communicator(mpi_comm))
{
set_functions_to_trigger_an_assert();
}
if (ida_mem)
IDAFree(&ida_mem);
#ifdef DEAL_II_WITH_MPI
- MPI_Comm_free(&communicator);
+ if (is_serial_vector<VectorType>::value == false)
+ {
+ const int ierr = MPI_Comm_free(&communicator);
+ AssertThrowMPI(ierr);
+ }
#endif
}
KINSOL<VectorType>::KINSOL(const AdditionalData &data, const MPI_Comm mpi_comm) :
data(data),
kinsol_mem(nullptr),
- communicator(Utilities::MPI::duplicate_communicator(mpi_comm))
+ communicator(is_serial_vector<VectorType>::value ?
+ MPI_COMM_SELF :
+ Utilities::MPI::duplicate_communicator(mpi_comm))
{
set_functions_to_trigger_an_assert();
}
if (kinsol_mem)
KINFree(&kinsol_mem);
#ifdef DEAL_II_WITH_MPI
- MPI_Comm_free(&communicator);
+ if (is_serial_vector<VectorType>::value == false)
+ {
+ const int ierr = MPI_Comm_free(&communicator);
+ AssertThrowMPI(ierr);
+ }
#endif
}
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2017 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.
+//
+//-----------------------------------------------------------
+
+// Make sure that we can initialize the SUNDIALS solvers
+// for serial vector types without initializing MPI
+
+#include "../tests.h"
+#include <deal.II/lac/vector.h>
+#include <deal.II/sundials/arkode.h>
+#include <deal.II/sundials/ida.h>
+#include <deal.II/sundials/kinsol.h>
+
+
+int main()
+{
+ initlog();
+
+ SUNDIALS::IDA<Vector<double>> ida_solver;
+ (void) ida_solver;
+ deallog << "IDA OK" << std::endl;
+
+ SUNDIALS::KINSOL<Vector<double>> kinsol_solver;
+ (void) kinsol_solver;
+ deallog << "KINSOL OK" << std::endl;
+
+ SUNDIALS::ARKode<Vector<double>> arkode_solver;
+ (void) arkode_solver;
+ deallog << "ARKODE OK" << std::endl;
+
+ return 0;
+}
+