Reference documentation for deal.II version GIT relicensing-439-g5fda5c893d 2024-04-20 06:50:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
trilinos_epetra_vector.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2016 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_trilinos_epetra_vector_h
16#define dealii_trilinos_epetra_vector_h
17
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_TRILINOS
22
26
31
32# include <Epetra_FEVector.h>
33
34# include <memory>
35
37
38namespace LinearAlgebra
39{
40 // Forward declaration
41 template <typename Number>
42 class ReadWriteVector;
43
51 namespace EpetraWrappers
52 {
58 {
59 public:
60 using value_type = double;
62 };
63
74 namespace internal
75 {
85 class VectorReference
86 {
87 private:
88 using value_type = VectorTraits::value_type;
89 using size_type = VectorTraits::size_type;
90
95 VectorReference(Vector &vector, const size_type index);
96
97 public:
101 VectorReference(const VectorReference &) = default;
102
114 const VectorReference &
115 operator=(const VectorReference &r) const;
116
120 VectorReference &
121 operator=(const VectorReference &r);
122
126 const VectorReference &
127 operator=(const value_type &s) const;
128
132 const VectorReference &
133 operator+=(const value_type &s) const;
134
138 const VectorReference &
139 operator-=(const value_type &s) const;
140
144 const VectorReference &
145 operator*=(const value_type &s) const;
146
150 const VectorReference &
151 operator/=(const value_type &s) const;
152
157 operator value_type() const;
158
162 DeclException1(ExcTrilinosError,
163 int,
164 << "An error with error number " << arg1
165 << " occurred while calling a Trilinos function");
166
167 /*
168 * Access to a an element that is not (locally-)owned.
169 *
170 * @ingroup Exceptions
171 */
173 ExcAccessToNonLocalElement,
174 size_type,
175 size_type,
176 size_type,
177 size_type,
178 << "You are trying to access element " << arg1
179 << " of a distributed vector, but this element is not stored "
180 << "on the current processor. Note: There are " << arg2
181 << " elements stored "
182 << "on the current processor from within the range [" << arg3 << ','
183 << arg4 << "] but Trilinos vectors need not store contiguous "
184 << "ranges on each processor, and not every element in "
185 << "this range may in fact be stored locally."
186 << "\n\n"
187 << "A common source for this kind of problem is that you "
188 << "are passing a 'fully distributed' vector into a function "
189 << "that needs read access to vector elements that correspond "
190 << "to degrees of freedom on ghost cells (or at least to "
191 << "'locally active' degrees of freedom that are not also "
192 << "'locally owned'). You need to pass a vector that has these "
193 << "elements as ghost entries.");
194
195 private:
199 Vector &vector;
200
204 const size_type index;
205
206 // Make the vector class a friend, so that it can create objects of the
207 // present type.
208 friend class ::LinearAlgebra::EpetraWrappers::Vector;
209 }; // class VectorReference
210
211 } // namespace internal
225 class Vector : public ReadVector<VectorTraits::value_type>,
226 public Subscriptor
227 {
228 public:
234
243 Vector();
244
249 Vector(const Vector &V);
250
257 explicit Vector(const IndexSet &parallel_partitioner,
258 const MPI_Comm communicator);
259
266 void
268 const MPI_Comm communicator,
269 const bool omit_zeroing_entries = false);
270
275 void
276 reinit(const Vector &V, const bool omit_zeroing_entries = false);
277
281 virtual void
284 ArrayView<double> &elements) const override;
285
315 Vector &
316 operator=(const Vector &V);
317
322 Vector &
323 operator=(const double s);
324
333 void
336 VectorOperation::values operation,
337 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
339
344 void
345 import(
347 const VectorOperation::values operation,
348 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
350 {
352 }
353
370 operator()(const size_type index);
371
380 operator()(const size_type index) const;
381
388 operator[](const size_type index);
389
396 operator[](const size_type index) const;
397
408 Vector &
409 operator*=(const double factor);
410
414 Vector &
415 operator/=(const double factor);
416
420 Vector &
421 operator+=(const Vector &V);
422
426 Vector &
427 operator-=(const Vector &V);
428
433 double
434 operator*(const Vector &V) const;
435
439 void
440 add(const double a);
441
446 void
447 add(const double a, const Vector &V);
448
453 void
454 add(const double a, const Vector &V, const double b, const Vector &W);
455
460 void
461 sadd(const double s, const double a, const Vector &V);
462
469 void
471
475 void
476 equ(const double a, const Vector &V);
477
481 bool
482 all_zero() const;
483
495 double
496 mean_value() const;
497
502 double
503 l1_norm() const;
504
509 double
510 l2_norm() const;
511
516 double
517 linfty_norm() const;
518
541 double
542 add_and_dot(const double a, const Vector &V, const Vector &W);
543
555 bool
556 has_ghost_elements() const;
557
562 virtual size_type
563 size() const override;
564
570 locally_owned_size() const;
571
576 get_mpi_communicator() const;
577
591
610 void
611 compress(const VectorOperation::values operation);
612
617 const Epetra_FEVector &
618 trilinos_vector() const;
619
624 Epetra_FEVector &
626
630 void
631 print(std::ostream &out,
632 const unsigned int precision = 3,
633 const bool scientific = true,
634 const bool across = true) const;
635
639 std::size_t
640 memory_consumption() const;
641
649
656
663 int,
664 << "An error with error number " << arg1
665 << " occurred while calling a Trilinos function");
666
667 private:
673 void
675 const MPI_Comm mpi_comm);
676
680 std::unique_ptr<Epetra_FEVector> vector;
681
686
691 std::shared_ptr<const CommunicationPattern> epetra_comm_pattern;
692
693 // Make the reference class a friend.
695 };
696
697# ifndef DOXYGEN
698
699 // VectorReference
700 namespace internal
701 {
702 inline VectorReference::VectorReference(Vector &vector,
703 const size_type index)
704 : vector(vector)
705 , index(index)
706 {}
707
708
709
710 inline const VectorReference &
711 VectorReference::operator=(const VectorReference &r) const
712 {
713 // as explained in the class
714 // documentation, this is not the copy
715 // operator. so simply pass on to the
716 // "correct" assignment operator
717 *this = static_cast<value_type>(r);
718
719 return *this;
720 }
721
722
723
724 inline VectorReference &
725 VectorReference::operator=(const VectorReference &r)
726 {
727 // as above
728 *this = static_cast<value_type>(r);
729
730 return *this;
731 }
732
733
734
735 inline const VectorReference &
736 VectorReference::operator=(const value_type &value) const
737 {
738 (*vector.vector)[0][index] = value;
739
740 return *this;
741 }
742
743
744
745 inline const VectorReference &
746 VectorReference::operator+=(const value_type &value) const
747 {
748 value_type new_value = static_cast<value_type>(*this) + value;
749 (*vector.vector)[0][index] = new_value;
750
751 return *this;
752 }
753
754
755
756 inline const VectorReference &
757 VectorReference::operator-=(const value_type &value) const
758 {
759 value_type new_value = static_cast<value_type>(*this) - value;
760 (*vector.vector)[0][index] = new_value;
761
762 return *this;
763 }
764
765
766
767 inline const VectorReference &
768 VectorReference::operator*=(const value_type &value) const
769 {
770 value_type new_value = static_cast<value_type>(*this) * value;
771 (*vector.vector)[0][index] = new_value;
772
773 return *this;
774 }
775
776
777
778 inline const VectorReference &
779 VectorReference::operator/=(const value_type &value) const
780 {
781 value_type new_value = static_cast<value_type>(*this) / value;
782 (*vector.vector)[0][index] = new_value;
783
784 return *this;
785 }
786 } // namespace internal
787
788# endif /* DOXYGEN */
789
790
791 inline internal::VectorReference
793 {
794 return internal::VectorReference(*this, index);
795 }
796
797 inline internal::VectorReference
799 {
800 return operator()(index);
801 }
802
803 inline Vector::value_type
804 Vector::operator[](const size_type index) const
805 {
806 return operator()(index);
807 }
808
809
810 inline bool
812 {
813 return false;
814 }
815 } // namespace EpetraWrappers
816} // namespace LinearAlgebra
817
818
822template <>
823struct is_serial_vector<LinearAlgebra::EpetraWrappers::Vector> : std::false_type
824{};
825
827
828#endif
829
830#endif
void equ(const double a, const Vector &V)
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)
value_type operator()(const size_type index) const
virtual size_type size() const override
reference operator[](const size_type index)
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)
std::shared_ptr< const CommunicationPattern > epetra_comm_pattern
void reinit(const IndexSet &parallel_partitioner, const MPI_Comm communicator, const bool omit_zeroing_entries=false)
reference operator()(const size_type index)
const internal::VectorReference const_reference
double add_and_dot(const double a, const Vector &V, const Vector &W)
bool has_ghost_elements() const
Number operator[](const size_type i) const
Number operator()(const size_type i) const
#define DEAL_II_DEPRECATED
Definition config.h:207
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcVectorTypeNotCompatible()
#define DeclException0(Exception0)
Definition exceptions.h:471
static ::ExceptionBase & ExcTrilinosError(int arg1)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition exceptions.h:585
static ::ExceptionBase & ExcDifferentParallelPartitioning()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
unsigned int global_dof_index
Definition types.h:81