18 #ifndef LIBMESH_TRILINOS_EPETRA_VECTOR_H 19 #define LIBMESH_TRILINOS_EPETRA_VECTOR_H 22 #include "libmesh/libmesh_common.h" 24 #ifdef LIBMESH_TRILINOS_HAVE_EPETRA 27 #include "libmesh/numeric_vector.h" 28 #include "libmesh/parallel.h" 31 #include "libmesh/ignore_warnings.h" 32 #include <Epetra_CombineMode.h> 33 #include <Epetra_Map.h> 34 #include <Epetra_MultiVector.h> 35 #include <Epetra_Vector.h> 36 #include <Epetra_MpiComm.h> 37 #include "libmesh/restore_warnings.h" 47 class Epetra_IntSerialDenseVector;
48 class Epetra_SerialDenseVector;
54 template <
typename T>
class SparseMatrix;
101 const std::vector<numeric_index_type> & ghost,
125 virtual void close ()
override;
130 virtual void clear () noexcept
override;
132 virtual void zero ()
override;
134 virtual std::unique_ptr<NumericVector<T>>
zero_clone ()
const override;
136 virtual std::unique_ptr<NumericVector<T>>
clone ()
const override;
140 const bool fast=
false,
144 const bool fast=
false,
149 const std::vector<numeric_index_type> & ghost,
150 const bool fast =
false,
154 const bool fast =
false)
override;
162 virtual Real min ()
const override;
164 virtual Real max ()
const override;
166 virtual T
sum ()
const override;
200 virtual void add (
const T s)
override;
213 const std::vector<numeric_index_type> & dof_indices)
override;
227 virtual void insert (
const T * v,
228 const std::vector<numeric_index_type> & dof_indices)
override;
230 virtual void scale (
const T factor)
override;
232 virtual void abs()
override;
236 virtual void localize (std::vector<T> & v_local)
const override;
241 const std::vector<numeric_index_type> & send_list)
const override;
243 virtual void localize (std::vector<T> & v_local,
244 const std::vector<numeric_index_type> & indices)
const override;
248 const std::vector<numeric_index_type> & send_list)
override;
257 const std::vector<numeric_index_type> & rows)
const override;
283 std::unique_ptr<Epetra_Map>
_map;
302 const double * values);
315 const Epetra_SerialDenseVector & values);
323 const double * values);
336 const Epetra_SerialDenseVector & values);
340 const int * numValuesPerID,
341 const double * values);
345 const int * numValuesPerID,
346 const double * values);
370 const double * values,
375 const int * numValuesPerID,
376 const double * values,
385 const double * values,
415 template <
typename T>
420 _destroy_vec_on_exit(true),
424 nonlocalIDs_(nullptr),
425 nonlocalElementSize_(nullptr),
427 allocatedNonlocalLength_(0),
428 nonlocalCoefs_(nullptr),
430 ignoreNonLocalEntries_(false)
437 template <
typename T>
443 _destroy_vec_on_exit(true),
447 nonlocalIDs_(nullptr),
448 nonlocalElementSize_(nullptr),
450 allocatedNonlocalLength_(0),
451 nonlocalCoefs_(nullptr),
453 ignoreNonLocalEntries_(false)
461 template <
typename T>
468 _destroy_vec_on_exit(true),
472 nonlocalIDs_(nullptr),
473 nonlocalElementSize_(nullptr),
475 allocatedNonlocalLength_(0),
476 nonlocalCoefs_(nullptr),
478 ignoreNonLocalEntries_(false)
480 this->
init(n, n_local,
false,
type);
486 template <
typename T>
491 _destroy_vec_on_exit(false),
495 nonlocalIDs_(nullptr),
496 nonlocalElementSize_(nullptr),
498 allocatedNonlocalLength_(0),
499 nonlocalCoefs_(nullptr),
501 ignoreNonLocalEntries_(false)
510 _map = std::make_unique<Epetra_Map>
511 (
_vec->GlobalLength(),
514 Epetra_MpiComm (this->
comm().
get()));
527 template <
typename T>
532 const std::vector<numeric_index_type> & ghost,
535 _destroy_vec_on_exit(true),
539 nonlocalIDs_(nullptr),
540 nonlocalElementSize_(nullptr),
542 allocatedNonlocalLength_(0),
543 nonlocalCoefs_(nullptr),
545 ignoreNonLocalEntries_(false)
547 this->
init(n, n_local, ghost,
false,
type);
563 template <
typename T>
572 template <
typename T>
603 _map = std::make_unique<Epetra_Map>
604 (
static_cast<int>(n),
607 Epetra_MpiComm (this->comm().get()));
609 _vec =
new Epetra_Vector(*_map);
611 myFirstID_ = _vec->Map().MinMyGID();
612 myNumIDs_ = _vec->Map().NumMyElements();
617 _vec->ExtractView(&myCoefs_, &dummy);
620 this->_is_closed =
true;
628 template <
typename T>
632 const std::vector<numeric_index_type> & ,
637 this->
init(n, n_local, fast, type);
642 template <
typename T>
648 this->
init(n,n,fast,type);
653 template <
typename T>
660 unsigned char global_last_edit = last_edit;
661 this->comm().max(global_last_edit);
664 if (global_last_edit == 1)
665 this->GlobalAssemble(Insert);
666 else if (global_last_edit == 2)
667 this->GlobalAssemble(Add);
671 this->_is_closed =
true;
677 template <
typename T>
684 if (this->_destroy_vec_on_exit)
699 template <
typename T>
706 _vec->PutScalar(0.0);
711 template <
typename T>
716 cloned_vector->
init(*
this);
717 return std::unique_ptr<NumericVector<T>>(cloned_vector);
722 template <
typename T>
727 cloned_vector->
init(*
this,
true);
728 *cloned_vector = *
this;
729 return std::unique_ptr<NumericVector<T>>(cloned_vector);
734 template <
typename T>
740 return _vec->GlobalLength();
745 template <
typename T>
751 return _vec->MyLength();
754 template <
typename T>
760 return _vec->Map().MinMyGID();
765 template <
typename T>
771 return _vec->Map().MaxMyGID()+1;
775 template <
typename T>
781 (i < this->last_local_index())) );
783 return (*_vec)[i-this->first_local_index()];
788 template <
typename T>
796 _vec->MinValue(&
value);
803 template <
typename T>
811 _vec->MaxValue(&
value);
818 template <
typename T>
826 std::swap(_vec, v.
_vec);
843 template <
typename T>
848 return std::numeric_limits<int>::max();
863 #endif // #ifdef LIBMESH_TRILINOS_HAVE_EPETRA 864 #endif // LIBMESH_TRILINOS_EPETRA_VECTOR_H
virtual NumericVector< T > & operator/=(const NumericVector< T > &v) override
Computes the component-wise division of this vector's entries by another's, .
bool closed()
Checks that the library has been closed.
This class provides a nice interface to the Trilinos Epetra_Vector object.
virtual T operator()(const numeric_index_type i) const override
std::unique_ptr< Epetra_Map > _map
Holds the distributed Map.
virtual void pointwise_mult(const NumericVector< T > &vec1, const NumericVector< T > &vec2) override
Computes (summation not implied) i.e.
int allocatedNonlocalLength_
virtual void reciprocal() override
Computes the component-wise reciprocal, .
virtual numeric_index_type last_local_index() const override
virtual void add(const numeric_index_type i, const T value) override
Adds value to the vector entry specified by i.
virtual void close() override
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
virtual std::unique_ptr< NumericVector< T > > zero_clone() const override
virtual numeric_index_type size() const =0
Epetra_Vector * _vec
Actual Epetra vector datatype to hold vector entries.
virtual void localize_to_one(std::vector< T > &v_local, const processor_id_type proc_id=0) const override
Creates a local copy of the global vector in v_local only on processor proc_id.
Provides a uniform interface to vector storage schemes for different linear algebra libraries...
const Parallel::Communicator & comm() const
virtual void insert(const T *v, const std::vector< numeric_index_type > &dof_indices) override
Inserts the entries of v in *this at the locations specified by v.
int inputNonlocalValues(int GID, int numValues, const double *values, bool accumulate)
bool _is_initialized
true once init() has been called.
virtual void conjugate() override
Negates the imaginary component of each entry in the vector.
virtual Real max() const override
The libMesh namespace provides an interface to certain functionality in the library.
virtual void zero() override
Set all entries to zero.
int * nonlocalElementSize_
virtual T sum() const override
virtual void init(const numeric_index_type n, const numeric_index_type n_local, const bool fast=false, const ParallelType ptype=AUTOMATIC)=0
Change the dimension of the vector to n.
int GlobalAssemble(Epetra_CombineMode mode=Add)
Gather any overlapping/shared data into the non-overlapping partitioning defined by the Map that was ...
virtual numeric_index_type size() const override
virtual void add_vector_transpose(const NumericVector< T > &v, const SparseMatrix< T > &A) override
Computes , i.e.
uint8_t processor_id_type
int ReplaceGlobalValues(int numIDs, const int *GIDs, const double *values)
Copy values into the vector overwriting any values that already exist for the specified indices...
int inputValues(int numIDs, const int *GIDs, const double *values, bool accumulate)
bool _destroy_vec_on_exit
This boolean value should only be set to false for the constructor which takes a Epetra Vec object...
void setIgnoreNonLocalEntries(bool flag)
Set whether or not non-local data values should be ignored.
unsigned char last_edit
Keep track of whether the last write operation on this vector was nothing (0) or a sum (1) or an add ...
dof_id_type numeric_index_type
bool _is_initialized
Flag that tells if init() has been called.
EpetraVector(const Parallel::Communicator &comm, const ParallelType type=AUTOMATIC)
Dummy-Constructor.
virtual NumericVector< T > & operator*=(const NumericVector< T > &v) override
Computes the component-wise multiplication of this vector's entries by another's, ...
virtual NumericVector< T > & operator+=(const NumericVector< T > &v) override
Adds v to *this, .
virtual Real l1_norm() const override
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
virtual void abs() override
Sets for each entry in the vector.
virtual NumericVector< T > & operator-=(const NumericVector< T > &v) override
Subtracts v from *this, .
int inputNonlocalValue(int GID, double value, bool accumulate)
virtual void clear() noexcept override
clear() is called from the destructor, so it should not throw.
bool ignoreNonLocalEntries_
virtual Real min() const override
virtual Real linfty_norm() const override
virtual void init(const numeric_index_type N, const numeric_index_type n_local, const bool fast=false, const ParallelType type=AUTOMATIC) override
Change the dimension of the vector to n.
virtual std::size_t max_allowed_id() const override
ParallelType _type
Type of vector.
void destroyNonlocalData()
ParallelType type() const
virtual void localize(std::vector< T > &v_local) const override
Creates a copy of the global vector in the local vector v_local.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void swap(NumericVector< T > &v)
Swaps the contents of this with v.
virtual T dot(const NumericVector< T > &v) const override
int SumIntoGlobalValues(int numIDs, const int *GIDs, const double *values)
Accumulate values into the vector, adding them to any values that already exist for the specified ind...
virtual Real l2_norm() const override
virtual numeric_index_type local_size() const =0
virtual void create_subvector(NumericVector< T > &subvector, const std::vector< numeric_index_type > &rows) const override
Fills in subvector from this vector using the indices in rows.
virtual numeric_index_type first_local_index() const override
void FEoperatorequals(const EpetraVector &source)
virtual void swap(NumericVector< T > &v) override
Swaps the contents of this with v.
virtual std::unique_ptr< NumericVector< T > > clone() const override
bool initialized()
Checks that library initialization has been done.
int * numeric_trilinos_cast(const numeric_index_type *p)
virtual void scale(const T factor) override
Scale each element of the vector by the given factor.
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices) override
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i]...
bool _is_closed
Flag which tracks whether the vector's values are consistent on all processors after insertion or add...
virtual numeric_index_type local_size() const override
EpetraVector & operator=(const EpetraVector &)=delete
ParallelType
Defines an enum for parallel data structure types.