Reference documentation for deal.II version GIT relicensing-218-g1f873f3929 2024-03-28 15:00: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 {
61 class Vector : public ReadVector<double>, public Subscriptor
62 {
63 public:
64 using value_type = double;
66 using real_type = double;
67
71 Vector();
72
77 Vector(const Vector &V);
78
85 explicit Vector(const IndexSet &parallel_partitioner,
86 const MPI_Comm communicator);
87
94 void
95 reinit(const IndexSet &parallel_partitioner,
96 const MPI_Comm communicator,
97 const bool omit_zeroing_entries = false);
98
103 void
104 reinit(const Vector &V, const bool omit_zeroing_entries = false);
105
109 virtual void
112 ArrayView<double> &elements) const override;
113
119 Vector &
120 operator=(const Vector &V);
121
126 Vector &
127 operator=(const double s);
128
137 void
140 VectorOperation::values operation,
141 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
142 &communication_pattern = {});
143
148 void
149 import(
151 const VectorOperation::values operation,
152 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
153 &communication_pattern = {})
154 {
155 import_elements(V, operation, communication_pattern);
156 }
157
161 Vector &
162 operator*=(const double factor);
163
167 Vector &
168 operator/=(const double factor);
169
173 Vector &
174 operator+=(const Vector &V);
175
179 Vector &
180 operator-=(const Vector &V);
181
186 double
187 operator*(const Vector &V) const;
188
192 void
193 add(const double a);
194
199 void
200 add(const double a, const Vector &V);
201
206 void
207 add(const double a, const Vector &V, const double b, const Vector &W);
208
213 void
214 sadd(const double s, const double a, const Vector &V);
215
222 void
223 scale(const Vector &scaling_factors);
224
228 void
229 equ(const double a, const Vector &V);
230
234 bool
235 all_zero() const;
236
240 double
241 mean_value() const;
242
247 double
248 l1_norm() const;
249
254 double
255 l2_norm() const;
256
261 double
262 linfty_norm() const;
263
286 double
287 add_and_dot(const double a, const Vector &V, const Vector &W);
292 bool
293 has_ghost_elements() const;
294
299 virtual size_type
300 size() const override;
301
307 locally_owned_size() const;
308
313 get_mpi_communicator() const;
314
328
339 void
340 compress(const VectorOperation::values operation);
341
346 const Epetra_FEVector &
347 trilinos_vector() const;
348
353 Epetra_FEVector &
355
359 void
360 print(std::ostream &out,
361 const unsigned int precision = 3,
362 const bool scientific = true,
363 const bool across = true) const;
364
368 std::size_t
369 memory_consumption() const;
370
376
383
390 int,
391 << "An error with error number " << arg1
392 << " occurred while calling a Trilinos function");
393
394 private:
400 void
401 create_epetra_comm_pattern(const IndexSet &source_index_set,
402 const MPI_Comm mpi_comm);
403
407 std::unique_ptr<Epetra_FEVector> vector;
408
413
418 std::shared_ptr<const CommunicationPattern> epetra_comm_pattern;
419 };
420
421
422 inline bool
424 {
425 return false;
426 }
427 } // namespace EpetraWrappers
428} // namespace LinearAlgebra
429
430
434template <>
435struct is_serial_vector<LinearAlgebra::EpetraWrappers::Vector> : std::false_type
436{};
437
439
440#endif
441
442#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)
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)
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: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)
static ::ExceptionBase & ExcDifferentParallelPartitioning()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
unsigned int global_dof_index
Definition types.h:81