18 #ifdef DEAL_II_WITH_TRILINOS
28 # include <boost/io/ios_state.hpp>
30 # include <Epetra_Export.h>
31 # include <Epetra_Import.h>
32 # include <Epetra_Vector.h>
54 vector.vector->Map().LID(
60 vector.vector->Map().MinMyGID(),
61 vector.vector->Map().MaxMyGID()));
64 return (*(vector.vector))[0][local_index];
76 , vector(new Epetra_FEVector(
83 const MPI_Comm &communicator)
86 reinit(parallel_partitioning, communicator);
112 const MPI_Comm &communicator)
122 vector = std::make_unique<Epetra_FEVector>(
131 const MPI_Comm &communicator)
134 reinit(local, ghost, communicator,
false);
144 Epetra_Map map(0, 0, Epetra_MpiComm(MPI_COMM_SELF));
147 vector = std::make_unique<Epetra_FEVector>(map);
155 const MPI_Comm &communicator,
160 const bool overlapping =
166 vector = std::make_unique<Epetra_FEVector>(map);
196 const bool omit_zeroing_entries,
197 const bool allow_different_maps)
204 if (allow_different_maps ==
false)
210 const Epetra_MpiComm *my_comm =
211 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Comm());
212 const Epetra_MpiComm *v_comm =
213 dynamic_cast<const Epetra_MpiComm *
>(&v.
vector->Comm());
214 const bool same_communicators =
215 my_comm !=
nullptr && v_comm !=
nullptr &&
216 my_comm->DataPtr() == v_comm->DataPtr();
217 if (!same_communicators ||
220 vector = std::make_unique<Epetra_FEVector>(v.
vector->Map());
225 else if (omit_zeroing_entries ==
false)
234 ierr =
vector->PutScalar(0.0);
247 Assert(omit_zeroing_entries ==
false,
249 "It is not possible to exchange data with the "
250 "option 'omit_zeroing_entries' set, which would not write "
256 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
258 const int ierr =
vector->Import(*v.
vector, data_exchange, Insert);
264 const Epetra_MpiComm *comm_ptr =
265 dynamic_cast<const Epetra_MpiComm *
>(&(v.
vector->Comm()));
290 size_type n_elements = 0, added_elements = 0, block_offset = 0;
292 n_elements += v.
block(block).local_size();
293 std::vector<TrilinosWrappers::types::int_type> global_ids(n_elements, -1);
298 v.
block(block).trilinos_partitioner());
300 global_ids[added_elements++] = glob_elements[i] + block_offset;
303 block_offset += v.
block(block).size();
307 Epetra_Map new_map(v.
size(),
311 v.
block(0).trilinos_partitioner().Comm());
313 auto actual_vec = std::make_unique<Epetra_FEVector>(new_map);
318 v.
block(block).trilinos_vector().ExtractCopy(entries, 0);
319 entries += v.
block(block).local_size();
322 if (import_data ==
true)
325 *actual_vec)) == v.
size(),
330 Epetra_Import data_exchange(
vector->Map(), actual_vec->Map());
332 const int ierr =
vector->Import(*actual_vec, data_exchange, Insert);
338 vector = std::move(actual_vec);
340 const Epetra_MpiComm *comm_ptr =
341 dynamic_cast<const Epetra_MpiComm *
>(&(
vector->Comm()));
355 const MPI_Comm &communicator,
356 const bool vector_writable)
360 if (vector_writable ==
false)
362 IndexSet parallel_partitioner = locally_owned_entries;
366 vector = std::make_unique<Epetra_FEVector>(map);
373 ExcMessage(
"A writable vector must not have ghost entries in "
374 "its parallel partitioning"));
376 if (
vector->Map().SameAs(map) ==
false)
377 vector = std::make_unique<Epetra_FEVector>(map);
380 const int ierr =
vector->PutScalar(0.);
385 IndexSet nonlocal_entries(ghost_entries);
389 Epetra_Map nonlocal_map =
392 std::make_unique<Epetra_MultiVector>(nonlocal_map, 1);
412 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner,
413 const bool vector_writable)
415 this->
reinit(partitioner->locally_owned_range(),
416 partitioner->ghost_indices(),
417 partitioner->get_mpi_communicator(),
427 ExcMessage(
"Vector is not constructed properly."));
431 const Epetra_MpiComm *my_comm =
432 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Comm());
433 const Epetra_MpiComm *v_comm =
434 dynamic_cast<const Epetra_MpiComm *
>(&v.
vector->Comm());
435 const bool same_communicators = my_comm !=
nullptr && v_comm !=
nullptr &&
436 my_comm->DataPtr() == v_comm->DataPtr();
463 if (same_communicators && v.
vector->Map().SameAs(
vector->Map()))
476 (v.
vector->Map().UniqueGIDs() ||
vector->Map().UniqueGIDs()))
508 template <
typename number>
534 "Cannot find exchange information!"));
536 ExcMessage(
"The input vector has overlapping data, "
537 "which is not allowed."));
543 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
544 const int ierr =
vector->Import(*v.
vector, data_exchange, Insert);
559 "Both vectors need to have the same size for import() to work!"));
569 (*this)[idx] = rwv[idx];
574 (*this)[idx] += rwv[idx];
602 "compress() can only be called with VectorOperation add, insert, or unknown"));
612 "The last operation on the Vector and the given last action in the compress() call do not agree!"));
620 const double double_mode = mode;
621 const Epetra_MpiComm *comm_ptr =
629 "Not all processors agree whether the last operation on "
630 "this vector was an addition or a set operation. This will "
631 "prevent the compress() operation from succeeding."));
638 ierr =
vector->GlobalAssemble(mode);
664 if (trilinos_i == -1)
670 vector->Map().MaxMyGID()));
673 value = (*vector)[0][trilinos_i];
683 if (allow_different_maps ==
false)
691 # if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0)
692 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
694 vector->Import(*v.
vector, data_exchange, Epetra_AddLocalAlso);
701 Epetra_MultiVector dummy(
vector->Map(), 1,
false);
702 Epetra_Import data_exchange(dummy.Map(), v.
vector->Map());
704 int ierr = dummy.Import(*v.
vector, data_exchange, Insert);
707 ierr =
vector->Update(1.0, dummy, 1.0);
724 if ((*(v.
vector))[0][i] != (*vector)[0][i])
737 return (!(*
this == v));
749 unsigned int flag = 0;
762 const Epetra_MpiComm *mpi_comm =
763 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Map().Comm());
766 return num_nonzero == 0;
778 unsigned int flag = 0;
791 const auto max_n_negative =
793 return max_n_negative == 0;
800 const unsigned int precision,
801 const bool scientific,
802 const bool across)
const
805 boost::io::ios_flags_saver restore_flags(out);
808 out.precision(precision);
810 out.setf(std::ios::scientific, std::ios::floatfield);
812 out.setf(std::ios::fixed, std::ios::floatfield);
816 auto global_id = [&](
const size_type index) {
819 out <<
"size:" <<
size() <<
" local_size:" <<
local_size() <<
" :"
822 out <<
"[" << global_id(i) <<
"]: " << (*(
vector))[0][i]
828 int leading_dimension;
829 int ierr =
vector->ExtractView(&val, &leading_dimension);
834 out <<
static_cast<double>(val[i]) <<
' ';
837 out <<
static_cast<double>(val[i]) << std::endl;
867 return sizeof(*this) +
874 # include "trilinos_vector.inst"
unsigned int n_blocks() const
BlockType & block(const unsigned int i)
size_type n_elements() const
void set_size(const size_type size)
void subtract_set(const IndexSet &other)
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
bool is_ascending_and_one_to_one(const MPI_Comm &communicator) const
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
const IndexSet & get_stored_elements() const
types::global_dof_index size_type
const Epetra_BlockMap & trilinos_partitioner() const
void add(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
const MPI_Comm & get_mpi_communicator() const
std::unique_ptr< Epetra_MultiVector > nonlocal_vector
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 local_size() const
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
IndexSet locally_owned_elements() const
void reinit(const Vector &v, const bool omit_zeroing_entries=false, const bool allow_different_maps=false)
bool has_ghost_elements() const
void compress(::VectorOperation::values operation)
bool operator!=(const Vector &v) const
bool operator==(const Vector &v) const
Epetra_CombineMode last_action
std::pair< size_type, size_type > local_range() const
Vector & operator=(const TrilinosScalar s)
void import(const LinearAlgebra::ReadWriteVector< double > &rwv, const VectorOperation::values operation)
bool is_non_negative() const
::types::global_dof_index size_type
std::size_t memory_consumption() const
const Epetra_CrsMatrix & trilinos_matrix() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcGhostsPresent()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcAccessToNonLocalElement(size_type arg1, size_type arg2, size_type arg3, size_type arg4)
#define AssertThrow(cond, exc)
void swap(MemorySpaceData< T, MemorySpace > &u, MemorySpaceData< T, MemorySpace > &v)
void swap(BlockVector &u, BlockVector &v)
TrilinosWrappers::types::int_type global_length(const Epetra_MultiVector &vector)
int gid(const Epetra_BlockMap &map, int i)
TrilinosWrappers::types::int_type * my_global_elements(const Epetra_BlockMap &map)
TrilinosWrappers::types::int64_type n_global_elements(const Epetra_BlockMap &map)
T sum(const T &t, const MPI_Comm &mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator)
T max(const T &t, const MPI_Comm &mpi_communicator)
const Epetra_Comm & comm_self()