* too many details. All the algorithms described below are implement in
* classes and functions in namespace parallel::distributed.
*
+ * One important aspect in parallel computations using MPI is that write
+ * access to matrix and vector elements requires a call to compress() after
+ * the operation is finished and before the object is used (for example read
+ * from). Also see @ref GlossCompress.
*
* <h4>Other resources</h4>
*
+
//-------------------------------------------------------------------------
// $Id$
// Version: $Name$
* ways:
*
* - You tell the object that you want to compress what operation is
- * intended. The TrilinosWrappers::VectorBase::compress() can take
- * such an additional argument.
+ * intended. This can be done using the VectorOperation argument in
+ * the various compress() functions.
* - You do a fake addition or set operation on the object in question. For
* example, you can add a zero to an element of the matrix or vector,
* which has no effect other than telling the object that the next
<ol>
+<li>
+Changed: unify the concept of compress() for all linear algebra
+objects. Introduce type VectorOperation to decide between
+add and insert. Implement also for serial vectors. Note:
+this breaks distributed::vector::compress(bool).
+<br>
+(Timo Heister, 2012/08/13)
+
<li>
Changed: Support for the METIS 4.x has been replaced with support for
METIS 5.x. Use <code>--with-metis=path/to/metis</code> to configure
#include <deal.II/lac/exceptions.h>
#include <deal.II/lac/full_matrix.h>
#include <deal.II/lac/matrix_iterator.h>
+#include <deal.II/lac/vector.h>
#include <cmath>
/**
* Call the compress() function on all
* the subblocks of the matrix.
+ *
+ *
+ * See @ref GlossCompress "Compressing
+ * distributed objects" for more
+ * information.
*/
- void compress ();
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
/**
* Multiply the entire matrix by a
template <class MatrixType>
inline
void
-BlockMatrixBase<MatrixType>::compress ()
+BlockMatrixBase<MatrixType>::compress (::dealii::VectorOperation::values operation)
{
for (unsigned int r=0; r<n_block_rows(); ++r)
for (unsigned int c=0; c<n_block_cols(); ++c)
- block(r,c).compress ();
+ block(r,c).compress (operation);
}
*/
~BlockVector ();
+ /**
+ * Call the compress() function on all
+ * the subblocks.
+ *
+ * This functionality only needs to be
+ * called if using MPI based vectors and
+ * exists in other objects for
+ * compatibility.
+ *
+ * See @ref GlossCompress "Compressing
+ * distributed objects" for more
+ * information.
+ */
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
+
/**
* Copy operator: fill all components of
* the vector with the given scalar
return *this;
}
+template <typename Number>
+inline
+void BlockVector<Number>::compress (::dealii::VectorOperation::values operation)
+{
+ for (unsigned int i=0; i<this->n_blocks();++i)
+ this->components[i].compress(operation);
+}
+
template <typename Number>
/**
* Call the compress() function on all
* the subblocks of the matrix.
+ *
+ * This functionality only needs to be
+ * called if using MPI based vectors and
+ * exists in other objects for
+ * compatibility.
+ *
+ * See @ref GlossCompress "Compressing
+ * distributed objects" for more
+ * information.
*/
- void compress ();
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
+
/**
* Access to a single block.
template <class VectorType>
inline
void
-BlockVectorBase<VectorType>::compress ()
+BlockVectorBase<VectorType>::compress (::dealii::VectorOperation::values operation)
{
for (unsigned int i=0; i<n_blocks(); ++i)
- block(i).compress ();
+ block(i).compress (operation);
}
/**
* This function copies the data that has
* accumulated in the data buffer for ghost
- * indices to the owning processor. If the
- * optional argument @p add_ghost_data is set
- * to true, the data is added into the
- * respective positions of the owning
- * processor, otherwise the new data
- * overwrites the old content in the host
- * vector. In that case, data coming from
- * different processor to the same target
- * entry should be identical. However, no
- * checking is performed, so it is the user's
- * responsibility to ensure consistency of
- * data.
+ * indices to the owning processor.
*
- * For the meaning of this argument, see the
- * entry on
- * @ref GlossCompress "Compressing distributed vectors and matrices"
+ * For the meaning of this argument,
+ * see the entry on @ref
+ * GlossCompress "Compressing
+ * distributed vectors and matrices"
* in the glossary.
*/
- void compress (const bool add_ghost_data = true);
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
+
/**
* Fills the data field for ghost indices with
template <typename Number>
inline
void
- Vector<Number>::compress (const bool add_ghost_data)
+ Vector<Number>::compress (::dealii::VectorOperation::values operation)
{
compress_start ();
- compress_finish (add_ghost_data);
+ if (operation == ::dealii::VectorOperation::insert)
+ compress_finish (false);
+ else
+ compress_finish (true);
}
# include <deal.II/base/subscriptor.h>
# include <deal.II/lac/full_matrix.h>
# include <deal.II/lac/exceptions.h>
+# include <deal.II/lac/vector.h>
# include <petscmat.h>
# include <deal.II/base/std_cxx1x/shared_ptr.h>
* for more information.
* more information.
*/
- void compress ();
-
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
/**
* Return the value of the entry
* (<i>i,j</i>). This may be an
for (unsigned int i=0; i<v.size(); ++i)
(*this)(i) = v(i);
- compress ();
+ compress (::dealii::VectorOperation::insert);
return *this;
}
for (unsigned int i=0; i<v.size(); ++i)
(*this)(i) = v(i);
- compress ();
+ compress (::dealii::VectorOperation::insert);
return *this;
}
# include <deal.II/base/subscriptor.h>
# include <deal.II/lac/exceptions.h>
+# include <deal.II/lac/vector.h>
# include <vector>
# include <utility>
*
* See @ref GlossCompress "Compressing distributed objects"
* for more information.
- * more information.
*/
- void compress ();
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
/**
* Set all components of the vector to
*/
IndexSet ghost_indices;
- /**
- * PETSc doesn't allow to mix additions
- * to matrix entries and overwriting
- * them (to make synchronisation of
- * parallel computations
- * simpler). Since the interface of the
- * existing classes don't support the
- * notion of not interleaving things,
- * we have to emulate this
- * ourselves. The way we do it is to,
- * for each access operation, store
- * whether it is an insertion or an
- * addition. If the previous one was of
- * different type, then we first have
- * to flush the PETSc buffers;
- * otherwise, we can simply go on.
- *
- * The following structure and variable
- * declare and store the previous
- * state.
- */
- struct LastAction
- {
- enum Values { none, insert, add };
- };
-
-
/**
* Store whether the last action was a
* write or add operation. This
* even though the vector object they
* refer to is constant.
*/
- mutable LastAction::Values last_action;
+ mutable ::dealii::VectorOperation::values last_action;
/**
* Make the reference class a friend.
const VectorReference &
VectorReference::operator = (const PetscScalar &value) const
{
- Assert ((vector.last_action == VectorBase::LastAction::insert)
+ Assert ((vector.last_action == VectorOperation::insert)
||
- (vector.last_action == VectorBase::LastAction::none),
- ExcWrongMode (VectorBase::LastAction::insert,
+ (vector.last_action == VectorOperation::unknown),
+ ExcWrongMode (VectorOperation::insert,
vector.last_action));
#ifdef PETSC_USE_64BIT_INDICES
= VecSetValues (vector, 1, &petsc_i, &value, INSERT_VALUES);
AssertThrow (ierr == 0, ExcPETScError(ierr));
- vector.last_action = VectorBase::LastAction::insert;
+ vector.last_action = VectorOperation::insert;
return *this;
}
const VectorReference &
VectorReference::operator += (const PetscScalar &value) const
{
- Assert ((vector.last_action == VectorBase::LastAction::add)
+ Assert ((vector.last_action == VectorOperation::add)
||
- (vector.last_action == VectorBase::LastAction::none),
- ExcWrongMode (VectorBase::LastAction::add,
+ (vector.last_action == VectorOperation::unknown),
+ ExcWrongMode (VectorOperation::add,
vector.last_action));
- vector.last_action = VectorBase::LastAction::add;
+ vector.last_action = VectorOperation::add;
// we have to do above actions in any
// case to be consistent with the MPI
const VectorReference &
VectorReference::operator -= (const PetscScalar &value) const
{
- Assert ((vector.last_action == VectorBase::LastAction::add)
+ Assert ((vector.last_action == VectorOperation::add)
||
- (vector.last_action == VectorBase::LastAction::none),
- ExcWrongMode (VectorBase::LastAction::add,
+ (vector.last_action == VectorOperation::unknown),
+ ExcWrongMode (VectorOperation::add,
vector.last_action));
- vector.last_action = VectorBase::LastAction::add;
+ vector.last_action = VectorOperation::add;
// we have to do above actions in any
// case to be consistent with the MPI
const VectorReference &
VectorReference::operator *= (const PetscScalar &value) const
{
- Assert ((vector.last_action == VectorBase::LastAction::insert)
+ Assert ((vector.last_action == VectorOperation::insert)
||
- (vector.last_action == VectorBase::LastAction::none),
- ExcWrongMode (VectorBase::LastAction::insert,
+ (vector.last_action == VectorOperation::unknown),
+ ExcWrongMode (VectorOperation::insert,
vector.last_action));
- vector.last_action = VectorBase::LastAction::insert;
+ vector.last_action = VectorOperation::insert;
// we have to do above actions in any
// case to be consistent with the MPI
const VectorReference &
VectorReference::operator /= (const PetscScalar &value) const
{
- Assert ((vector.last_action == VectorBase::LastAction::insert)
+ Assert ((vector.last_action == VectorOperation::insert)
||
- (vector.last_action == VectorBase::LastAction::none),
- ExcWrongMode (VectorBase::LastAction::insert,
+ (vector.last_action == VectorOperation::unknown),
+ ExcWrongMode (VectorOperation::insert,
vector.last_action));
- vector.last_action = VectorBase::LastAction::insert;
+ vector.last_action = VectorOperation::insert;
// we have to do above actions in any
// case to be consistent with the MPI
void reinit (const ::dealii::BlockSparseMatrix<double> &deal_ii_sparse_matrix,
const double drop_tolerance=1e-13);
- /**
- * This function calls the compress()
- * command of all matrices after
- * the assembly is
- * completed. Note that all MPI
- * processes need to call this
- * command (whereas the individual
- * assembly routines will most probably
- * only be called on each processor
- * individually) before any
- * can complete it.
- */
- void compress ();
-
/**
* Returns the state of the
* matrix, i.e., whether
*/
~BlockVector ();
+ /**
+ * use compress(VectorOperation) instead
+ *
+ * @deprecated
+ *
+ * See @ref GlossCompress "Compressing
+ * distributed objects" for more
+ * information.
+ */
+ void compress (const Epetra_CombineMode last_action);
+
+ /**
+ * so it is not hidden
+ */
+ using BlockVectorBase<Vector>::compress;
+
/**
* Copy operator: fill all
* components of the vector that
}
+ inline
+ void
+ BlockVector::compress (const Epetra_CombineMode last_action)
+ {
+ ::dealii::VectorOperation::values last_action_;
+ if (last_action == Add)
+ last_action_ = ::dealii::VectorOperation::add;
+ else if (last_action == Insert)
+ last_action_ = ::dealii::VectorOperation::insert;
+ else
+ AssertThrow(false, ExcNotImplemented());
+
+ this->compress(last_action_);
+ }
+
inline
void
BlockVector::swap (BlockVector &v)
void import_nonlocal_data_for_fe (const TrilinosWrappers::BlockSparseMatrix &m,
const BlockVector &v);
- /**
- * Compress the underlying
- * representation of the Trilinos
- * object, i.e. flush the buffers
- * of the vector object if it has
- * any. This function is
- * necessary after writing into a
- * vector element-by-element and
- * before anything else can be
- * done on it.
- *
- * The (defaulted) argument can
- * be used to specify the
- * compress mode
- * (<code>Add</code> or
- * <code>Insert</code>) in case
- * the vector has not been
- * written to since the last
- * time this function was
- * called. The argument is
- * ignored if the vector has
- * been added or written to
- * since the last time
- * compress() was called.
- *
- * See @ref GlossCompress "Compressing distributed objects"
- * for more information.
- * more information.
- */
- void compress (const Epetra_CombineMode last_action = Zero);
+
+ /**
+ * use compress(VectorOperation) instead
+ *
+ * @deprecated
+ *
+ * See @ref GlossCompress "Compressing
+ * distributed objects" for more
+ * information.
+ */
+ void compress (const Epetra_CombineMode last_action);
+
+ /**
+ * so it is not hidden
+ */
+ using BlockVectorBase<Vector>::compress;
+
/**
* Returns the state of the
+ inline
+ void
+ BlockVector::compress (const Epetra_CombineMode last_action)
+ {
+ ::dealii::VectorOperation::values last_action_;
+ if (last_action == Add)
+ last_action_ = ::dealii::VectorOperation::add;
+ else if (last_action == Insert)
+ last_action_ = ::dealii::VectorOperation::insert;
+ else
+ AssertThrow(false, ExcNotImplemented());
+
+ this->compress(last_action_);
+ }
+
+
+
template <typename Number>
BlockVector &
BlockVector::operator = (const ::dealii::BlockVector<Number> &v)
* See @ref GlossCompress "Compressing distributed objects"
* for more information.
*/
- void compress ();
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
/**
* Set the element (<i>i,j</i>)
inline
void
- SparseMatrix::compress ()
+ SparseMatrix::compress (::dealii::VectorOperation::values operation)
{
// flush buffers
int ierr;
*
* See @ref GlossCompress "Compressing distributed objects"
* for more information.
- * more information.
*/
- void compress (const Epetra_CombineMode last_action = Zero);
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown);
+
+ /**
+ * @deprecated
+ */
+ void compress (const Epetra_CombineMode last_action);
/**
* Returns the state of the
inline
void
- VectorBase::compress (const Epetra_CombineMode given_last_action)
+ VectorBase::compress (const Epetra_CombineMode last_action)
{
- Epetra_CombineMode mode = (last_action != Zero) ?
- last_action : given_last_action;
+ ::dealii::VectorOperation::values last_action_;
+ if (last_action == Add)
+ last_action_ = ::dealii::VectorOperation::add;
+ else if (last_action == Insert)
+ last_action_ = ::dealii::VectorOperation::insert;
+ else
+ AssertThrow(false, ExcNotImplemented());
+
+ compress(last_action_);
+ }
+
+
+
+ inline
+ void
+ VectorBase::compress (::dealii::VectorOperation::values given_last_action)
+ {
+ //Select which mode to send to
+ //Trilinos. Note that we use last_action
+ //if available and ignore what the user
+ //tells us to detect wrongly mixed
+ //operations. Typically given_last_action
+ //is only used on machines that do not
+ //execute an operation (because they have
+ //no own cells for example).
+ Epetra_CombineMode mode = last_action;
+ if (last_action == Zero)
+ {
+ if (given_last_action==::dealii::VectorOperation::add)
+ mode = Add;
+ else if (given_last_action==::dealii::VectorOperation::insert)
+ mode = Insert;
+ }
+
#ifdef DEBUG
# ifdef DEAL_II_COMPILER_SUPPORTS_MPI
// check that every process has decided
*@{
*/
+/**
+ * This enum keeps track of the current operation in parallel linear algebra
+ * objects like Vectors and Matrices.
+ *
+ * It is used in the various compress() functions. They also exist in serial
+ * codes for compatibility and are empty there.
+ *
+ * See @ref GlossCompress "Compressing distributed objects" for more
+ * information.
+ */
+struct VectorOperation
+{
+ enum values { unknown, insert, add };
+};
+
/**
* Numerical vector of data. For this class there are different types
* this class, it is immaterial and thus
* an empty function.
*/
- void compress () const;
+ void compress (::dealii::VectorOperation::values operation
+ =::dealii::VectorOperation::unknown) const;
/**
* Change the dimension of the vector to
template <typename Number>
inline
void
-Vector<Number>::compress () const
+Vector<Number>::compress (::dealii::VectorOperation::values) const
{}
= function_values_scalar[fe_index][dof_to_rep_index_table[fe_index][i]];
}
}
+ vec.compress(::dealii::VectorOperation::insert);
}
namespace distributed
{
- namespace
- {
- template <class VECTOR>
- void compress_vector_insert(VECTOR & vec)
- {
- vec.compress();
- }
-
-#ifdef DEAL_II_USE_TRILINOS
- void compress_vector_insert(TrilinosWrappers::Vector & vec)
- {
- vec.compress(Insert);
- }
-
- void compress_vector_insert(TrilinosWrappers::BlockVector & vec)
- {
- for (unsigned int i=0;i<vec.n_blocks();++i)
- vec.block(i).compress(Insert);
- }
-
- void compress_vector_insert(TrilinosWrappers::MPI::Vector & vec)
- {
- vec.compress(Insert);
- }
-
- void compress_vector_insert(TrilinosWrappers::MPI::BlockVector & vec)
- {
- for (unsigned int i=0;i<vec.n_blocks();++i)
- vec.block(i).compress(Insert);
- }
-#endif
- }
-
-
-
template<int dim, typename VECTOR, class DH>
SolutionTransfer<dim, VECTOR, DH>::SolutionTransfer(const DH &dof)
:
std_cxx1x::_3,
std_cxx1x::ref(all_out)));
+
for (typename std::vector<VECTOR*>::iterator it=all_out.begin();
it !=all_out.end();
++it)
- compress_vector_insert(*(*it));
+ (*it)->compress(::dealii::VectorOperation::insert);
input_vectors.clear();
}
// constraints, so we need to explicitly
// state, that the others are doing an
// insert here:
- vec.compress (Insert);
+ vec.compress (::dealii::VectorOperation::insert);
}
it->entries[i].second);
vec(it->line) = new_value;
}
- vec.block(block).compress(Insert);
+ vec.block(block).compress(::dealii::VectorOperation::insert);
}
}
void
- MatrixBase::compress ()
+ MatrixBase::compress (::dealii::VectorOperation::values operation)
{
// flush buffers
int ierr;
VectorBase::VectorBase ()
:
ghosted(false),
- last_action (LastAction::none),
+ last_action (::dealii::VectorOperation::unknown),
attained_ownership(true)
{}
Subscriptor (),
ghosted(v.ghosted),
ghost_indices(v.ghost_indices),
- last_action (LastAction::none),
+ last_action (::dealii::VectorOperation::unknown),
attained_ownership(true)
{
int ierr = VecDuplicate (v.vector, &vector);
Subscriptor (),
vector(v),
ghosted(false),
- last_action (LastAction::none),
+ last_action (::dealii::VectorOperation::unknown),
attained_ownership(false)
{}
void
- VectorBase::compress ()
+ VectorBase::compress (::dealii::VectorOperation::values operation)
{
// note that one may think that we only need to do something
- // if in fact the state is anything but LastAction::none. but
+ // if in fact the state is anything but last_action::unknown. but
// that's not true: one frequently gets into situations where
// only one processor (or a subset of processors) actually writes
// something into a vector, but we still need to call
- // VecAssemblyBegin/End on all processors, even those that are
- // still in LastAction::none state
+ // VecAssemblyBegin/End on all processors.
int ierr;
ierr = VecAssemblyBegin(vector);
AssertThrow (ierr == 0, ExcPETScError(ierr));
// reset the last action field to indicate that we're back to
// a pristine state
- last_action = VectorBase::LastAction::none;
+ last_action = ::dealii::VectorOperation::unknown;
}
std::size_t
VectorBase::memory_consumption () const
{
- std::size_t mem = sizeof(Vec)+sizeof(LastAction::Values)
+ std::size_t mem = sizeof(Vec)+sizeof(last_action)
+MemoryConsumption::memory_consumption(ghosted)
+MemoryConsumption::memory_consumption(ghost_indices);
const PetscScalar *values,
const bool add_values)
{
- Assert ((last_action == (add_values ?
- LastAction::add :
- LastAction::insert))
+ ::dealii::VectorOperation::values action = (add_values ?
+ ::dealii::VectorOperation::add :
+ ::dealii::VectorOperation::insert);
+ Assert ((last_action == action)
||
- (last_action == LastAction::none),
- internal::VectorReference::ExcWrongMode ((add_values ?
- LastAction::add :
- LastAction::insert),
+ (last_action == ::dealii::VectorOperation::unknown),
+ internal::VectorReference::ExcWrongMode (action,
last_action));
// VecSetValues complains if we
// set the mode here, independent of whether we have actually
// written elements or whether the list was empty
- last_action = (add_values ?
- VectorBase::LastAction::add :
- VectorBase::LastAction::insert);
+ last_action = action;
}
- void
- BlockSparseMatrix::compress()
- {
- for (unsigned int r=0; r<this->n_block_rows(); ++r)
- for (unsigned int c=0; c<this->n_block_cols(); ++c)
- this->block(r,c).compress();
- }
-
void
- void
- BlockVector::compress (const Epetra_CombineMode last_action)
- {
- for (unsigned int i=0; i<n_blocks(); ++i)
- components[i].compress(last_action);
- }
-
-
-
void BlockVector::print (std::ostream &out,
const unsigned int precision,
const bool scientific,