16 #include <deal.II/lac/trilinos_vector.h> 18 #ifdef DEAL_II_WITH_TRILINOS 20 # include <deal.II/base/mpi.h> 21 # include <deal.II/base/std_cxx14/memory.h> 23 # include <deal.II/lac/trilinos_index_access.h> 24 # include <deal.II/lac/trilinos_parallel_block_vector.h> 25 # include <deal.II/lac/trilinos_sparse_matrix.h> 27 # include <boost/io/ios_state.hpp> 29 # include <Epetra_Export.h> 30 # include <Epetra_Import.h> 31 # include <Epetra_Vector.h> 36 DEAL_II_NAMESPACE_OPEN
42 VectorReference::operator TrilinosScalar()
const 50 const TrilinosWrappers::types::int_type local_index =
51 vector.vector->Map().LID(
52 static_cast<TrilinosWrappers::types::int_type>(index));
57 vector.vector->Map().MinMyGID(),
58 vector.vector->Map().MaxMyGID()));
61 return (*(vector.vector))[0][local_index];
72 , vector(new Epetra_FEVector(
73 Epetra_Map(0, 0, 0,
Utilities::Trilinos::comm_self())))
79 const MPI_Comm &communicator)
82 reinit(parallel_partitioning, communicator);
91 vector = std_cxx14::make_unique<Epetra_FEVector>(*v.
vector);
108 const MPI_Comm &communicator)
112 static_cast<size_type
>(
118 vector = std_cxx14::make_unique<Epetra_FEVector>(
127 const MPI_Comm &communicator)
130 reinit(local, ghost, communicator,
false);
140 # ifdef DEAL_II_WITH_MPI 141 Epetra_Map map(0, 0, Epetra_MpiComm(MPI_COMM_SELF));
143 Epetra_Map map(0, 0, Epetra_SerialComm());
147 vector = std_cxx14::make_unique<Epetra_FEVector>(map);
155 const MPI_Comm &communicator,
160 const bool overlapping =
166 vector = std_cxx14::make_unique<Epetra_FEVector>(map);
183 const size_type n_elements_global =
196 const bool omit_zeroing_entries,
197 const bool allow_different_maps)
204 if (allow_different_maps ==
false)
210 # ifdef DEAL_II_WITH_MPI 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();
219 const bool same_communicators =
true;
221 if (!same_communicators ||
224 vector = std_cxx14::make_unique<Epetra_FEVector>(v.
vector->Map());
229 else if (omit_zeroing_entries ==
false)
238 ierr =
vector->PutScalar(0.0);
251 Assert(omit_zeroing_entries ==
false,
253 "It is not possible to exchange data with the " 254 "option 'omit_zeroing_entries' set, which would not write " 260 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
262 const int ierr =
vector->Import(*v.
vector, data_exchange, Insert);
267 # if defined(DEBUG) && defined(DEAL_II_WITH_MPI) 268 const Epetra_MpiComm *comm_ptr =
269 dynamic_cast<const Epetra_MpiComm *
>(&(v.
vector->Comm()));
271 const size_type n_elements_global =
294 size_type n_elements = 0, added_elements = 0, block_offset = 0;
295 for (size_type block = 0; block < v.
n_blocks(); ++block)
296 n_elements += v.
block(block).local_size();
297 std::vector<TrilinosWrappers::types::int_type> global_ids(n_elements, -1);
298 for (size_type block = 0; block < v.
n_blocks(); ++block)
300 TrilinosWrappers::types::int_type *glob_elements =
302 v.
block(block).vector_partitioner());
303 for (size_type i = 0; i < v.
block(block).local_size(); ++i)
304 global_ids[added_elements++] = glob_elements[i] + block_offset;
307 block_offset += v.
block(block).size();
311 Epetra_Map new_map(v.
size(),
315 v.
block(0).vector_partitioner().Comm());
317 auto actual_vec = std_cxx14::make_unique<Epetra_FEVector>(new_map);
319 TrilinosScalar *entries = (*actual_vec)[0];
320 for (size_type block = 0; block < v.
n_blocks(); ++block)
322 v.
block(block).trilinos_vector().ExtractCopy(entries, 0);
323 entries += v.
block(block).local_size();
326 if (import_data ==
true)
329 *actual_vec)) == v.
size(),
334 Epetra_Import data_exchange(
vector->Map(), actual_vec->Map());
336 const int ierr =
vector->Import(*actual_vec, data_exchange, Insert);
342 vector = std::move(actual_vec);
343 # if defined(DEBUG) && defined(DEAL_II_WITH_MPI) 344 const Epetra_MpiComm *comm_ptr =
345 dynamic_cast<const Epetra_MpiComm *
>(&(
vector->Comm()));
347 const size_type n_elements_global =
359 const MPI_Comm &communicator,
360 const bool vector_writable)
364 if (vector_writable ==
false)
366 IndexSet parallel_partitioner = locally_owned_entries;
370 vector = std_cxx14::make_unique<Epetra_FEVector>(map);
377 ExcMessage(
"A writable vector must not have ghost entries in " 378 "its parallel partitioning"));
380 if (
vector->Map().SameAs(map) ==
false)
381 vector = std_cxx14::make_unique<Epetra_FEVector>(map);
384 const int ierr =
vector->PutScalar(0.);
389 IndexSet nonlocal_entries(ghost_entries);
393 Epetra_Map nonlocal_map =
396 std_cxx14::make_unique<Epetra_MultiVector>(nonlocal_map, 1);
405 const size_type n_elements_global =
418 ExcMessage(
"Vector is not constructed properly."));
422 # ifdef DEAL_II_WITH_MPI 423 const Epetra_MpiComm *my_comm =
424 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Comm());
425 const Epetra_MpiComm *v_comm =
426 dynamic_cast<const Epetra_MpiComm *
>(&v.
vector->Comm());
427 const bool same_communicators = my_comm !=
nullptr && v_comm !=
nullptr &&
428 my_comm->DataPtr() == v_comm->DataPtr();
450 const bool same_communicators =
true;
457 if (same_communicators && v.
vector->Map().SameAs(
vector->Map()))
470 (v.
vector->Map().UniqueGIDs() ||
vector->Map().UniqueGIDs()))
478 vector = std_cxx14::make_unique<Epetra_FEVector>(*v.
vector);
503 template <
typename number>
515 for (size_type i = local_range.first; i < local_range.second; ++i)
516 (*vector)[0][i - local_range.first] = v(i);
529 "Cannot find exchange information!"));
531 ExcMessage(
"The input vector has overlapping data, " 532 "which is not allowed."));
538 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
539 const int ierr =
vector->Import(*v.
vector, data_exchange, Insert);
567 "compress() can only be called with VectorOperation add, insert, or unknown"));
577 "The last operation on the Vector and the given last action in the compress() call do not agree!"));
582 # ifdef DEAL_II_WITH_MPI 586 double double_mode = mode;
587 const Epetra_MpiComm *comm_ptr =
594 "Not all processors agree whether the last operation on " 595 "this vector was an addition or a set operation. This will " 596 "prevent the compress() operation from succeeding."));
604 ierr =
vector->GlobalAssemble(mode);
624 TrilinosWrappers::types::int_type trilinos_i =
vector->Map().LID(
625 static_cast<TrilinosWrappers::types::int_type>(index));
626 TrilinosScalar value = 0.;
630 if (trilinos_i == -1)
636 vector->Map().MaxMyGID()));
639 value = (*vector)[0][trilinos_i];
649 if (allow_different_maps ==
false)
657 # if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0) 658 Epetra_Import data_exchange(
vector->Map(), v.
vector->Map());
660 vector->Import(*v.
vector, data_exchange, Epetra_AddLocalAlso);
667 Epetra_MultiVector dummy(
vector->Map(), 1,
false);
668 Epetra_Import data_exchange(dummy.Map(), v.
vector->Map());
670 int ierr = dummy.Import(*v.
vector, data_exchange, Insert);
673 ierr =
vector->Update(1.0, dummy, 1.0);
690 if ((*(v.
vector))[0][i] != (*vector)[0][i])
703 return (!(*
this == v));
713 TrilinosScalar * start_ptr = (*vector)[0];
714 const TrilinosScalar *ptr = start_ptr, *eptr = start_ptr +
local_size();
715 unsigned int flag = 0;
726 # ifdef DEAL_II_WITH_MPI 729 const Epetra_MpiComm *mpi_comm =
730 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Map().Comm());
733 return num_nonzero == 0;
744 # ifdef DEAL_II_WITH_MPI 754 TrilinosScalar *start_ptr;
755 int leading_dimension;
756 int ierr =
vector->ExtractView(&start_ptr, &leading_dimension);
763 const TrilinosScalar *ptr = start_ptr, *eptr = start_ptr +
size();
782 const unsigned int precision,
783 const bool scientific,
784 const bool across)
const 787 boost::io::ios_flags_saver restore_flags(out);
790 out.precision(precision);
792 out.setf(std::ios::scientific, std::ios::floatfield);
794 out.setf(std::ios::fixed, std::ios::floatfield);
798 auto global_id = [&](
const size_type index) {
800 static_cast<TrilinosWrappers::types::int_type
>(index));
802 out <<
"size:" <<
size() <<
" local_size:" <<
local_size() <<
" :" 805 out <<
"[" << global_id(i) <<
"]: " << (*(
vector))[0][i]
811 int leading_dimension;
812 int ierr =
vector->ExtractView(&val, &leading_dimension);
816 for (size_type i = 0; i <
size(); ++i)
817 out << static_cast<double>(val[i]) <<
' ';
819 for (size_type i = 0; i <
size(); ++i)
820 out << static_cast<double>(val[i]) << std::endl;
850 return sizeof(*this) +
852 (
sizeof(double) +
sizeof(TrilinosWrappers::types::int_type));
861 # include "trilinos_vector.inst" 865 DEAL_II_NAMESPACE_CLOSE
867 #endif // DEAL_II_WITH_TRILINOS Vector & operator=(const TrilinosScalar s)
std::unique_ptr< Epetra_MultiVector > nonlocal_vector
bool is_non_negative() const
const Epetra_CrsMatrix & trilinos_matrix() const
static::ExceptionBase & ExcIO()
TrilinosWrappers::types::int_type global_length(const Epetra_MultiVector &vector)
bool is_ascending_and_one_to_one(const MPI_Comm &communicator) const
void reinit(const Vector &v, const bool omit_zeroing_entries=false, const bool allow_different_maps=false)
void add(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
size_type n_elements() const
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
std::pair< size_type, size_type > local_range() const
unsigned int n_blocks() const
#define AssertThrow(cond, exc)
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static::ExceptionBase & ExcMessage(std::string arg1)
bool operator!=(const Vector &v) const
reference operator()(const size_type index)
T sum(const T &t, const MPI_Comm &mpi_communicator)
TrilinosWrappers::types::int_type * my_global_elements(const Epetra_BlockMap &map)
void subtract_set(const IndexSet &other)
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
std::unique_ptr< Epetra_FEVector > vector
void compress(::VectorOperation::values operation)
const Epetra_Map & vector_partitioner() const
void import_nonlocal_data_for_fe(const ::TrilinosWrappers::SparseMatrix &matrix, const Vector &vector)
bool has_ghost_elements() const
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
static::ExceptionBase & ExcGhostsPresent()
TrilinosWrappers::types::int_type n_global_elements(const Epetra_BlockMap &map)
std::size_t memory_consumption() const
void swap(Vector< Number > &u, Vector< Number > &v)
void set_size(const size_type size)
bool operator==(const Vector &v) const
size_type local_size() const
static::ExceptionBase & ExcAccessToNonLocalElement(size_type arg1, size_type arg2, size_type arg3, size_type arg4)
static::ExceptionBase & ExcNotImplemented()
Epetra_CombineMode last_action
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator)
static::ExceptionBase & ExcTrilinosError(int arg1)
BlockType & block(const unsigned int i)
static::ExceptionBase & ExcInternalError()