16 #ifndef dealii__trilinos_vector_base_h 17 #define dealii__trilinos_vector_base_h 20 #include <deal.II/base/config.h> 22 #ifdef DEAL_II_WITH_TRILINOS 24 #include <deal.II/base/utilities.h> 25 # include <deal.II/base/std_cxx11/shared_ptr.h> 26 # include <deal.II/base/subscriptor.h> 27 # include <deal.II/lac/exceptions.h> 28 # include <deal.II/lac/vector.h> 29 # include <deal.II/base/mpi.h> 35 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
36 # include "Epetra_ConfigDefs.h" 37 # ifdef DEAL_II_WITH_MPI // only if MPI is installed 39 # include "Epetra_MpiComm.h" 41 # include "Epetra_SerialComm.h" 43 # include "Epetra_FEVector.h" 44 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
46 DEAL_II_NAMESPACE_OPEN
49 template <
typename number>
class Vector;
82 typedef ::types::global_dof_index
size_type;
100 VectorReference (VectorBase &vector,
101 const size_type index);
116 const VectorReference &
117 operator = (
const VectorReference &r)
const;
123 operator = (
const VectorReference &r);
128 const VectorReference &
129 operator = (
const TrilinosScalar &s)
const;
134 const VectorReference &
135 operator += (
const TrilinosScalar &s)
const;
140 const VectorReference &
141 operator -= (
const TrilinosScalar &s)
const;
146 const VectorReference &
147 operator *= (
const TrilinosScalar &s)
const;
152 const VectorReference &
153 operator /= (
const TrilinosScalar &s)
const;
159 operator TrilinosScalar ()
const;
166 <<
"An error with error number " << arg1
167 <<
" occurred while calling a Trilinos function");
178 const size_type index;
184 friend class ::TrilinosWrappers::VectorBase;
231 typedef TrilinosScalar real_type;
232 typedef ::types::global_dof_index size_type;
233 typedef value_type *iterator;
234 typedef const value_type *const_iterator;
235 typedef internal::VectorReference reference;
236 typedef const internal::VectorReference const_reference;
273 const bool omit_zeroing_entries =
false);
299 bool is_compressed () const DEAL_II_DEPRECATED;
314 operator = (const TrilinosScalar s);
332 template <typename Number>
334 operator = (const ::
Vector<Number> &v);
353 size_type size () const;
366 size_type local_size () const;
389 std::pair<size_type, size_type> local_range () const;
398 bool in_local_range (const size_type index) const;
413 IndexSet locally_owned_elements () const;
422 bool has_ghost_elements() const;
428 TrilinosScalar operator * (const
VectorBase &vec) const;
433 real_type norm_sqr () const;
438 TrilinosScalar mean_value () const;
445 TrilinosScalar minimal_value () const DEAL_II_DEPRECATED;
450 TrilinosScalar min () const;
455 TrilinosScalar max () const;
460 real_type l1_norm () const;
466 real_type l2_norm () const;
472 real_type lp_norm (const TrilinosScalar p) const;
477 real_type linfty_norm () const;
494 TrilinosScalar add_and_dot (const TrilinosScalar a,
503 bool all_zero () const;
510 bool is_non_negative () const;
529 operator () (const size_type index);
541 operator () (const size_type index) const;
549 operator [] (const size_type index);
557 operator [] (const size_type index) const;
565 void extract_subvector_to (const
std::vector<size_type> &indices,
566 std::vector<TrilinosScalar> &values) const;
572 template <typename ForwardIterator, typename OutputIterator>
573 void extract_subvector_to (ForwardIterator indices_begin,
574 const ForwardIterator indices_end,
575 OutputIterator values_begin) const;
587 TrilinosScalar el (const size_type index) const DEAL_II_DEPRECATED;
605 const_iterator begin () const;
617 const_iterator end () const;
633 void set (const
std::vector<size_type> &indices,
634 const
std::vector<TrilinosScalar> &values);
640 void set (const
std::vector<size_type> &indices,
641 const ::
Vector<TrilinosScalar> &values);
648 void set (const size_type n_elements,
649 const size_type *indices,
650 const TrilinosScalar *values);
656 void add (const
std::vector<size_type> &indices,
657 const
std::vector<TrilinosScalar> &values);
663 void add (const
std::vector<size_type> &indices,
664 const ::
Vector<TrilinosScalar> &values);
671 void add (const size_type n_elements,
672 const size_type *indices,
673 const TrilinosScalar *values);
678 VectorBase &operator *= (const TrilinosScalar factor);
683 VectorBase &operator /= (const TrilinosScalar factor);
699 void add (const TrilinosScalar s);
714 const
bool allow_different_maps = false);
719 void add (const TrilinosScalar a,
725 void add (const TrilinosScalar a,
727 const TrilinosScalar b,
734 void sadd (const TrilinosScalar s,
740 void sadd (const TrilinosScalar s,
741 const TrilinosScalar a,
749 void sadd (const TrilinosScalar s,
750 const TrilinosScalar a,
752 const TrilinosScalar b,
761 void sadd (const TrilinosScalar s,
762 const TrilinosScalar a,
764 const TrilinosScalar b,
766 const TrilinosScalar c,
774 void scale (const
VectorBase &scaling_factors);
779 void equ (const TrilinosScalar a,
787 void equ (const TrilinosScalar a,
789 const TrilinosScalar b,
816 const Epetra_MultiVector &trilinos_vector () const;
822 Epetra_FEVector &trilinos_vector ();
828 const Epetra_Map &vector_partitioner () const;
836 void print (const
char *format = 0) const DEAL_II_DEPRECATED;
845 void print (
std::ostream &out,
846 const
unsigned int precision = 3,
847 const
bool scientific = true,
848 const
bool across = true) const;
868 std::
size_t memory_consumption () const;
874 const MPI_Comm &get_mpi_communicator () const;
887 << "An error with error number " << arg1
888 << " occurred while calling a Trilinos function");
894 size_type, size_type, size_type, size_type,
895 << "You tried to access element " << arg1
896 << " of a distributed vector, but this element is not stored "
897 << "on the current processor. Note: There are "
898 << arg2 << " elements stored "
899 << "on the current processor from within the range "
900 << arg3 << " through " << arg4
901 << " but Trilinos vectors need not store contiguous "
902 << "ranges on each processor, and not every element in "
903 << "this range may in fact be stored locally.");
918 Epetra_CombineMode last_action;
944 std_cxx11::shared_ptr<Epetra_MultiVector> nonlocal_vector;
956 friend class MPI::Vector;
984 VectorReference::VectorReference (
VectorBase &vector,
985 const size_type index)
993 const VectorReference &
994 VectorReference::operator = (
const VectorReference &r)
const 1000 *
this =
static_cast<TrilinosScalar
> (r);
1009 VectorReference::operator = (
const VectorReference &r)
1012 *
this =
static_cast<TrilinosScalar
> (r);
1019 const VectorReference &
1020 VectorReference::operator = (
const TrilinosScalar &value)
const 1022 vector.
set (1, &index, &value);
1029 const VectorReference &
1030 VectorReference::operator += (
const TrilinosScalar &value)
const 1032 vector.
add (1, &index, &value);
1039 const VectorReference &
1040 VectorReference::operator -= (
const TrilinosScalar &value)
const 1042 TrilinosScalar new_value = -value;
1043 vector.
add (1, &index, &new_value);
1050 const VectorReference &
1051 VectorReference::operator *= (
const TrilinosScalar &value)
const 1053 TrilinosScalar new_value =
static_cast<TrilinosScalar
>(*this) * value;
1054 vector.
set (1, &index, &new_value);
1061 const VectorReference &
1062 VectorReference::operator /= (
const TrilinosScalar &value)
const 1064 TrilinosScalar new_value =
static_cast<TrilinosScalar
>(*this) / value;
1065 vector.
set (1, &index, &new_value);
1085 std::pair<size_type, size_type> range = local_range();
1087 return ((index >= range.first) && (index < range.second));
1096 Assert(owned_elements.size()==size(),
1097 ExcMessage(
"The locally owned elements have not been properly initialized!" 1098 " This happens for example if this object has been initialized" 1099 " with exactly one overlapping IndexSet."));
1100 return owned_elements;
1115 internal::VectorReference
1118 return internal::VectorReference (*
this, index);
1124 internal::VectorReference
1127 return operator() (index);
1135 return operator() (index);
1142 std::vector<TrilinosScalar> &values)
const 1144 for (size_type i = 0; i < indices.size(); ++i)
1145 values[i] =
operator()(indices[i]);
1150 template <
typename ForwardIterator,
typename OutputIterator>
1153 const ForwardIterator indices_end,
1154 OutputIterator values_begin)
const 1156 while (indices_begin != indices_end)
1158 *values_begin = operator()(*indices_begin);
1167 VectorBase::iterator
1170 return (*vector)[0];
1176 VectorBase::iterator
1179 return (*vector)[0]+local_size();
1185 VectorBase::const_iterator
1188 return (*vector)[0];
1194 VectorBase::const_iterator
1197 return (*vector)[0]+local_size();
1205 const bool omit_zeroing_entries)
1207 Assert (vector.get() != 0,
1208 ExcMessage(
"Vector has not been constructed properly."));
1210 if (omit_zeroing_entries ==
false ||
1212 vector.reset (
new Epetra_FEVector(*v.
vector));
1215 nonlocal_vector.reset(
new Epetra_MultiVector(v.
nonlocal_vector->Map(), 1));
1228 int ierr = vector->PutScalar(s);
1231 if (nonlocal_vector.get() != 0)
1233 ierr = nonlocal_vector->PutScalar(0.);
1245 const std::vector<TrilinosScalar> &values)
1251 Assert (indices.size() == values.size(),
1254 set (indices.size(), &indices[0], &values[0]);
1262 const ::Vector<TrilinosScalar> &values)
1268 Assert (indices.size() == values.size(),
1271 set (indices.size(), &indices[0], values.begin());
1279 const size_type *indices,
1280 const TrilinosScalar *values)
1286 if (last_action == Add)
1288 const int ierr = vector->GlobalAssemble(Add);
1292 if (last_action != Insert)
1293 last_action = Insert;
1295 for (size_type i=0; i<n_elements; ++i)
1297 const size_type row = indices[i];
1298 const TrilinosWrappers::types::int_type local_row = vector->Map().LID(static_cast<TrilinosWrappers::types::int_type>(row));
1299 if (local_row != -1)
1300 (*vector)[0][local_row] = values[i];
1303 const int ierr = vector->ReplaceGlobalValues (1,
1304 (
const TrilinosWrappers::types::int_type *)(&row),
1322 const std::vector<TrilinosScalar> &values)
1327 Assert (indices.size() == values.size(),
1330 add (indices.size(), &indices[0], &values[0]);
1338 const ::Vector<TrilinosScalar> &values)
1343 Assert (indices.size() == values.size(),
1346 add (indices.size(), &indices[0], values.begin());
1354 const size_type *indices,
1355 const TrilinosScalar *values)
1361 if (last_action != Add)
1363 if (last_action == Insert)
1365 const int ierr = vector->GlobalAssemble(Insert);
1371 for (size_type i=0; i<n_elements; ++i)
1373 const size_type row = indices[i];
1374 const TrilinosWrappers::types::int_type local_row = vector->Map().LID(static_cast<TrilinosWrappers::types::int_type>(row));
1375 if (local_row != -1)
1376 (*vector)[0][local_row] += values[i];
1377 else if (nonlocal_vector.get() == 0)
1379 const int ierr = vector->SumIntoGlobalValues (1,
1380 (
const TrilinosWrappers::types::int_type *)(&row),
1389 const TrilinosWrappers::types::int_type my_row = nonlocal_vector->Map().LID(static_cast<TrilinosWrappers::types::int_type>(row));
1391 ExcMessage(
"Attempted to write into off-processor vector entry " 1392 "that has not be specified as being writable upon " 1394 (*nonlocal_vector)[0][my_row] += values[i];
1403 VectorBase::size_type
1406 #ifndef DEAL_II_WITH_64BIT_INDICES 1407 return (size_type) (vector->Map().MaxAllGID() + 1 - vector->Map().MinAllGID());
1409 return (size_type) (vector->Map().MaxAllGID64() + 1 - vector->Map().MinAllGID64());
1416 VectorBase::size_type
1419 return (size_type) vector->Map().NumMyElements();
1425 std::pair<VectorBase::size_type, VectorBase::size_type>
1428 #ifndef DEAL_II_WITH_64BIT_INDICES 1429 const TrilinosWrappers::types::int_type begin = vector->Map().MinMyGID();
1430 const TrilinosWrappers::types::int_type end = vector->Map().MaxMyGID()+1;
1432 const TrilinosWrappers::types::int_type begin = vector->Map().MinMyGID64();
1433 const TrilinosWrappers::types::int_type end = vector->Map().MaxMyGID64()+1;
1436 Assert (end-begin == vector->Map().NumMyElements(),
1437 ExcMessage (
"This function only makes sense if the elements that this " 1438 "vector stores on the current processor form a contiguous range. " 1439 "This does not appear to be the case for the current vector."));
1441 return std::make_pair (begin, end);
1451 ExcDifferentParallelPartitioning());
1454 TrilinosScalar result;
1456 const int ierr = vector->Dot(*(vec.
vector), &result);
1465 VectorBase::real_type
1468 const TrilinosScalar d = l2_norm();
1480 TrilinosScalar mean;
1481 const int ierr = vector->MeanValue (&mean);
1502 TrilinosScalar min_value;
1503 const int ierr = vector->MinValue (&min_value);
1515 TrilinosScalar max_value;
1516 const int ierr = vector->MaxValue (&max_value);
1525 VectorBase::real_type
1531 const int ierr = vector->Norm1 (&d);
1540 VectorBase::real_type
1546 const int ierr = vector->Norm2 (&d);
1555 VectorBase::real_type
1560 TrilinosScalar norm = 0;
1561 TrilinosScalar sum=0;
1562 const size_type n_local = local_size();
1566 for (size_type i=0; i<n_local; i++)
1567 sum += std::pow(std::fabs((*vector)[0][i]), p);
1569 norm = std::pow(sum, static_cast<TrilinosScalar>(1./p));
1577 VectorBase::real_type
1586 const int ierr = vector->NormInf (&d);
1617 const int ierr = vector->Scale(a);
1631 const TrilinosScalar factor = 1./a;
1635 const int ierr = vector->Scale(factor);
1650 ExcDifferentParallelPartitioning());
1652 const int ierr = vector->Update (1.0, *(v.
vector), 1.0);
1667 ExcDifferentParallelPartitioning());
1669 const int ierr = vector->Update (-1.0, *(v.
vector), 1.0);
1686 size_type n_local = local_size();
1687 for (size_type i=0; i<n_local; i++)
1688 (*vector)[0][i] += s;
1706 const int ierr = vector->Update(a, *(v.
vector), 1.);
1716 const TrilinosScalar b,
1730 const int ierr = vector->Update(a, *(v.
vector), b, *(w.
vector), 1.);
1754 Assert (this->vector->Map().SameAs(v.
vector->Map())==
true,
1755 ExcDifferentParallelPartitioning());
1756 const int ierr = vector->Update(1., *(v.
vector), s);
1771 const TrilinosScalar a,
1786 Assert (this->vector->Map().SameAs(v.
vector->Map())==
true,
1787 ExcDifferentParallelPartitioning());
1788 const int ierr = vector->Update(a, *(v.
vector), s);
1796 this->add(tmp,
true);
1805 const TrilinosScalar a,
1807 const TrilinosScalar b,
1826 Assert (this->vector->Map().SameAs(v.
vector->Map())==
true,
1827 ExcDifferentParallelPartitioning());
1828 Assert (this->vector->Map().SameAs(w.
vector->Map())==
true,
1829 ExcDifferentParallelPartitioning());
1830 const int ierr = vector->Update(a, *(v.
vector), b, *(w.
vector), s);
1835 this->sadd( s, a, v);
1836 this->sadd(1., b, w);
1845 const TrilinosScalar a,
1847 const TrilinosScalar b,
1849 const TrilinosScalar c,
1872 Assert (this->vector->Map().SameAs(v.
vector->Map())==
true,
1873 ExcDifferentParallelPartitioning());
1874 Assert (this->vector->Map().SameAs(w.
vector->Map())==
true,
1875 ExcDifferentParallelPartitioning());
1876 Assert (this->vector->Map().SameAs(x.
vector->Map())==
true,
1877 ExcDifferentParallelPartitioning());
1882 const int ierr = vector->Update(a, *(v.
vector), b, *(w.
vector), s);
1885 const int jerr = vector->Update(c, *(x.
vector), 1.);
1891 this->sadd( s, a, v);
1892 this->sadd(1., b, w);
1893 this->sadd(1., c, x);
1909 const int ierr = vector->Multiply (1.0, *(factors.
vector), *vector, 0.0);
1926 if (vector->Map().SameAs(v.
vector->Map())==
false)
1928 this->sadd(0., a, v);
1933 int ierr = vector->Update(a, *v.
vector, 0.0);
1953 const int ierr = vector->ReciprocalMultiply(1.0, *(w.
vector),
1962 const Epetra_MultiVector &
1965 return static_cast<const Epetra_MultiVector &
>(*vector);
1983 return static_cast<const Epetra_Map &
>(vector->Map());
1992 static MPI_Comm comm;
1994 #ifdef DEAL_II_WITH_MPI 1996 const Epetra_MpiComm *mpi_comm
1997 =
dynamic_cast<const Epetra_MpiComm *
>(&vector->Map().Comm());
1998 comm = mpi_comm->Comm();
2002 comm = MPI_COMM_SELF;
2017 DEAL_II_NAMESPACE_CLOSE
2019 #endif // DEAL_II_WITH_TRILINOS static::ExceptionBase & ExcTrilinosError(int arg1)
types::global_dof_index size_type
size_type local_size() const
bool is_compressed() const 1
const MPI_Comm & get_mpi_communicator() const
reference operator()(const size_type index)
real_type l2_norm() const
IndexSet locally_owned_elements() const
TrilinosScalar max() const
TrilinosScalar operator*(const VectorBase &vec) const
void sadd(const TrilinosScalar s, const VectorBase &V)
void swap(Vector &u, Vector &v)
std_cxx11::shared_ptr< Epetra_MultiVector > nonlocal_vector
void scale(const VectorBase &scaling_factors)
reference operator[](const size_type index)
void reinit(const VectorBase &v, const bool omit_zeroing_entries=false)
real_type l1_norm() const
#define AssertThrow(cond, exc)
void set(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
VectorBase & operator+=(const VectorBase &V)
void add(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
bool has_ghost_elements() const
static::ExceptionBase & ExcMessage(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
VectorBase & operator/=(const TrilinosScalar factor)
std_cxx11::shared_ptr< Epetra_FEVector > vector
TrilinosScalar minimal_value() const 1
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DeclException0(Exception0)
TrilinosScalar value_type
const Epetra_Map & vector_partitioner() const
VectorBase & operator-=(const VectorBase &V)
static::ExceptionBase & ExcGhostsPresent()
TrilinosScalar mean_value() const
void ratio(const VectorBase &a, const VectorBase &b) 1
void extract_subvector_to(const std::vector< size_type > &indices, std::vector< TrilinosScalar > &values) const
TrilinosScalar min() const
VectorBase & operator=(const TrilinosScalar s)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
real_type linfty_norm() const
const Epetra_MultiVector & trilinos_vector() const
bool in_local_range(const size_type index) const
TrilinosScalar add_and_dot(const TrilinosScalar a, const VectorBase &V, const VectorBase &W)
real_type lp_norm(const TrilinosScalar p) const
VectorBase & operator*=(const TrilinosScalar factor)
#define AssertIsFinite(number)
std::pair< size_type, size_type > local_range() const
real_type norm_sqr() const
void equ(const TrilinosScalar a, const VectorBase &V)