18 #ifdef DEAL_II_WITH_TRILINOS
26 # include <boost/io/ios_state.hpp>
28 # include <Epetra_Import.h>
29 # include <Epetra_Map.h>
30 # include <Epetra_MpiComm.h>
39 namespace EpetraWrappers
45 # ifdef DEAL_II_HAVE_CXX20
46 static_assert(concepts::is_vector_space_vector<Vector>);
51 , vector(new Epetra_FEVector(
59 , vector(new Epetra_FEVector(
V.trilinos_vector()))
65 const MPI_Comm communicator)
67 , vector(new Epetra_FEVector(
68 parallel_partitioner.make_trilinos_map(communicator, false)))
75 const MPI_Comm communicator,
76 const bool omit_zeroing_entries)
78 Epetra_Map input_map =
80 if (
vector->Map().SameAs(input_map) ==
false)
81 vector = std::make_unique<Epetra_FEVector>(input_map);
82 else if (omit_zeroing_entries ==
false)
84 const int ierr =
vector->PutScalar(0.);
95 reinit(
V.locally_owned_elements(),
96 V.get_mpi_communicator(),
97 omit_zeroing_entries);
109 const auto &map =
vector.Map();
111 for (
unsigned int i = 0; i < indices.
size(); ++i)
114 const auto trilinos_i =
116 elements[i] =
vector[0][trilinos_i];
129 if (
vector->Map().SameAs(
V.trilinos_vector().Map()))
133 if (
size() ==
V.size())
135 Epetra_Import data_exchange(
vector->Map(),
136 V.trilinos_vector().Map());
139 vector->Import(
V.trilinos_vector(), data_exchange, Insert);
144 vector = std::make_unique<Epetra_FEVector>(
V.trilinos_vector());
157 const int ierr =
vector->PutScalar(s);
170 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
171 &communication_pattern)
175 if (communication_pattern ==
nullptr)
181 V.get_stored_elements().size()) ||
185 V.get_stored_elements(),
186 dynamic_cast<const Epetra_MpiComm &
>(
vector->Comm()).Comm());
192 std::dynamic_pointer_cast<const CommunicationPattern>(
193 communication_pattern);
197 std::string(
"The communication pattern is not of type ") +
198 "LinearAlgebra::EpetraWrappers::CommunicationPattern."));
204 Epetra_FEVector source_vector(import_map.TargetMap());
205 double *
values = source_vector.Values();
209 vector->Export(source_vector, import_map, Insert);
211 vector->Export(source_vector, import_map, Add);
213 vector->Export(source_vector, import_map, Epetra_Max);
215 vector->Export(source_vector, import_map, Epetra_Min);
238 *
this *= 1. / factor;
249 if (
vector->Map().SameAs(
V.trilinos_vector().Map()))
251 const int ierr =
vector->Update(1.,
V.trilinos_vector(), 1.);
260 # if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0)
261 Epetra_Import data_exchange(
vector->Map(),
V.trilinos_vector().Map());
262 const int ierr =
vector->Import(
V.trilinos_vector(),
264 Epetra_AddLocalAlso);
271 Epetra_MultiVector dummy(
vector->Map(), 1,
false);
272 Epetra_Import data_exchange(dummy.Map(),
V.trilinos_vector().Map());
274 int ierr = dummy.Import(
V.trilinos_vector(), data_exchange, Insert);
277 ierr =
vector->Update(1.0, dummy, 1.0);
307 const int ierr =
vector->Dot(
V.trilinos_vector(), &result);
320 const unsigned local_size(
vector->MyLength());
321 for (
unsigned int i = 0; i < local_size; ++i)
334 const int ierr =
vector->Update(a,
V.trilinos_vector(), 1.);
391 if (
vector->Map().SameAs(
V.trilinos_vector().Map()) ==
false)
392 this->
sadd(0., a,
V);
396 int ierr =
vector->Update(a,
V.trilinos_vector(), 0.);
409 double *start_ptr = (*vector)[0];
410 const double *ptr = start_ptr, *eptr = start_ptr +
vector->MyLength();
411 unsigned int flag = 0;
423 const Epetra_MpiComm *mpi_comm =
424 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Map().Comm());
428 return num_nonzero == 0;
499 # ifndef DEAL_II_WITH_64BIT_INDICES
500 return vector->GlobalLength();
502 return vector->GlobalLength64();
511 return vector->MyLength();
519 const Epetra_MpiComm *epetra_comm =
520 dynamic_cast<const Epetra_MpiComm *
>(&(
vector->Comm()));
522 return epetra_comm->GetMpiComm();
533 if (
vector->Map().LinearMap())
535 # ifndef DEAL_II_WITH_64BIT_INDICES
539 vector->Map().MaxMyGID64() + 1);
542 else if (
vector->Map().NumMyElements() > 0)
545 # ifndef DEAL_II_WITH_64BIT_INDICES
546 unsigned int *vector_indices =
547 reinterpret_cast<unsigned int *
>(
vector->Map().MyGlobalElements());
552 is.
add_indices(vector_indices, vector_indices + n_indices);
566 const Epetra_FEVector &
584 const unsigned int precision,
585 const bool scientific,
586 const bool across)
const
589 boost::io::ios_flags_saver restore_flags(out);
594 int leading_dimension;
595 int ierr =
vector->ExtractView(&val, &leading_dimension);
599 out.precision(precision);
601 out.setf(std::ios::scientific, std::ios::floatfield);
603 out.setf(std::ios::fixed, std::ios::floatfield);
606 for (
int i = 0; i <
vector->MyLength(); ++i)
607 out << val[i] <<
' ';
609 for (
int i = 0; i <
vector->MyLength(); ++i)
610 out << val[i] << std::endl;
623 return sizeof(*this) +
632 const MPI_Comm mpi_comm)
Epetra_Map make_trilinos_map(const MPI_Comm communicator=MPI_COMM_WORLD, const bool overlapping=false) const
void add_range(const size_type begin, const size_type end)
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
Vector & operator-=(const Vector &V)
void equ(const double a, const Vector &V)
Vector & operator=(const Vector &V)
MPI_Comm get_mpi_communicator() const
void compress(const VectorOperation::values operation)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
std::unique_ptr< Epetra_FEVector > vector
const Epetra_FEVector & trilinos_vector() const
void import_elements(const ReadWriteVector< double > &V, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
void sadd(const double s, const double a, const Vector &V)
Vector & operator/=(const double factor)
virtual size_type size() const override
types::global_dof_index size_type
size_type locally_owned_size() const
::IndexSet locally_owned_elements() const
Vector & operator+=(const Vector &V)
virtual void extract_subvector_to(const ArrayView< const types::global_dof_index > &indices, ArrayView< double > &elements) const override
void scale(const Vector &scaling_factors)
void create_epetra_comm_pattern(const IndexSet &source_index_set, const MPI_Comm mpi_comm)
double mean_value() const
Vector & operator*=(const double factor)
::IndexSet source_stored_elements
double operator*(const Vector &V) const
std::shared_ptr< const CommunicationPattern > epetra_comm_pattern
std::size_t memory_consumption() const
void reinit(const IndexSet ¶llel_partitioner, const MPI_Comm communicator, const bool omit_zeroing_entries=false)
double linfty_norm() const
double add_and_dot(const double a, const Vector &V, const Vector &W)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcZero()
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDifferentParallelPartitioning()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertIsFinite(number)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim >>> &Du)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
T sum(const T &t, const MPI_Comm mpi_communicator)
const Epetra_Comm & comm_self()
void copy(const T *begin, const T *end, U *dest)