#include <vector>
-//TODO: How does a matrix or vector get its communicator if zero-initialized and later reinit'ed?
namespace PETScWrappers
{
* this flag has been set. The default
* value of this flag is @p{false}.
*/
- SparseMatrix (const MPI_Comm &communicator,
- const unsigned int m,
- const unsigned int n,
- const unsigned int local_rows,
- const unsigned int n_nonzero_per_row,
- const bool is_symmetric = false);
+ SparseMatrix (const MPI_Comm &communicator,
+ const unsigned int m,
+ const unsigned int n,
+ const unsigned int local_rows,
+ const unsigned int n_nonzero_per_row,
+ const bool is_symmetric = false);
/**
* Initialize a rectangular matrix with
SparseMatrix (const MPI_Comm &communicator,
const unsigned int m,
const unsigned int n,
- const unsigned int local_rows,
+ const unsigned int local_rows,
const std::vector<unsigned int> &row_lengths,
const bool is_symmetric = false);
* the same argument list as the
* present function.
*/
- void reinit (const MPI_Comm &communicator,
- const unsigned int m,
- const unsigned int n,
- const unsigned int local_rows,
- const unsigned int n_nonzero_per_row,
- const bool is_symmetric = false);
+ void reinit (const MPI_Comm &communicator,
+ const unsigned int m,
+ const unsigned int n,
+ const unsigned int local_rows,
+ const unsigned int n_nonzero_per_row,
+ const bool is_symmetric = false);
/**
* Throw away the present matrix and
*/
void do_reinit (const unsigned int m,
const unsigned int n,
- const unsigned int local_rows,
+ const unsigned int local_rows,
const std::vector<unsigned int> &row_lengths,
const bool is_symmetric = false);
/**
* Constructor. Set dimension to
- * @p{n} and initialize all
+ * @p n and initialize all
* elements with zero.
*
+ * @arg local_size denotes the size
+ * of the chunk that shall be stored
+ * on the present process.
+ *
+ * @arg communicator denotes the MPI
+ * communicator over which the
+ * different parts of the vector
+ * shall communicate
+ *
* The constructor is made explicit
* to avoid accidents like this:
* @p{v=0;}. Presumably, the user
* vector is replaced by one of
* length zero.
*/
- explicit Vector (const unsigned int n,
- const unsigned int local_size,
- const MPI_Comm &communicator);
+ explicit Vector (const MPI_Comm &communicator,
+ const unsigned int n,
+ const unsigned int local_size);
/**
* Copy-constructor from deal.II
* vectors. Sets the dimension to that
* of the given vector, and copies all
* elements.
+ *
+ * @arg local_size denotes the size
+ * of the chunk that shall be stored
+ * on the present process.
+ *
+ * @arg communicator denotes the MPI
+ * communicator over which the
+ * different parts of the vector
+ * shall communicate
*/
template <typename Number>
- explicit Vector (const ::Vector<Number> &v,
- const unsigned int local_size,
- const MPI_Comm &communicator);
+ explicit Vector (const MPI_Comm &communicator,
+ const ::Vector<Number> &v,
+ const unsigned int local_size);
/**
* Copy-constructor the
* values from a PETSc wrapper vector
* class.
+ *
+ * @arg local_size denotes the size
+ * of the chunk that shall be stored
+ * on the present process.
+ *
+ * @arg communicator denotes the MPI
+ * communicator over which the
+ * different parts of the vector
+ * shall communicate
*/
- explicit Vector (const VectorBase &v,
- const unsigned int local_size,
- const MPI_Comm &communicator);
+ explicit Vector (const MPI_Comm &communicator,
+ const VectorBase &v,
+ const unsigned int local_size);
/**
* Copy the given vector. Resize the
* elements are left an unspecified
* state.
*/
- void reinit (const unsigned int N,
- const unsigned int local_size,
- const MPI_Comm &communicator,
- const bool fast = false);
+ void reinit (const MPI_Comm &communicator,
+ const unsigned int N,
+ const unsigned int local_size,
+ const bool fast = false);
/**
* Change the dimension to that of
template <typename number>
- Vector::Vector (const ::Vector<number> &v,
- const unsigned int local_size,
- const MPI_Comm &communicator)
+ Vector::Vector (const MPI_Comm &communicator,
+ const ::Vector<number> &v,
+ const unsigned int local_size)
:
communicator (communicator)
{
- Vector::Vector (const unsigned int n,
- const unsigned int local_size,
- const MPI_Comm &communicator)
+ Vector::Vector (const MPI_Comm &communicator,
+ const unsigned int n,
+ const unsigned int local_size)
:
communicator (communicator)
{
- Vector::Vector (const VectorBase &v,
- const unsigned int local_size,
- const MPI_Comm &communicator)
+ Vector::Vector (const MPI_Comm &communicator,
+ const VectorBase &v,
+ const unsigned int local_size)
:
communicator (communicator)
{
void
- Vector::reinit (const unsigned int n,
+ Vector::reinit (const MPI_Comm &comm,
+ const unsigned int n,
const unsigned int local_sz,
- const MPI_Comm &comm,
const bool fast)
{
communicator = comm;
Assert (local_size <= n, ExcIndexRange (local_size, 0, n));
const int ierr
- = VecCreateMPI (PETSC_COMM_SELF, local_size, n, &vector);
+ = VecCreateMPI (communicator, local_size, PETSC_DETERMINE,
+ &vector);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ Assert (size() == n, ExcInternalError());
}
}