15#ifndef dealii_trilinos_vector_h
16#define dealii_trilinos_vector_h
21#ifdef DEAL_II_WITH_TRILINOS
33# include <Epetra_ConfigDefs.h>
34# include <Epetra_FEVector.h>
35# include <Epetra_LocalMap.h>
36# include <Epetra_Map.h>
37# include <Epetra_MpiComm.h>
50 template <
typename Number>
51 class ReadWriteVector;
102 VectorReference(MPI::Vector &vector,
const size_type index);
108 VectorReference(
const VectorReference &) =
default;
121 const VectorReference &
122 operator=(
const VectorReference &r)
const;
128 operator=(
const VectorReference &r);
133 const VectorReference &
139 const VectorReference &
145 const VectorReference &
151 const VectorReference &
157 const VectorReference &
171 <<
"An error with error number " << arg1
172 <<
" occurred while calling a Trilinos function");
187 friend class ::TrilinosWrappers::MPI::Vector;
194# ifndef DEAL_II_WITH_64BIT_INDICES
199 gid(
const Epetra_BlockMap &map,
int i)
208 gid(
const Epetra_BlockMap &map,
int i)
495 template <
typename Number>
497 const ::Vector<Number> &v,
629 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner,
695 template <
typename Number>
718 const ::TrilinosWrappers::SparseMatrix &matrix,
794 std::pair<size_type, size_type>
997 std::vector<TrilinosScalar> &values)
const;
1033 template <
typename ForwardIterator,
typename OutputIterator>
1086 set(
const std::vector<size_type> &indices,
1087 const std::vector<TrilinosScalar> &values);
1094 set(
const std::vector<size_type> &indices,
1095 const ::Vector<TrilinosScalar> &values);
1112 add(
const std::vector<size_type> &indices,
1113 const std::vector<TrilinosScalar> &values);
1120 add(
const std::vector<size_type> &indices,
1121 const ::Vector<TrilinosScalar> &values);
1231 const Epetra_MultiVector &
1256 print(std::ostream &out,
1257 const unsigned int precision = 3,
1258 const bool scientific =
true,
1259 const bool across =
true)
const;
1300 <<
"An error with error number " <<
arg1
1301 <<
" occurred while calling a Trilinos function");
1312 <<
"You are trying to access element " <<
arg1
1313 <<
" of a distributed vector, but this element is not stored "
1314 <<
"on the current processor. Note: There are " <<
arg2
1315 <<
" elements stored "
1316 <<
"on the current processor from within the range [" <<
arg3 <<
','
1317 <<
arg4 <<
"] but Trilinos vectors need not store contiguous "
1318 <<
"ranges on each processor, and not every element in "
1319 <<
"this range may in fact be stored locally."
1321 <<
"A common source for this kind of problem is that you "
1322 <<
"are passing a 'fully distributed' vector into a function "
1323 <<
"that needs read access to vector elements that correspond "
1324 <<
"to degrees of freedom on ghost cells (or at least to "
1325 <<
"'locally active' degrees of freedom that are not also "
1326 <<
"'locally owned'). You need to pass a vector that has these "
1327 <<
"elements as ghost entries.");
1401 inline VectorReference::VectorReference(
MPI::Vector &vector,
1402 const size_type index)
1408 inline const VectorReference &
1409 VectorReference::operator=(
const VectorReference &r)
const
1422 inline VectorReference &
1423 VectorReference::operator=(
const VectorReference &r)
1432 inline const VectorReference &
1435 vector.set(1, &index, &value);
1441 inline const VectorReference &
1444 vector.add(1, &index, &value);
1450 inline const VectorReference &
1454 vector.add(1, &index, &new_value);
1460 inline const VectorReference &
1464 vector.set(1, &index, &new_value);
1470 inline const VectorReference &
1474 vector.set(1, &index, &new_value);
1486 return ((index >=
range.first) && (index <
range.second));
1496 "The locally owned elements have not been properly initialized!"
1497 " This happens for example if this object has been initialized"
1498 " with exactly one overlapping IndexSet."));
1518 inline internal::VectorReference
1526 inline internal::VectorReference
1544 std::vector<TrilinosScalar> &values)
const
1546 for (
size_type i = 0; i < indices.size(); ++i)
1547 values[i] =
operator()(indices[i]);
1557 for (
unsigned int i = 0; i < indices.
size(); ++i)
1560 elements[i] = (*this)[indices[i]];
1566 template <
typename ForwardIterator,
typename OutputIterator>
1569 const ForwardIterator indices_end,
1570 OutputIterator values_begin)
const
1615 Vector::set(
const std::vector<size_type> &indices,
1616 const std::vector<TrilinosScalar> &values)
1624 set(indices.size(), indices.data(),
values.data());
1630 Vector::set(
const std::vector<size_type> &indices,
1631 const ::Vector<TrilinosScalar> &values)
1639 set(indices.size(), indices.data(),
values.begin());
1646 const size_type *indices,
1662 for (
size_type i = 0; i < n_elements; ++i)
1671 const int ierr =
vector->ReplaceGlobalValues(1, &row, &values[i]);
1686 Vector::add(
const std::vector<size_type> &indices,
1687 const std::vector<TrilinosScalar> &values)
1694 add(indices.size(), indices.data(),
values.data());
1700 Vector::add(
const std::vector<size_type> &indices,
1701 const ::Vector<TrilinosScalar> &values)
1708 add(indices.size(), indices.data(),
values.begin());
1715 const size_type *indices,
1732 for (
size_type i = 0; i < n_elements; ++i)
1741 const int ierr =
vector->SumIntoGlobalValues(
1758 "Attempted to write into off-processor vector entry "
1759 "that has not be specified as being writable upon "
1772# ifndef DEAL_II_WITH_64BIT_INDICES
1773 return vector->Map().MaxAllGID() + 1 -
vector->Map().MinAllGID();
1775 return vector->Map().MaxAllGID64() + 1 -
vector->Map().MinAllGID64();
1789 inline std::pair<Vector::size_type, Vector::size_type>
1792# ifndef DEAL_II_WITH_64BIT_INDICES
1795 vector->Map().MaxMyGID() + 1;
1798 vector->Map().MinMyGID64();
1800 vector->Map().MaxMyGID64() + 1;
1806 "This function only makes sense if the elements that this "
1807 "vector stores on the current processor form a contiguous range. "
1808 "This does not appear to be the case for the current vector."));
1998 const int ierr =
vector->Update(1.0, *(v.vector), 1.0);
2013 const int ierr =
vector->Update(-1.0, *(v.vector), 1.0);
2046 const int ierr =
vector->Update(a, *(v.vector), 1.);
2067 const int ierr =
vector->Update(a, *(v.vector), b, *(
w.vector), 1.);
2089 Assert(this->
vector->Map().SameAs(v.vector->Map()) ==
true,
2091 const int ierr =
vector->Update(1., *(v.vector), s);
2120 Assert(this->
vector->Map().SameAs(v.vector->Map()) ==
true,
2122 const int ierr =
vector->Update(a, *(v.vector), s);
2130 this->
add(tmp,
true);
2159 if (
vector->Map().SameAs(v.vector->Map()) ==
false)
2161 this->
sadd(0., a, v);
2166 int ierr =
vector->Update(a, *v.vector, 0.0);
2175 inline const Epetra_MultiVector &
2178 return static_cast<const Epetra_MultiVector &
>(*vector);
2183 inline Epetra_FEVector &
2191 inline const Epetra_BlockMap &
2202 const Epetra_MpiComm *mpi_comm =
2203 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Map().Comm());
2204 return mpi_comm->Comm();
2209 template <
typename number>
2211 const ::Vector<number> &v,
2248 namespace LinearOperatorImplementation
2261 template <
typename Matrix>
2265 bool omit_zeroing_entries)
2267 v.
reinit(matrix.locally_owned_range_indices(),
2268 matrix.get_mpi_communicator(),
2269 omit_zeroing_entries);
2272 template <
typename Matrix>
2276 bool omit_zeroing_entries)
2278 v.
reinit(matrix.locally_owned_domain_indices(),
2279 matrix.get_mpi_communicator(),
2280 omit_zeroing_entries);
size_type n_elements() const
real_type l1_norm() const
Vector & operator/=(const TrilinosScalar factor)
void compress(VectorOperation::values operation)
Vector(const IndexSet ¶llel_partitioning, const ::Vector< Number > &v, const MPI_Comm communicator=MPI_COMM_WORLD)
void sadd(const TrilinosScalar s, const Vector &V)
TrilinosScalar mean_value() const
void add(const size_type n_elements, const size_type *indices, const TrilinosScalar *values)
void add(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
std::unique_ptr< Epetra_MultiVector > nonlocal_vector
void add(const TrilinosScalar s)
void import_elements(const LinearAlgebra::ReadWriteVector< double > &rwv, const VectorOperation::values operation)
void update_ghost_values() const
MPI_Comm get_mpi_communicator() const
const Epetra_BlockMap & trilinos_partitioner() const
reference operator()(const size_type index)
void import_nonlocal_data_for_fe(const ::TrilinosWrappers::SparseMatrix &matrix, const Vector &vector)
std::unique_ptr< Epetra_FEVector > vector
size_type size() const override
bool in_local_range(const size_type index) const
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
friend class internal::VectorReference
real_type lp_norm(const TrilinosScalar p) const
IndexSet locally_owned_elements() const
virtual void extract_subvector_to(const ArrayView< const size_type > &indices, ArrayView< TrilinosScalar > &elements) const override
Vector & operator-=(const Vector &V)
real_type l2_norm() const
Vector & operator+=(const Vector &V)
const Epetra_MultiVector & trilinos_vector() const
const internal::VectorReference const_reference
void reinit(const Vector &v, const bool omit_zeroing_entries=false, const bool allow_different_maps=false)
bool has_ghost_elements() const
std::pair< size_type, size_type > local_range() const
real_type norm_sqr() const
void extract_subvector_to(const std::vector< size_type > &indices, std::vector< TrilinosScalar > &values) const
bool operator!=(const Vector &v) const
~Vector() override=default
Epetra_FEVector & trilinos_vector()
const_iterator begin() const
real_type linfty_norm() const
internal::VectorReference reference
TrilinosScalar min() const
void set(const std::vector< size_type > &indices, const ::Vector< TrilinosScalar > &values)
void set(const size_type n_elements, const size_type *indices, const TrilinosScalar *values)
void scale(const Vector &scaling_factors)
void equ(const TrilinosScalar a, const Vector &V)
bool operator==(const Vector &v) const
void sadd(const TrilinosScalar s, const TrilinosScalar a, const Vector &V)
Epetra_CombineMode last_action
reference operator[](const size_type index)
const value_type * const_iterator
TrilinosScalar value_type
size_type locally_owned_size() const
TrilinosScalar operator[](const size_type index) const
TrilinosScalar add_and_dot(const TrilinosScalar a, const Vector &V, const Vector &W)
Vector & operator=(const TrilinosScalar s)
Vector & operator*=(const TrilinosScalar factor)
bool is_non_negative() const
::types::global_dof_index size_type
TrilinosScalar operator*(const Vector &vec) const
void add(const std::vector< size_type > &indices, const ::Vector< TrilinosScalar > &values)
void extract_subvector_to(ForwardIterator indices_begin, const ForwardIterator indices_end, OutputIterator values_begin) const
std::size_t memory_consumption() const
void add(const TrilinosScalar a, const Vector &V, const TrilinosScalar b, const Vector &W)
void add(const TrilinosScalar a, const Vector &V)
void set(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
void swap(Vector &u, Vector &v)
TrilinosScalar max() const
Vector & operator=(const ::Vector< Number > &v)
const_iterator end() const
typename numbers::NumberTraits< Number >::real_type real_type
bool has_ghost_elements() const
const value_type * const_iterator
virtual size_type size() const override
size_type locally_owned_size() const
static void reinit_domain_vector(const Matrix &matrix, TrilinosWrappers::MPI::Vector &v, bool omit_zeroing_entries)
static void reinit_range_vector(const Matrix &matrix, TrilinosWrappers::MPI::Vector &v, bool omit_zeroing_entries)
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DeclException0(Exception0)
static ::ExceptionBase & ExcGhostsPresent()
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
static ::ExceptionBase & ExcDifferentParallelPartitioning()
static ::ExceptionBase & ExcAccessToNonLocalElement(size_type arg1, size_type arg2, size_type arg3, size_type arg4)
#define Assert(cond, exc)
#define AssertIsFinite(number)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcTrilinosError(int arg1)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
types::global_dof_index size_type
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &Du)
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
int gid(const Epetra_BlockMap &map, int i)
T sum(const T &t, const MPI_Comm mpi_communicator)
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)
unsigned int global_dof_index