This gives a common interface to all vectors for use by FEValues.
#include <deal.II/lac/block_indices.h>
#include <deal.II/lac/exceptions.h>
+#include <deal.II/lac/read_vector.h>
#include <deal.II/lac/vector.h>
#include <deal.II/lac/vector_operation.h>
* @ref GlossBlockLA "Block (linear algebra)"
*/
template <class VectorType>
-class BlockVectorBase : public Subscriptor
+class BlockVectorBase : public Subscriptor,
+ public ReadVector<typename VectorType::value_type>
{
public:
/**
* Return dimension of the vector. This is the sum of the dimensions of all
* components.
*/
- std::size_t
- size() const;
+ size_type
+ size() const override;
/**
* Return local dimension of the vector. This is the sum of the local
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<OtherNumber> & values) const;
+ virtual void
+ extract_subvector_to(const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<value_type> &entries) const override;
+
/**
* Instead of getting individual elements of a vector via operator(),
* this function allows getting a whole set of elements at once. In
template <class VectorType>
-inline std::size_t
+inline typename BlockVectorBase<VectorType>::size_type
BlockVectorBase<VectorType>::size() const
{
return block_indices.total_size();
+template <typename VectorType>
+inline void
+BlockVectorBase<VectorType>::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<value_type> & entries) const
+{
+ AssertDimension(indices.size(), entries.size());
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ {
+ AssertIndexRange(indices[i], size());
+ entries[i] = (*this)[indices[i]];
+ }
+}
+
+
+
template <typename VectorType>
template <typename ForwardIterator, typename OutputIterator>
inline void
#include <deal.II/base/partitioner.h>
#include <deal.II/base/subscriptor.h>
+#include <deal.II/lac/read_vector.h>
#include <deal.II/lac/vector_operation.h>
#include <deal.II/lac/vector_space_vector.h>
#include <deal.II/lac/vector_type_traits.h>
*/
template <typename Number, typename MemorySpace = MemorySpace::Host>
class Vector : public ::dealii::LinearAlgebra::VectorSpaceVector<Number>,
+ public ::dealii::ReadVector<Number>,
public Subscriptor
{
public:
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<OtherNumber> & values) const;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> &elements) const override;
+
/**
* Instead of getting individual elements of a vector via operator(),
* this function allows getting a whole set of elements at once. In
+ template <typename Number, typename MemorySpaceType>
+ void
+ Vector<Number, MemorySpaceType>::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> & elements) const
+ {
+ AssertDimension(indices.size(), elements.size());
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ {
+ AssertIndexRange(indices[i], size());
+ elements[i] = (*this)[indices[i]];
+ }
+ }
+
+
+
template <typename Number, typename MemorySpaceType>
bool
Vector<Number, MemorySpaceType>::all_zero() const
*
* @ingroup PETScWrappers
*/
- class VectorBase : public Subscriptor
+ class VectorBase : public ReadVector<PetscScalar>, public Subscriptor
{
public:
/**
* Return the global dimension of the vector.
*/
size_type
- size() const;
+ size() const override;
/**
* Return the local dimension of the vector, i.e. the number of elements
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<PetscScalar> & values) const;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<PetscScalar> & elements) const override;
+
/**
* Instead of getting individual elements of a vector via operator(),
* this function allows getting a whole set of elements at once. In
extract_subvector_to(indices.begin(), indices.end(), values.begin());
}
+ inline void
+ VectorBase::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<PetscScalar> & elements) const
+ {
+ AssertDimension(indices.size(), elements.size());
+ extract_subvector_to(indices.begin(), indices.end(), elements.begin());
+ }
+
+
template <typename ForwardIterator, typename OutputIterator>
inline void
VectorBase::extract_subvector_to(const ForwardIterator indices_begin,
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_lac_read_vector_h
+#define dealii_lac_read_vector_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/array_view.h>
+#include <deal.II/base/types.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * @addtogroup Vectors
+ * @{
+ */
+
+/**
+ * Base class for providing read-only access to vector elements.
+ *
+ * deal.II supports a large number of vector classes, including both its own
+ * serial and parallel vectors as well as vector classes from external
+ * libraries like PETSc and Trilinos. ReadVector is a common base class for
+ * all vector classes and defines a minimal interface for efficiently
+ * accessing vector elements.
+ */
+template <typename Number>
+class ReadVector
+{
+public:
+ using size_type = types::global_dof_index;
+
+ /**
+ * Return the size of the vector.
+ */
+ virtual size_type
+ size() const = 0;
+
+ /**
+ * Extract a subset of the vector specified by @p indices into the output
+ * array @p elements.
+ */
+ virtual void
+ extract_subvector_to(const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> &elements) const = 0;
+};
+
+/** @} */
+
+DEAL_II_NAMESPACE_CLOSE
+#endif
#include <deal.II/base/template_constraints.h>
#include <deal.II/base/types.h>
+#include <deal.II/lac/read_vector.h>
#include <deal.II/lac/vector_operation.h>
#include <cstdlib>
* get the first index of the largest range.
*/
template <typename Number>
- class ReadWriteVector : public Subscriptor
+ class ReadWriteVector : public Subscriptor, public ReadVector<Number>
{
public:
/**
* number returned by the current function.
*/
size_type
- size() const;
+ size() const override;
/**
* This function returns the number of elements stored. It is smaller or
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<Number2> & values) const;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> & entries) const override;
+
/**
* Instead of getting individual elements of a vector via operator(),
* this function allows getting a whole set of elements at once. In
+ template <typename Number>
+ void
+ ReadWriteVector<Number>::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> & entries) const
+ {
+ AssertDimension(indices.size(), entries.size());
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ entries[i] = (*this)[indices[i]];
+ }
+
+
+
template <typename Number>
template <typename ForwardIterator, typename OutputIterator>
inline void
# include <deal.II/base/mpi_stub.h>
# include <deal.II/base/subscriptor.h>
+# include <deal.II/lac/read_vector.h>
# include <deal.II/lac/trilinos_epetra_communication_pattern.h>
# include <deal.II/lac/vector_operation.h>
# include <deal.II/lac/vector_space_vector.h>
* @ingroup TrilinosWrappers
* @ingroup Vectors
*/
- class Vector : public VectorSpaceVector<double>, public Subscriptor
+ class Vector : public VectorSpaceVector<double>,
+ public ReadVector<double>,
+ public Subscriptor
{
public:
+ using value_type = double;
+ using size_type = types::global_dof_index;
+
/**
* Constructor. Create a vector of dimension zero.
*/
reinit(const VectorSpaceVector<double> &V,
const bool omit_zeroing_entries = false) override;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<double> &elements) const override;
+
/**
* Copy function. This function takes a Vector and copies all the
* elements. The Vector will have the same parallel distribution as @p
# include <deal.II/base/mpi_stub.h>
# include <deal.II/base/subscriptor.h>
+# include <deal.II/lac/read_vector.h>
# include <deal.II/lac/trilinos_tpetra_communication_pattern.h>
# include <deal.II/lac/vector_operation.h>
# include <deal.II/lac/vector_space_vector.h>
* @ingroup Vectors
*/
template <typename Number>
- class Vector : public VectorSpaceVector<Number>, public Subscriptor
+ class Vector : public VectorSpaceVector<Number>,
+ public ReadVector<Number>,
+ public Subscriptor
{
public:
using value_type = Number;
reinit(const VectorSpaceVector<Number> &V,
const bool omit_zeroing_entries = false) override;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> &elements) const override;
+
/**
* Copy function. This function takes a Vector and copies all the
* elements. The Vector will have the same parallel distribution as @p
+ template <typename Number>
+ void
+ Vector::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> & elements) const
+ {
+ AssertDimension(indices.size(), elements.size());
+ const auto &vector = trilinos_vector();
+ const auto &map = vector.Map();
+
+# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
+ auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
+ Tpetra::Access::ReadOnly);
+# else
+ vector.template sync<Kokkos::HostSpace>();
+ auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
+# endif
+ auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
+
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ {
+ AssertIndexRange(indices[i], size());
+ const trilinos_i = map->getLocalElement(
+ static_cast<TrilinosWrappers::types::int_type>(indices[i]));
+ elements[i] = vector_1d(trilinos_i);
+ }
+ }
+
+
+
template <typename Number>
Vector<Number> &
Vector<Number>::operator=(const Vector<Number> &V)
# include <deal.II/base/subscriptor.h>
# include <deal.II/lac/exceptions.h>
+# include <deal.II/lac/read_vector.h>
# include <deal.II/lac/vector.h>
# include <deal.II/lac/vector_operation.h>
# include <deal.II/lac/vector_type_traits.h>
* @ingroup TrilinosWrappers
* @ingroup Vectors
*/
- class Vector : public Subscriptor
+ class Vector : public Subscriptor, public ReadVector<TrilinosScalar>
{
public:
/**
* Return the global dimension of the vector.
*/
size_type
- size() const;
+ size() const override;
/**
* Return the local dimension of the vector, i.e. the number of elements
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<TrilinosScalar> & values) const;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(const ArrayView<const size_type> &indices,
+ ArrayView<TrilinosScalar> &elements) const override;
+
/**
* Instead of getting individual elements of a vector via operator(),
* this function allows getting a whole set of elements at once. In
+ inline void
+ Vector::extract_subvector_to(const ArrayView<const size_type> &indices,
+ ArrayView<TrilinosScalar> &elements) const
+ {
+ AssertDimension(indices.size(), elements.size());
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ {
+ AssertIndexRange(indices[i], size());
+ elements[i] = (*this)[indices[i]];
+ }
+ }
+
+
+
template <typename ForwardIterator, typename OutputIterator>
inline void
Vector::extract_subvector_to(ForwardIterator indices_begin,
#include <deal.II/base/numbers.h>
#include <deal.II/base/subscriptor.h>
+#include <deal.II/lac/read_vector.h>
#include <deal.II/lac/vector_operation.h>
#include <deal.II/lac/vector_type_traits.h>
* in the manual).
*/
template <typename Number>
-class Vector : public Subscriptor
+class Vector : public Subscriptor, public ReadVector<Number>
{
public:
/**
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<OtherNumber> & values) const;
+ /**
+ * Extract a range of elements all at once.
+ */
+ virtual void
+ extract_subvector_to(const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> &elements) const override;
+
/**
* Instead of getting individual elements of a vector via operator(),
* this function allows getting a whole set of elements at once. In
/**
* Return dimension of the vector.
*/
- size_type
- size() const;
+ virtual size_type
+ size() const override;
/**
* Return local dimension of the vector. Since this vector does not support
+template <typename Number>
+void
+Vector<Number>::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<Number> & elements) const
+{
+ AssertDimension(indices.size(), elements.size());
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ {
+ AssertIndexRange(indices[i], size());
+ elements[i] = (*this)[indices[i]];
+ }
+}
+
+
+
template <typename Number>
Vector<Number> &
Vector<Number>::operator+=(const Vector<Number> &v)
+ void
+ Vector::extract_subvector_to(
+ const ArrayView<const types::global_dof_index> &indices,
+ ArrayView<double> & elements) const
+ {
+ AssertDimension(indices.size(), elements.size());
+ const auto &vector = trilinos_vector();
+ const auto &map = vector.Map();
+
+ for (unsigned int i = 0; i < indices.size(); ++i)
+ {
+ AssertIndexRange(indices[i], size());
+ const auto trilinos_i =
+ map.LID(static_cast<TrilinosWrappers::types::int_type>(indices[i]));
+ elements[i] = vector[0][trilinos_i];
+ }
+ }
+
+
+
Vector &
Vector::operator=(const Vector &V)
{