16 #include <deal.II/base/std_cxx14/memory.h> 18 #include <deal.II/lac/trilinos_epetra_vector.h> 20 #ifdef DEAL_II_WITH_TRILINOS 22 # ifdef DEAL_II_WITH_MPI 24 # include <deal.II/base/index_set.h> 26 # include <deal.II/lac/read_write_vector.h> 28 # include <boost/io/ios_state.hpp> 30 # include <Epetra_Import.h> 31 # include <Epetra_Map.h> 32 # include <Epetra_MpiComm.h> 37 DEAL_II_NAMESPACE_OPEN
41 namespace EpetraWrappers
44 : vector(new Epetra_FEVector(
45 Epetra_Map(0, 0, 0,
Utilities::Trilinos::comm_self())))
58 const MPI_Comm &communicator)
59 :
vector(new Epetra_FEVector(
60 parallel_partitioner.make_trilinos_map(communicator, false)))
67 const MPI_Comm &communicator,
68 const bool omit_zeroing_entries)
70 Epetra_Map input_map =
72 if (
vector->Map().SameAs(input_map) ==
false)
73 vector = std_cxx14::make_unique<Epetra_FEVector>(input_map);
74 else if (omit_zeroing_entries ==
false)
76 const int ierr =
vector->PutScalar(0.);
86 const bool omit_zeroing_entries)
89 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
97 omit_zeroing_entries);
115 Epetra_Import data_exchange(
vector->Map(),
138 const int ierr =
vector->PutScalar(s);
151 std::shared_ptr<const CommunicationPatternBase> communication_pattern)
155 if (communication_pattern ==
nullptr)
166 dynamic_cast<const Epetra_MpiComm &
>(
vector->Comm()).Comm());
173 communication_pattern);
177 std::string(
"The communication pattern is not of type ") +
178 "LinearAlgebra::EpetraWrappers::CommunicationPattern."));
184 Epetra_FEVector source_vector(
import.TargetMap());
185 double * values = source_vector.Values();
186 std::copy(V.
begin(), V.
end(), values);
189 vector->Export(source_vector,
import, Insert);
191 vector->Export(source_vector,
import, Add);
193 vector->Export(source_vector,
import, Epetra_Max);
195 vector->Export(source_vector,
import, Epetra_Min);
218 *
this *= 1. / factor;
229 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
246 # if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0) 247 Epetra_Import data_exchange(
vector->Map(),
251 Epetra_AddLocalAlso);
258 Epetra_MultiVector dummy(
vector->Map(), 1,
false);
259 Epetra_Import data_exchange(dummy.Map(),
266 ierr =
vector->Update(1.0, dummy, 1.0);
290 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
314 const unsigned local_size(
vector->MyLength());
315 for (
unsigned int i = 0; i < local_size; ++i)
325 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
348 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
351 Assert(dynamic_cast<const Vector *>(&W) !=
nullptr,
365 const int ierr =
vector->Update(
379 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
396 Assert(dynamic_cast<const Vector *>(&scaling_factors) !=
nullptr,
400 const Vector &down_scaling_factors =
401 dynamic_cast<const Vector &
>(scaling_factors);
405 const int ierr =
vector->Multiply(1.0,
419 Assert(dynamic_cast<const Vector *>(&V) !=
nullptr,
426 this->
sadd(0., a, V);
443 double * start_ptr = (*vector)[0];
444 const double *ptr = start_ptr, *eptr = start_ptr +
vector->MyLength();
445 unsigned int flag = 0;
457 const Epetra_MpiComm *mpi_comm =
458 dynamic_cast<const Epetra_MpiComm *
>(&
vector->Map().Comm());
462 return num_nonzero == 0;
472 int ierr =
vector->MeanValue(&mean_value);
485 int ierr =
vector->Norm1(&norm);
498 int ierr =
vector->Norm2(&norm);
511 int ierr =
vector->NormInf(&norm);
535 # ifndef DEAL_II_WITH_64BIT_INDICES 536 return vector->GlobalLength();
538 return vector->GlobalLength64();
547 const Epetra_MpiComm *epetra_comm =
548 dynamic_cast<const Epetra_MpiComm *
>(&(
vector->Comm()));
550 return epetra_comm->GetMpiComm();
561 if (
vector->Map().LinearMap())
563 # ifndef DEAL_II_WITH_64BIT_INDICES 567 vector->Map().MaxMyGID64() + 1);
570 else if (
vector->Map().NumMyElements() > 0)
572 const size_type n_indices =
vector->Map().NumMyElements();
573 # ifndef DEAL_II_WITH_64BIT_INDICES 574 unsigned int *vector_indices =
575 (
unsigned int *)
vector->Map().MyGlobalElements();
577 size_type *vector_indices =
578 (size_type *)
vector->Map().MyGlobalElements64();
580 is.
add_indices(vector_indices, vector_indices + n_indices);
589 const Epetra_FEVector &
607 const unsigned int precision,
608 const bool scientific,
609 const bool across)
const 612 boost::io::ios_flags_saver restore_flags(out);
617 int leading_dimension;
618 int ierr =
vector->ExtractView(&val, &leading_dimension);
622 out.precision(precision);
624 out.setf(std::ios::scientific, std::ios::floatfield);
626 out.setf(std::ios::fixed, std::ios::floatfield);
629 for (
int i = 0; i <
vector->MyLength(); ++i)
630 out << val[i] <<
' ';
632 for (
int i = 0; i <
vector->MyLength(); ++i)
633 out << val[i] << std::endl;
646 return sizeof(*this) +
648 (
sizeof(double) +
sizeof(TrilinosWrappers::types::int_type));
655 const MPI_Comm &mpi_comm)
666 DEAL_II_NAMESPACE_CLOSE
virtual double l1_norm() const override
static::ExceptionBase & ExcDifferentParallelPartitioning()
virtual bool all_zero() const override
static::ExceptionBase & ExcIO()
const Epetra_FEVector & trilinos_vector() const
virtual ::IndexSet locally_owned_elements() const override
virtual double l2_norm() const override
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
virtual size_type size() const override
virtual void scale(const VectorSpaceVector< double > &scaling_factors) override
#define AssertThrow(cond, exc)
virtual std::size_t memory_consumption() const override
std::unique_ptr< Epetra_FEVector > vector
virtual void equ(const double a, const VectorSpaceVector< double > &V) override
static::ExceptionBase & ExcMessage(std::string arg1)
T sum(const T &t, const MPI_Comm &mpi_communicator)
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void create_epetra_comm_pattern(const IndexSet &source_index_set, const MPI_Comm &mpi_comm)
virtual Vector & operator*=(const double factor) override
virtual void import(const ReadWriteVector< double > &V, VectorOperation::values operation, std::shared_ptr< const CommunicationPatternBase > communication_pattern=std::shared_ptr< const CommunicationPatternBase >()) override
Vector & operator=(const Vector &V)
static::ExceptionBase & ExcTrilinosError(int arg1)
static::ExceptionBase & ExcVectorTypeNotCompatible()
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
void reinit(const IndexSet ¶llel_partitioner, const MPI_Comm &communicator, const bool omit_zeroing_entries=false)
std::shared_ptr< const CommunicationPattern > epetra_comm_pattern
virtual void sadd(const double s, const double a, const VectorSpaceVector< double > &V) override
void add_range(const size_type begin, const size_type end)
virtual Vector & operator+=(const VectorSpaceVector< double > &V) override
const IndexSet & get_stored_elements() const
virtual Vector & operator/=(const double factor) override
virtual void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const override
virtual double mean_value() const override
MPI_Comm get_mpi_communicator() const
virtual Vector & operator-=(const VectorSpaceVector< double > &V) override
static::ExceptionBase & ExcNotImplemented()
virtual double operator*(const VectorSpaceVector< double > &V) const override
virtual void add(const double a) override
static::ExceptionBase & ExcZero()
#define AssertIsFinite(number)
virtual double linfty_norm() const override
static::ExceptionBase & ExcInternalError()
virtual double add_and_dot(const double a, const VectorSpaceVector< double > &V, const VectorSpaceVector< double > &W) override
::IndexSet source_stored_elements