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(
59 vector.vector->Map().NumMyElements(),
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);
113 const MPI_Comm communicator)
123 vector = std::make_unique<Epetra_FEVector>(
132 const MPI_Comm communicator)
135 reinit(local, ghost, communicator,
false);
145 Epetra_Map map(0, 0, Epetra_MpiComm(MPI_COMM_SELF));
148 vector = std::make_unique<Epetra_FEVector>(map);
156 const MPI_Comm communicator,
161 const bool overlapping =
167 vector = std::make_unique<Epetra_FEVector>(map);
197 const bool omit_zeroing_entries,
198 const bool allow_different_maps)
205 if (allow_different_maps ==
false)
211 const Epetra_MpiComm *my_comm =
212 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Comm());
213 const Epetra_MpiComm *v_comm =
214 dynamic_cast<const Epetra_MpiComm *
>(&v.
vector->Comm());
215 const bool same_communicators =
216 my_comm !=
nullptr && v_comm !=
nullptr &&
217 my_comm->DataPtr() == v_comm->DataPtr();
218 if (!same_communicators ||
221 vector = std::make_unique<Epetra_FEVector>(v.
vector->Map());
226 else if (omit_zeroing_entries ==
false)
235 ierr =
vector->PutScalar(0.0);
248 Assert(omit_zeroing_entries ==
false,
250 "It is not possible to exchange data with the "
251 "option 'omit_zeroing_entries' set, which would not write "
257 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
259 const int ierr =
vector->Import(*v.
vector, data_exchange, Insert);
265 const Epetra_MpiComm *comm_ptr =
266 dynamic_cast<const Epetra_MpiComm *
>(&(v.
vector->Comm()));
291 size_type n_elements = 0, added_elements = 0, block_offset = 0;
293 n_elements += v.
block(block).vector->Map().NumMyElements();
294 std::vector<TrilinosWrappers::types::int_type> global_ids(n_elements, -1);
299 v.
block(block).trilinos_partitioner());
300 size_type vector_size = v.
block(block).vector->Map().NumMyElements();
301 for (
size_type i = 0; i < vector_size; ++i)
302 global_ids[added_elements++] = glob_elements[i] + block_offset;
305 block_offset += v.
block(block).size();
309 Epetra_Map new_map(v.
size(),
313 v.
block(0).trilinos_partitioner().Comm());
315 auto actual_vec = std::make_unique<Epetra_FEVector>(new_map);
320 v.
block(block).trilinos_vector().ExtractCopy(entries, 0);
321 entries += v.
block(block).vector->Map().NumMyElements();
324 if (import_data ==
true)
327 *actual_vec)) == v.
size(),
332 Epetra_Import data_exchange(
vector->Map(), actual_vec->Map());
334 const int ierr =
vector->Import(*actual_vec, data_exchange, Insert);
340 vector = std::move(actual_vec);
342 const Epetra_MpiComm *comm_ptr =
343 dynamic_cast<const Epetra_MpiComm *
>(&(
vector->Comm()));
357 const MPI_Comm communicator,
358 const bool vector_writable)
362 if (vector_writable ==
false)
364 IndexSet parallel_partitioner = locally_owned_entries;
368 vector = std::make_unique<Epetra_FEVector>(map);
375 ExcMessage(
"A writable vector must not have ghost entries in "
376 "its parallel partitioning"));
378 if (
vector->Map().SameAs(map) ==
false)
379 vector = std::make_unique<Epetra_FEVector>(map);
382 const int ierr =
vector->PutScalar(0.);
387 IndexSet nonlocal_entries(ghost_entries);
391 Epetra_Map nonlocal_map =
394 std::make_unique<Epetra_MultiVector>(nonlocal_map, 1);
414 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner,
415 const bool make_ghosted,
416 const bool vector_writable)
420 Assert(partitioner->ghost_indices_initialized(),
421 ExcMessage(
"You asked to create a ghosted vector, but the "
422 "partitioner does not provide ghost indices."));
424 this->
reinit(partitioner->locally_owned_range(),
425 partitioner->ghost_indices(),
426 partitioner->get_mpi_communicator(),
431 this->
reinit(partitioner->locally_owned_range(),
432 partitioner->get_mpi_communicator());
442 ExcMessage(
"Vector is not constructed properly."));
446 const Epetra_MpiComm *my_comm =
447 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Comm());
448 const Epetra_MpiComm *v_comm =
449 dynamic_cast<const Epetra_MpiComm *
>(&v.
vector->Comm());
450 const bool same_communicators = my_comm !=
nullptr && v_comm !=
nullptr &&
451 my_comm->DataPtr() == v_comm->DataPtr();
478 if (same_communicators && v.
vector->Map().SameAs(
vector->Map()))
491 (v.
vector->Map().UniqueGIDs() ||
vector->Map().UniqueGIDs()))
524 template <
typename number>
550 "Cannot find exchange information!"));
552 ExcMessage(
"The input vector has overlapping data, "
553 "which is not allowed."));
559 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
560 const int ierr =
vector->Import(*v.
vector, data_exchange, Insert);
575 "Both vectors need to have the same size for import_elements() to work!"));
585 (*this)[idx] = rwv[idx];
590 (*this)[idx] += rwv[idx];
618 "compress() can only be called with VectorOperation add, insert, or unknown"));
628 "The last operation on the Vector and the given last action in the compress() call do not agree!"));
636 const double double_mode = mode;
637 const Epetra_MpiComm *comm_ptr =
645 "Not all processors agree whether the last operation on "
646 "this vector was an addition or a set operation. This will "
647 "prevent the compress() operation from succeeding."));
654 ierr =
vector->GlobalAssemble(mode);
680 if (trilinos_i == -1)
684 vector->Map().NumMyElements(),
686 vector->Map().MaxMyGID()));
689 value = (*vector)[0][trilinos_i];
699 if (allow_different_maps ==
false)
707 # if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0)
708 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
710 vector->Import(*v.
vector, data_exchange, Epetra_AddLocalAlso);
717 Epetra_MultiVector dummy(
vector->Map(), 1,
false);
718 Epetra_Import data_exchange(dummy.Map(), v.
vector->Map());
720 int ierr = dummy.Import(*v.
vector, data_exchange, Insert);
723 ierr =
vector->Update(1.0, dummy, 1.0);
735 if (
vector->Map().NumMyElements() != v.
vector->Map().NumMyElements())
739 for (
size_type i = 0; i < vector_size; ++i)
753 return (!(*
this == v));
765 *eptr = start_ptr +
vector->Map().NumMyElements();
766 unsigned int flag = 0;
779 const Epetra_MpiComm *mpi_comm =
780 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Map().Comm());
783 return num_nonzero == 0;
795 *eptr = start_ptr +
vector->Map().NumMyElements();
796 unsigned int flag = 0;
809 const auto max_n_negative =
811 return max_n_negative == 0;
818 const unsigned int precision,
819 const bool scientific,
820 const bool across)
const
823 boost::io::ios_flags_saver restore_flags(out);
826 out.precision(precision);
828 out.setf(std::ios::scientific, std::ios::floatfield);
830 out.setf(std::ios::fixed, std::ios::floatfield);
833 if (
size() != vector_size)
835 auto global_id = [&](
const size_type index) {
838 out <<
"size:" <<
size()
839 <<
" locally_owned_size:" <<
vector->Map().NumMyElements() <<
" :"
841 for (
size_type i = 0; i < vector_size; ++i)
842 out <<
"[" << global_id(i) <<
"]: " << (*(
vector))[0][i]
848 int leading_dimension;
849 int ierr =
vector->ExtractView(&val, &leading_dimension);
854 out <<
static_cast<double>(val[i]) <<
' ';
857 out <<
static_cast<double>(val[i]) << std::endl;
887 return sizeof(*this) +
888 this->
vector->Map().NumMyElements() *
894 # include "trilinos_vector.inst"
virtual size_type size() const override
unsigned int n_blocks() const
BlockType & block(const unsigned int i)
bool is_ascending_and_one_to_one(const MPI_Comm communicator) const
size_type n_elements() const
void set_size(const size_type size)
Epetra_Map make_trilinos_map(const MPI_Comm communicator=MPI_COMM_WORLD, const bool overlapping=false) const
void subtract_set(const IndexSet &other)
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
size_type size() const override
const IndexSet & get_stored_elements() const
types::global_dof_index size_type
void compress(VectorOperation::values operation)
const Epetra_BlockMap & trilinos_partitioner() const
void add(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
std::unique_ptr< Epetra_MultiVector > nonlocal_vector
void import_elements(const LinearAlgebra::ReadWriteVector< double > &rwv, const VectorOperation::values operation)
MPI_Comm get_mpi_communicator() 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
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
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)
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)
T max(const T &t, const MPI_Comm mpi_communicator)
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm mpi_communicator)
const Epetra_Comm & comm_self()