*/
using size_type = types::global_dof_index;
+
/**
* The AdditionalData struct contains data for controlling the MUMPS
* execution.
* 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.
* a contains the actual values of the matrix. a[k] is the value of the matrix
* entry (i,j) if irn[k] == i and jcn[k] == j.
*/
- double *a;
+ std::unique_ptr<double[]> a;
/**
* The right-hand side vector. This is the right hand side of the system.
/**
* irn contains the row indices of the non-zero entries of the matrix.
*/
- int *irn;
+ std::unique_ptr<int[]> irn;
/**
* jcn contains the column indices of the non-zero entries of the matrix.
*/
- int *jcn;
+ std::unique_ptr<int[]> jcn;
/**
* The number of rows of the matrix. The matrix is square.
id.sym = 0;
// Use MPI_COMM_WORLD as communicator
- id.comm_fortran = -987654;
+ id.comm_fortran = (MUMPS_INT)MPI_Comm_c2f(MPI_COMM_WORLD);
dmumps_c(&id);
if (additional_data.output_details == false)
SparseDirectMUMPS::~SparseDirectMUMPS()
{
+ // MUMPS destructor
id.job = -2;
dmumps_c(&id);
-
- // Do some cleaning
- if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
- {
- delete[] a;
- delete[] irn;
- delete[] jcn;
- }
}
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)
nnz = matrix.n_actually_nonzero_elements();
// representation of the matrix
- a = new double[nnz];
+ a = std::make_unique<double[]>(nnz);
// 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 = new int[nnz];
- jcn = new int[nnz];
+ irn = std::make_unique<int[]>(nnz);
+ jcn = std::make_unique<int[]>(nnz);
size_type n_non_zero_elements = 0;
}
id.n = n;
id.nnz = n_non_zero_elements;
- id.irn = irn;
- id.jcn = jcn;
- id.a = a;
+ id.irn = irn.get();
+ id.jcn = jcn.get();
+ id.a = a.get();
}
}
dmumps_c(&id);
}
-
-
void
SparseDirectMUMPS::vmult(Vector<double> &dst, const Vector<double> &src) const
{
// and that the matrix has at least one nonzero element:
Assert(nnz != 0, ExcNotInitialized());
-
Assert(n == dst.size(), ExcMessage("Destination vector has the wrong size."));
Assert(n == src.size(), ExcMessage("Source vector has the wrong size."));