]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add MPI communicator to SparseDirectMUMPS constructor. Use integers from MUMPS 18475/head
authorMarco Feder <marco.feder@sissa.it>
Sat, 24 May 2025 16:07:41 +0000 (18:07 +0200)
committerMarco Feder <marco.feder@sissa.it>
Sat, 24 May 2025 17:29:55 +0000 (19:29 +0200)
include/deal.II/lac/sparse_direct.h
source/lac/sparse_direct.cc

index 236230f1ae19428485eb4ba423bf2f051843f86c..c1ee01d1c738b2e7532f98723e27cefc7a2f14ce 100644 (file)
@@ -48,6 +48,16 @@ namespace types
 #else
   using suitesparse_index = long int;
 #endif
+
+#ifdef DEAL_II_WITH_MUMPS
+  using mumps_index = MUMPS_INT;
+  using mumps_nnz   = MUMPS_INT8;
+#else
+  using mumps_index       = int;
+  using mumps_nnz         = std::size_t;
+#endif
+
+
 } // namespace types
 
 /**
@@ -475,12 +485,12 @@ public:
         : blr_ucfs(blr_ucfs)
         , lowrank_threshold(lowrank_threshold)
       {}
-
       /**
        * If true, the Block Low-Rank approximation is used with the UCFS
        * algorithm, an algorithm with higher compression than the standard one.
        */
       bool blr_ucfs;
+
       /**
        * Threshold for the low-rank truncation of the blocks.
        */
@@ -510,11 +520,13 @@ public:
      * If true, the MUMPS solver will print out error statistics.
      */
     bool error_statistics;
+
     /*
      * If true, the MUMPS solver will use the symmetric factorization. This is
      * only possible if the matrix is symmetric.
      */
     bool symmetric;
+
     /*
      * If true, the MUMPS solver will use the positive definite factorization.
      * This is only possible if the matrix is symmetric and positive definite.
@@ -536,7 +548,8 @@ public:
   /**
    * Constructor, takes <tt>AdditionalData</tt> to control MUMPS execution.
    */
-  SparseDirectMUMPS(const AdditionalData &additional_data = AdditionalData());
+  SparseDirectMUMPS(const AdditionalData &additional_data = AdditionalData(),
+                    const MPI_Comm       &communicator    = MPI_COMM_WORLD);
 
   /**
    * Destructor.
@@ -590,6 +603,7 @@ public:
 private:
 #ifdef DEAL_II_WITH_MUMPS
   mutable DMUMPS_STRUC_C id;
+
 #endif // DEAL_II_WITH_MUMPS
 
   /**
@@ -606,12 +620,12 @@ private:
   /**
    * irn contains the row indices of the non-zero entries of the matrix.
    */
-  std::unique_ptr<int[]> irn;
+  std::unique_ptr<types::mumps_index[]> irn;
 
   /**
    * jcn contains the column indices of the non-zero entries of the matrix.
    */
-  std::unique_ptr<int[]> jcn;
+  std::unique_ptr<types::mumps_index[]> jcn;
 
   /**
    * The number of rows of the matrix. The matrix is square.
@@ -621,7 +635,7 @@ private:
   /**
    * The number of non-zero entries in the matrix.
    */
-  std::size_t nnz;
+  types::mumps_nnz nnz;
 
   /**
    * This function hands over to MUMPS the system's <tt>matrix</tt>.
@@ -646,6 +660,11 @@ private:
    * Struct that holds the additional data for the MUMPS solver.
    */
   AdditionalData additional_data;
+
+  /**
+   * MPI_Comm object for the MUMPS solver.
+   */
+  const MPI_Comm mpi_communicator;
 };
 
 DEAL_II_NAMESPACE_CLOSE
index f4d6e6377f35c6b2f6edba82d0a4fba19fa160b2..b4b695b75cfe3836984dbb77c46d19590cee787d 100644 (file)
@@ -861,9 +861,11 @@ SparseDirectUMFPACK::n() const
 
 #ifdef DEAL_II_WITH_MUMPS
 
-SparseDirectMUMPS::SparseDirectMUMPS(const AdditionalData &data)
+SparseDirectMUMPS::SparseDirectMUMPS(const AdditionalData &data,
+                                     const MPI_Comm       &communicator)
+  : additional_data(data)
+  , mpi_communicator(communicator)
 {
-  this->additional_data = data;
   // Initialize MUMPS instance:
   id.job = -1;
   id.par = 1;
@@ -883,8 +885,7 @@ SparseDirectMUMPS::SparseDirectMUMPS(const AdditionalData &data)
   else
     id.sym = 0;
 
-  // Use MPI_COMM_WORLD as communicator
-  id.comm_fortran = (MUMPS_INT)MPI_Comm_c2f(MPI_COMM_WORLD);
+  id.comm_fortran = (MUMPS_INT)MPI_Comm_c2f(mpi_communicator);
   dmumps_c(&id);
 
   if (additional_data.output_details == false)
@@ -926,12 +927,9 @@ void
 SparseDirectMUMPS::initialize_matrix(const Matrix &matrix)
 {
   Assert(matrix.n() == matrix.m(), ExcMessage("Matrix needs to be square."));
-  // Here we should be checking if the matrix is respecting the symmetry given
-  //(I.E. sym = 0 for non-symmetric matrix, sym = 1 for posdef matrix, sym = 2
-  // for general symmetric).
 
   // Hand over matrix to MUMPS as centralized assembled matrix
-  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+  if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
     {
       // Set number of unknowns
       n = matrix.n();
@@ -945,8 +943,8 @@ SparseDirectMUMPS::initialize_matrix(const Matrix &matrix)
       // matrix indices pointing to the row and column dimensions
       // respectively of the matrix representation above (a): ie. a[k] is
       // the matrix element (irn[k], jcn[k])
-      irn = std::make_unique<int[]>(nnz);
-      jcn = std::make_unique<int[]>(nnz);
+      irn = std::make_unique<MUMPS_INT[]>(nnz);
+      jcn = std::make_unique<MUMPS_INT[]>(nnz);
 
       size_type n_non_zero_elements = 0;
 
@@ -1001,7 +999,7 @@ SparseDirectMUMPS::copy_rhs_to_mumps(const Vector<double> &new_rhs) const
   Assert(n == new_rhs.size(),
          ExcMessage("Matrix size and rhs length must be equal."));
 
-  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+  if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
     {
       rhs.resize(n);
       for (size_type i = 0; i < n; ++i)
@@ -1022,7 +1020,7 @@ SparseDirectMUMPS::copy_solution(Vector<double> &vector) const
          ExcMessage("Class not initialized with a rhs vector."));
 
   // Copy solution into the given vector
-  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+  if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
     {
       for (size_type i = 0; i < n; ++i)
         vector(i) = rhs[i];

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.