From e400eb2eaffe290cc1b608eaf15f414d5570b6b7 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 3 Jul 2023 14:41:56 -0600 Subject: [PATCH] Get rid of VectorSpaceVector in the LinearAlgebra::TpetraWrappers::Vector class. --- include/deal.II/lac/trilinos_tpetra_vector.h | 132 ++++++++------- .../lac/trilinos_tpetra_vector.templates.h | 157 ++++++------------ 2 files changed, 120 insertions(+), 169 deletions(-) diff --git a/include/deal.II/lac/trilinos_tpetra_vector.h b/include/deal.II/lac/trilinos_tpetra_vector.h index 3ad3cd01f2..6cd715bd1a 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.h @@ -95,8 +95,7 @@ namespace LinearAlgebra { /** * This class implements a wrapper to the Trilinos distributed vector - * class Tpetra::Vector. This class is derived from the - * LinearAlgebra::VectorSpaceVector class and requires Trilinos to be + * class Tpetra::Vector. This class requires Trilinos to be * compiled with MPI support. * * Tpetra uses Kokkos for thread-parallelism and chooses the execution and @@ -114,14 +113,12 @@ namespace LinearAlgebra * @ingroup Vectors */ template - class Vector : public VectorSpaceVector, - public ReadVector, - public Subscriptor + class Vector : public ReadVector, public Subscriptor { public: using value_type = Number; - - using size_type = typename VectorSpaceVector::size_type; + using real_type = typename numbers::NumberTraits::real_type; + using size_type = types::global_dof_index; /** * Constructor. Create a vector of dimension zero. @@ -158,9 +155,8 @@ namespace LinearAlgebra * Change the dimension to that of the vector @p V. The elements of @p V are not * copied. */ - virtual void - reinit(const VectorSpaceVector &V, - const bool omit_zeroing_entries = false) override; + void + reinit(const Vector &V, const bool omit_zeroing_entries = false); /** * Extract a range of elements all at once. @@ -182,8 +178,8 @@ namespace LinearAlgebra * Sets all elements of the vector to the scalar @p s. This operation is * only allowed if @p s is equal to zero. */ - virtual Vector & - operator=(const Number s) override; + Vector & + operator=(const Number s); /** * Imports all the elements present in the vector's IndexSet from the @@ -193,22 +189,22 @@ namespace LinearAlgebra * communication pattern is used multiple times. This can be used to * improve performance. */ - virtual void + void import_elements( const ReadWriteVector &V, VectorOperation::values operation, std::shared_ptr - communication_pattern = {}) override; + communication_pattern = {}); /** * @deprecated Use import_elements() instead. */ DEAL_II_DEPRECATED - virtual void + void import(const ReadWriteVector &V, VectorOperation::values operation, std::shared_ptr - communication_pattern = {}) override + communication_pattern = {}) { import_elements(V, operation, communication_pattern); } @@ -216,65 +212,63 @@ namespace LinearAlgebra /** * Multiply the entire vector by a fixed factor. */ - virtual Vector & - operator*=(const Number factor) override; + Vector & + operator*=(const Number factor); /** * Divide the entire vector by a fixed factor. */ - virtual Vector & - operator/=(const Number factor) override; + Vector & + operator/=(const Number factor); /** * Add the vector @p V to the present one. */ - virtual Vector & - operator+=(const VectorSpaceVector &V) override; + Vector & + operator+=(const Vector &V); /** * Subtract the vector @p V from the present one. */ - virtual Vector & - operator-=(const VectorSpaceVector &V) override; + Vector & + operator-=(const Vector &V); /** * Return the scalar product of two vectors. The vectors need to have the * same layout. */ - virtual Number - operator*(const VectorSpaceVector &V) const override; + Number + operator*(const Vector &V) const; /** * Add @p a to all components. Note that @p is a scalar not a vector. */ - virtual void - add(const Number a) override; + void + add(const Number a); /** * Simple addition of a multiple of a vector, i.e. *this += * a*V. The vectors need to have the same layout. */ - virtual void - add(const Number a, const VectorSpaceVector &V) override; + void + add(const Number a, const Vector &V); /** * Multiple addition of multiple of a vector, i.e. *this> += * a*V+b*W. The vectors need to have the same layout. */ - virtual void - add(const Number a, - const VectorSpaceVector &V, - const Number b, - const VectorSpaceVector &W) override; + void + add(const Number a, + const Vector &V, + const Number b, + const Vector &W); /** * Scaling and simple addition of a multiple of a vector, i.e. *this * = s*(*this)+a*V. */ - virtual void - sadd(const Number s, - const Number a, - const VectorSpaceVector &V) override; + void + sadd(const Number s, const Number a, const Vector &V); /** * Scale each element of this vector by the corresponding element in the @@ -282,47 +276,47 @@ namespace LinearAlgebra * (and immediate re-assignment) by a diagonal scaling matrix. The * vectors need to have the same layout. */ - virtual void - scale(const VectorSpaceVector &scaling_factors) override; + void + scale(const Vector &scaling_factors); /** * Assignment *this = a*V. */ - virtual void - equ(const Number a, const VectorSpaceVector &V) override; + void + equ(const Number a, const Vector &V); /** * Return whether the vector contains only elements with value zero. */ - virtual bool - all_zero() const override; + bool + all_zero() const; /** * Return the mean value of the element of this vector. */ - virtual Number - mean_value() const override; + Number + mean_value() const; /** * Return the l1 norm of the vector (i.e., the sum of the * absolute values of all entries among all processors). */ - virtual typename LinearAlgebra::VectorSpaceVector::real_type - l1_norm() const override; + real_type + l1_norm() const; /** * Return the l2 norm of the vector (i.e., the square root of * the sum of the square of all entries among all processors). */ - virtual typename LinearAlgebra::VectorSpaceVector::real_type - l2_norm() const override; + real_type + l2_norm() const; /** * Return the maximum norm of the vector (i.e., the maximum absolute value * among all entries and among all processors). */ - virtual typename LinearAlgebra::VectorSpaceVector::real_type - linfty_norm() const override; + real_type + linfty_norm() const; /** * Performs a combined operation of a vector addition and a subsequent @@ -346,10 +340,10 @@ namespace LinearAlgebra * implemented as * $\left=\sum_i v_i \bar{w_i}$. */ - virtual Number - add_and_dot(const Number a, - const VectorSpaceVector &V, - const VectorSpaceVector &W) override; + Number + add_and_dot(const Number a, + const Vector &V, + const Vector &W); /** * This function always returns false and is present only for backward * compatibility. @@ -388,9 +382,21 @@ namespace LinearAlgebra * vec.locally_owned_elements() == complete_index_set(vec.size()) * @endcode */ - virtual ::dealii::IndexSet - locally_owned_elements() const override; + ::dealii::IndexSet + locally_owned_elements() const; + /** + * Compress the underlying representation of the Trilinos object, i.e. + * flush the buffers of the vector object if it has any. This function is + * necessary after writing into a vector element-by-element and before + * anything else can be done on it. + * + * See + * @ref GlossCompress "Compressing distributed objects" + * for more information. + */ + void + compress(const VectorOperation::values operation); /** * Return a const reference to the underlying Trilinos * Tpetra::Vector class. @@ -408,17 +414,17 @@ namespace LinearAlgebra /** * Prints the vector to the output stream @p out. */ - virtual void + void print(std::ostream & out, const unsigned int precision = 3, const bool scientific = true, - const bool across = true) const override; + const bool across = true) const; /** * Return the memory consumption of this class in bytes. */ - virtual std::size_t - memory_consumption() const override; + std::size_t + memory_consumption() const; /** * The vectors have different partitioning, i.e. their IndexSet objects diff --git a/include/deal.II/lac/trilinos_tpetra_vector.templates.h b/include/deal.II/lac/trilinos_tpetra_vector.templates.h index e1f1dc06b5..fc8fb5c057 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h @@ -101,18 +101,11 @@ namespace LinearAlgebra template void - Vector::reinit(const VectorSpaceVector &V, - const bool omit_zeroing_entries) + Vector::reinit(const Vector &V, + const bool omit_zeroing_entries) { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); - - reinit(down_V.locally_owned_elements(), - down_V.get_mpi_communicator(), + reinit(V.locally_owned_elements(), + V.get_mpi_communicator(), omit_zeroing_entries); } @@ -296,34 +289,26 @@ namespace LinearAlgebra template Vector & - Vector::operator+=(const VectorSpaceVector &V) + Vector::operator+=(const Vector &V) { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); // If the maps are the same we can update right away. - if (vector->getMap()->isSameAs(*(down_V.trilinos_vector().getMap()))) + if (vector->getMap()->isSameAs(*(V.trilinos_vector().getMap()))) { - vector->update(1., down_V.trilinos_vector(), 1.); + vector->update(1., V.trilinos_vector(), 1.); } else { - Assert(this->size() == down_V.size(), - ExcDimensionMismatch(this->size(), down_V.size())); + Assert(this->size() == V.size(), + ExcDimensionMismatch(this->size(), V.size())); // TODO: Tpetra doesn't have a combine mode that also updates local // elements, maybe there is a better workaround. Tpetra::Vector dummy( vector->getMap(), false); Tpetra::Import data_exchange( - down_V.trilinos_vector().getMap(), dummy.getMap()); + V.trilinos_vector().getMap(), dummy.getMap()); - dummy.doImport(down_V.trilinos_vector(), - data_exchange, - Tpetra::INSERT); + dummy.doImport(V.trilinos_vector(), data_exchange, Tpetra::INSERT); vector->update(1.0, dummy, 1.0); } @@ -335,7 +320,7 @@ namespace LinearAlgebra template Vector & - Vector::operator-=(const VectorSpaceVector &V) + Vector::operator-=(const Vector &V) { this->add(-1., V); @@ -346,20 +331,14 @@ namespace LinearAlgebra template Number - Vector::operator*(const VectorSpaceVector &V) const + Vector::operator*(const Vector &V) const { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); - Assert(this->size() == down_V.size(), - ExcDimensionMismatch(this->size(), down_V.size())); - Assert(vector->getMap()->isSameAs(*down_V.trilinos_vector().getMap()), + Assert(this->size() == V.size(), + ExcDimensionMismatch(this->size(), V.size())); + Assert(vector->getMap()->isSameAs(*V.trilinos_vector().getMap()), ExcDifferentParallelPartitioning()); - return vector->dot(down_V.trilinos_vector()); + return vector->dot(V.trilinos_vector()); } @@ -397,68 +376,45 @@ namespace LinearAlgebra template void - Vector::add(const Number a, const VectorSpaceVector &V) + Vector::add(const Number a, const Vector &V) { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); AssertIsFinite(a); - Assert(vector->getMap()->isSameAs(*(down_V.trilinos_vector().getMap())), + Assert(vector->getMap()->isSameAs(*(V.trilinos_vector().getMap())), ExcDifferentParallelPartitioning()); - vector->update(a, down_V.trilinos_vector(), 1.); + vector->update(a, V.trilinos_vector(), 1.); } template void - Vector::add(const Number a, - const VectorSpaceVector &V, - const Number b, - const VectorSpaceVector &W) - { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - // Check that casting will work. - Assert(dynamic_cast *>(&W) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); - // Downcast W. If fails, throws an exception. - const Vector &down_W = dynamic_cast &>(W); - Assert(vector->getMap()->isSameAs(*(down_V.trilinos_vector().getMap())), + Vector::add(const Number a, + const Vector &V, + const Number b, + const Vector &W) + { + Assert(vector->getMap()->isSameAs(*(V.trilinos_vector().getMap())), ExcDifferentParallelPartitioning()); - Assert(vector->getMap()->isSameAs(*(down_W.trilinos_vector().getMap())), + Assert(vector->getMap()->isSameAs(*(W.trilinos_vector().getMap())), ExcDifferentParallelPartitioning()); AssertIsFinite(a); AssertIsFinite(b); - vector->update( - a, down_V.trilinos_vector(), b, down_W.trilinos_vector(), 1.); + vector->update(a, V.trilinos_vector(), b, W.trilinos_vector(), 1.); } template void - Vector::sadd(const Number s, - const Number a, - const VectorSpaceVector &V) + Vector::sadd(const Number s, + const Number a, + const Vector &V) { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - *this *= s; - // Downcast V. It fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); - Vector tmp(down_V); + + Vector tmp(V); tmp *= a; *this += tmp; } @@ -467,45 +423,28 @@ namespace LinearAlgebra template void - Vector::scale(const VectorSpaceVector &scaling_factors) + Vector::scale(const Vector &scaling_factors) { - // Check that casting will work. - Assert(dynamic_cast *>(&scaling_factors) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast scaling_factors. If fails, throws an exception. - const Vector &down_scaling_factors = - dynamic_cast &>(scaling_factors); Assert(vector->getMap()->isSameAs( - *(down_scaling_factors.trilinos_vector().getMap())), + *(scaling_factors.trilinos_vector().getMap())), ExcDifferentParallelPartitioning()); - vector->elementWiseMultiply(1., - *down_scaling_factors.vector, - *vector, - 0.); + vector->elementWiseMultiply(1., *scaling_factors.vector, *vector, 0.); } template void - Vector::equ(const Number a, const VectorSpaceVector &V) + Vector::equ(const Number a, const Vector &V) { - // Check that casting will work. - Assert(dynamic_cast *>(&V) != nullptr, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast &>(V); // If we don't have the same map, copy. - if (vector->getMap()->isSameAs(*down_V.trilinos_vector().getMap()) == - false) + if (vector->getMap()->isSameAs(*V.trilinos_vector().getMap()) == false) this->sadd(0., a, V); else { // Otherwise, just update - vector->update(a, down_V.trilinos_vector(), 0.); + vector->update(a, V.trilinos_vector(), 0.); } } @@ -554,7 +493,7 @@ namespace LinearAlgebra template - typename LinearAlgebra::VectorSpaceVector::real_type + typename Vector::real_type Vector::l1_norm() const { return vector->norm1(); @@ -563,7 +502,7 @@ namespace LinearAlgebra template - typename LinearAlgebra::VectorSpaceVector::real_type + typename Vector::real_type Vector::l2_norm() const { return vector->norm2(); @@ -572,7 +511,7 @@ namespace LinearAlgebra template - typename LinearAlgebra::VectorSpaceVector::real_type + typename Vector::real_type Vector::linfty_norm() const { return vector->normInf(); @@ -582,9 +521,9 @@ namespace LinearAlgebra template Number - Vector::add_and_dot(const Number a, - const VectorSpaceVector &V, - const VectorSpaceVector &W) + Vector::add_and_dot(const Number a, + const Vector &V, + const Vector &W) { this->add(a, V); @@ -653,6 +592,12 @@ namespace LinearAlgebra + template + void + Vector::compress(const VectorOperation::values /*operation*/) + {} + + template const Tpetra::Vector & Vector::trilinos_vector() const -- 2.39.5