Reference documentation for deal.II version GIT 574c7c8486 2023-09-22 10:20: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\}}\)
trilinos_epetra_vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2015 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_trilinos_epetra_vector_h
17 #define dealii_trilinos_epetra_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_TRILINOS
23 
24 # include <deal.II/base/index_set.h>
25 # include <deal.II/base/mpi_stub.h>
27 
28 # include <deal.II/lac/read_vector.h>
32 
33 # include <Epetra_FEVector.h>
34 
35 # include <memory>
36 
38 
39 namespace LinearAlgebra
40 {
41  // Forward declaration
42  template <typename Number>
43  class ReadWriteVector;
44 
52  namespace EpetraWrappers
53  {
62  class Vector : public ReadVector<double>, public Subscriptor
63  {
64  public:
65  using value_type = double;
67  using real_type = double;
68 
72  Vector();
73 
78  Vector(const Vector &V);
79 
86  explicit Vector(const IndexSet &parallel_partitioner,
87  const MPI_Comm communicator);
88 
95  void
96  reinit(const IndexSet &parallel_partitioner,
97  const MPI_Comm communicator,
98  const bool omit_zeroing_entries = false);
99 
104  void
105  reinit(const Vector &V, const bool omit_zeroing_entries = false);
106 
110  virtual void
113  ArrayView<double> &elements) const override;
114 
120  Vector &
121  operator=(const Vector &V);
122 
127  Vector &
128  operator=(const double s);
129 
138  void
140  const ReadWriteVector<double> &V,
141  VectorOperation::values operation,
142  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
143  &communication_pattern = {});
144 
149  void
150  import(
151  const ReadWriteVector<double> &V,
152  const VectorOperation::values operation,
153  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
154  &communication_pattern = {})
155  {
156  import_elements(V, operation, communication_pattern);
157  }
158 
162  Vector &
163  operator*=(const double factor);
164 
168  Vector &
169  operator/=(const double factor);
170 
174  Vector &
175  operator+=(const Vector &V);
176 
180  Vector &
181  operator-=(const Vector &V);
182 
187  double
188  operator*(const Vector &V) const;
189 
193  void
194  add(const double a);
195 
200  void
201  add(const double a, const Vector &V);
202 
207  void
208  add(const double a, const Vector &V, const double b, const Vector &W);
209 
214  void
215  sadd(const double s, const double a, const Vector &V);
216 
223  void
224  scale(const Vector &scaling_factors);
225 
229  void
230  equ(const double a, const Vector &V);
231 
235  bool
236  all_zero() const;
237 
241  double
242  mean_value() const;
243 
248  double
249  l1_norm() const;
250 
255  double
256  l2_norm() const;
257 
262  double
263  linfty_norm() const;
264 
287  double
288  add_and_dot(const double a, const Vector &V, const Vector &W);
293  bool
294  has_ghost_elements() const;
295 
300  virtual size_type
301  size() const override;
302 
307  size_type
308  locally_owned_size() const;
309 
313  MPI_Comm
314  get_mpi_communicator() const;
315 
327  ::IndexSet
328  locally_owned_elements() const;
329 
340  void
341  compress(const VectorOperation::values operation);
342 
347  const Epetra_FEVector &
348  trilinos_vector() const;
349 
354  Epetra_FEVector &
355  trilinos_vector();
356 
360  void
361  print(std::ostream &out,
362  const unsigned int precision = 3,
363  const bool scientific = true,
364  const bool across = true) const;
365 
369  std::size_t
370  memory_consumption() const;
371 
377 
384 
391  int,
392  << "An error with error number " << arg1
393  << " occurred while calling a Trilinos function");
394 
395  private:
401  void
402  create_epetra_comm_pattern(const IndexSet &source_index_set,
403  const MPI_Comm mpi_comm);
404 
408  std::unique_ptr<Epetra_FEVector> vector;
409 
414 
419  std::shared_ptr<const CommunicationPattern> epetra_comm_pattern;
420  };
421 
422 
423  inline bool
425  {
426  return false;
427  }
428  } // namespace EpetraWrappers
429 } // namespace LinearAlgebra
430 
431 
435 template <>
437 {};
438 
440 
441 #endif
442 
443 #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)
Vector & operator/=(const double factor)
virtual size_type size() const override
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)
Vector & operator*=(const double factor)
std::shared_ptr< const CommunicationPattern > epetra_comm_pattern
void reinit(const IndexSet &parallel_partitioner, const MPI_Comm communicator, const bool omit_zeroing_entries=false)
double add_and_dot(const double a, const Vector &V, const Vector &W)
#define DEAL_II_DEPRECATED
Definition: config.h:178
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define DeclException0(Exception0)
Definition: exceptions.h:467
static ::ExceptionBase & ExcVectorTypeNotCompatible()
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcDifferentParallelPartitioning()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:512
static const char V
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
unsigned int global_dof_index
Definition: types.h:82