#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
/**
: 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.
*/
* 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.
/**
* 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.
private:
#ifdef DEAL_II_WITH_MUMPS
mutable DMUMPS_STRUC_C id;
+
#endif // DEAL_II_WITH_MUMPS
/**
/**
* 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.
/**
* 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>.
* 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
#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;
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)
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();
// 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;
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)
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];