Reference documentation for deal.II version 9.0.0
|
#include <deal.II/lac/read_write_vector.h>
Classes | |
class | FunctorTemplate |
Public Types | |
typedef Number | value_type |
Public Member Functions | |
1: Basic Object-handling | |
ReadWriteVector () | |
ReadWriteVector (const ReadWriteVector< Number > &in_vector) | |
ReadWriteVector (const size_type size) | |
ReadWriteVector (const IndexSet &locally_stored_indices) | |
~ReadWriteVector ()=default | |
virtual void | reinit (const size_type size, const bool omit_zeroing_entries=false) |
template<typename Number2 > | |
void | reinit (const ReadWriteVector< Number2 > &in_vector, const bool omit_zeroing_entries=false) |
virtual void | reinit (const IndexSet &locally_stored_indices, const bool omit_zeroing_entries=false) |
void | reinit (const TrilinosWrappers::MPI::Vector &trilinos_vec) |
template<typename Functor > | |
void | apply (const Functor &func) |
void | swap (ReadWriteVector< Number > &v) |
ReadWriteVector< Number > & | operator= (const ReadWriteVector< Number > &in_vector) |
template<typename Number2 > | |
ReadWriteVector< Number > & | operator= (const ReadWriteVector< Number2 > &in_vector) |
ReadWriteVector< Number > & | operator= (const Number s) |
void | import (const distributed::Vector< Number > &vec, VectorOperation::values operation, const std::shared_ptr< const CommunicationPatternBase > &communication_pattern=std::shared_ptr< const CommunicationPatternBase >()) |
void | import (const PETScWrappers::MPI::Vector &petsc_vec, VectorOperation::values operation, const std::shared_ptr< const CommunicationPatternBase > &communication_pattern=std::shared_ptr< const CommunicationPatternBase >()) |
void | import (const TrilinosWrappers::MPI::Vector &trilinos_vec, VectorOperation::values operation, std::shared_ptr< const CommunicationPatternBase > communication_pattern=std::shared_ptr< const CommunicationPatternBase >()) |
void | import (const EpetraWrappers::Vector &epetra_vec, VectorOperation::values operation, std::shared_ptr< const CommunicationPatternBase > communication_pattern=std::shared_ptr< const CommunicationPatternBase >()) |
void | import (const CUDAWrappers::Vector< Number > &cuda_vec, VectorOperation::values operation, std::shared_ptr< const CommunicationPatternBase > communication_pattern=std::shared_ptr< const CommunicationPatternBase >()) |
size_type | size () const |
size_type | n_elements () const |
const IndexSet & | get_stored_elements () const |
iterator | begin () |
const_iterator | begin () const |
iterator | end () |
const_iterator | end () const |
2: Data-Access | |
Number | operator() (const size_type global_index) const |
Number & | operator() (const size_type global_index) |
Number | operator[] (const size_type global_index) const |
Number & | operator[] (const size_type global_index) |
template<typename Number2 > | |
void | extract_subvector_to (const std::vector< size_type > &indices, std::vector< Number2 > &values) const |
template<typename ForwardIterator , typename OutputIterator > | |
void | extract_subvector_to (ForwardIterator indices_begin, const ForwardIterator indices_end, OutputIterator values_begin) const |
Number | local_element (const size_type local_index) const |
Number & | local_element (const size_type local_index) |
3: Modification of vectors | |
template<typename Number2 > | |
void | add (const std::vector< size_type > &indices, const std::vector< Number2 > &values) |
template<typename Number2 > | |
void | add (const std::vector< size_type > &indices, const ReadWriteVector< Number2 > &values) |
template<typename Number2 > | |
void | add (const size_type n_elements, const size_type *indices, const Number2 *values) |
void | print (std::ostream &out, const unsigned int precision=3, const bool scientific=true) const |
std::size_t | memory_consumption () const |
Public Member Functions inherited from Subscriptor | |
Subscriptor () | |
Subscriptor (const Subscriptor &) | |
Subscriptor (Subscriptor &&) noexcept | |
virtual | ~Subscriptor () |
Subscriptor & | operator= (const Subscriptor &) |
Subscriptor & | operator= (Subscriptor &&) noexcept |
void | subscribe (const char *identifier=nullptr) const |
void | unsubscribe (const char *identifier=nullptr) const |
unsigned int | n_subscriptions () const |
void | list_subscribers () const |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
Protected Member Functions | |
void | import (const Epetra_MultiVector &multivector, const IndexSet &locally_owned_elements, VectorOperation::values operation, const MPI_Comm &mpi_comm, const std::shared_ptr< const CommunicationPatternBase > &communication_pattern) |
unsigned int | global_to_local (const types::global_dof_index global_index) const |
void | resize_val (const size_type new_allocated_size) |
EpetraWrappers::CommunicationPattern | create_epetra_comm_pattern (const IndexSet &source_index_set, const MPI_Comm &mpi_comm) |
Protected Attributes | |
IndexSet | stored_elements |
IndexSet | source_stored_elements |
std::shared_ptr< CommunicationPatternBase > | comm_pattern |
std::unique_ptr< Number[], decltype(free) * > | values |
std::shared_ptr< ::parallel::internal::TBBPartitioner > | thread_loop_partitioner |
Friends | |
template<typename Number2 > | |
class | ReadWriteVector |
Additional Inherited Members | |
Static Public Member Functions inherited from Subscriptor | |
static ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
ReadWriteVector is intended to represent vectors in \({\mathbb R}^N\) for which it stores all or a subset of elements. The latter case in important in parallel computations, where \(N\) may be so large that no processor can actually all elements of a solution vector, but where this is also not necessary: one typically only has to store the values of degrees of freedom that live on cells that are locally owned plus potentially those degrees of freedom that live on ghost cells.
This class allows to access individual elements to be read or written. However, it does not allow global operations such as taking the norm. ReadWriteVector can be used to read and write elements in vectors derived from VectorSpaceVector such as TrilinosWrappers::MPI::Vector and PETScWrappers::MPI::Vector.
Most of the time, one will simply read from or write into a vector of the current class using the global numbers of these degrees of freedom. This is done using operator() or operator[] which call global_to_local() to transform the global index into a local one. In such cases, it is clear that one can only access elements of the vector that the current object indeed stores.
However, it is also possible to access elements in the order in which they are stored by the current object. In other words, one is not interested in accessing elements with their global indices, but instead using an enumeration that only takes into account the elements that are actually stored. This is facilitated by the local_element() function. To this end, it is necessary to know in which order the current class stores its element. The elements of all the consecutive ranges are stored in ascending order of the first index of each range. The function largest_range_starting_index() of IndexSet can be used to get the first index of the largest range.
Definition at line 43 of file la_parallel_vector.h.
typedef Number LinearAlgebra::ReadWriteVector< Number >::value_type |
Declare standard types used in all containers. These types parallel those in the C++
standard libraries vector<...>
class.
Definition at line 132 of file read_write_vector.h.
LinearAlgebra::ReadWriteVector< Number >::ReadWriteVector | ( | ) |
Empty constructor.
LinearAlgebra::ReadWriteVector< Number >::ReadWriteVector | ( | const ReadWriteVector< Number > & | in_vector | ) |
Copy constructor.
|
explicit |
Construct a vector given the size, the stored elements have their index in [0,size).
|
explicit |
Construct a vector whose stored elements indices are given by the IndexSet locally_stored_indices
.
|
default |
Destructor.
|
virtual |
Set the global size of the vector to size
. The stored elements have their index in [0,size).
If the flag omit_zeroing_entries
is set to false, the memory will be initialized with zero, otherwise the memory will be untouched (and the user must make sure to fill it with reasonable data before using it).
Reimplemented in LinearAlgebra::Vector< Number >.
void LinearAlgebra::ReadWriteVector< Number >::reinit | ( | const ReadWriteVector< Number2 > & | in_vector, |
const bool | omit_zeroing_entries = false |
||
) |
Uses the same IndexSet as the one of the input vector in_vector
and allocates memory for this vector.
If the flag omit_zeroing_entries
is set to false, the memory will be initialized with zero, otherwise the memory will be untouched (and the user must make sure to fill it with reasonable data before using it).
|
virtual |
Initializes the vector. The indices are specified by locally_stored_indices
.
If the flag omit_zeroing_entries
is set to false, the memory will be initialized with zero, otherwise the memory will be untouched (and the user must make sure to fill it with reasonable data before using it). locally_stored_indices.
Reimplemented in LinearAlgebra::Vector< Number >.
void LinearAlgebra::ReadWriteVector< Number >::reinit | ( | const TrilinosWrappers::MPI::Vector & | trilinos_vec | ) |
Initialize this ReadWriteVector by supplying access to all locally available entries in the given ghosted or non-ghosted vector.
trilinos_vec
.This function is mainly written for backwards-compatibility to get element access to a ghosted TrilinosWrappers::MPI::Vector inside the library.
void LinearAlgebra::ReadWriteVector< Number >::apply | ( | const Functor & | func | ) |
Apply the functor func
to each element of the vector. The functor should look like
void LinearAlgebra::ReadWriteVector< Number >::swap | ( | ReadWriteVector< Number > & | v | ) |
Swap the contents of this vector and the other vector v
. One could do this operation with a temporary variable and copying over the data elements, but this function is significantly more efficient since it only swaps the pointers to the data of the two vectors and therefore does not need to allocate temporary storage and move data around.
This function is analogous to the swap
function of all C++ standard containers. Also, there is a global function swap(u,v)
that simply calls u.swap(v)
, again in analogy to standard functions.
ReadWriteVector<Number>& LinearAlgebra::ReadWriteVector< Number >::operator= | ( | const ReadWriteVector< Number > & | in_vector | ) |
Copies the data and the IndexSet of the input vector in_vector
.
ReadWriteVector<Number>& LinearAlgebra::ReadWriteVector< Number >::operator= | ( | const ReadWriteVector< Number2 > & | in_vector | ) |
Copies the data and the IndexSet of the input vector in_vector
.
ReadWriteVector<Number>& LinearAlgebra::ReadWriteVector< Number >::operator= | ( | const Number | s | ) |
Sets all elements of the vector to the scalar s
. This operation is only allowed if s
is equal to zero.
void LinearAlgebra::ReadWriteVector< Number >::import | ( | const distributed::Vector< Number > & | vec, |
VectorOperation::values | operation, | ||
const std::shared_ptr< const CommunicationPatternBase > & | communication_pattern = std::shared_ptr< const CommunicationPatternBase >() |
||
) |
Imports all the elements present in the vector's IndexSet from the input vector vec
. VectorOperation::values operation
is used to decide if the elements in V
should be added to the current vector or replace the current elements. The last parameter can be used if the same communication pattern is used multiple times. This can be used to improve performance.
void LinearAlgebra::ReadWriteVector< Number >::import | ( | const PETScWrappers::MPI::Vector & | petsc_vec, |
VectorOperation::values | operation, | ||
const std::shared_ptr< const CommunicationPatternBase > & | communication_pattern = std::shared_ptr< const CommunicationPatternBase >() |
||
) |
Imports all the elements present in the vector's IndexSet from the input vector petsc_vec
. VectorOperation::values operation
is used to decide if the elements in V
should be added to the current vector or replace the current elements. The last parameter can be used if the same communication pattern is used multiple times. This can be used to improve performance.
void LinearAlgebra::ReadWriteVector< Number >::import | ( | const TrilinosWrappers::MPI::Vector & | trilinos_vec, |
VectorOperation::values | operation, | ||
std::shared_ptr< const CommunicationPatternBase > | communication_pattern = std::shared_ptr< const CommunicationPatternBase >() |
||
) |
Imports all the elements present in the vector's IndexSet from the input vector trilinos_vec
. VectorOperation::values operation
is used to decide if the elements in V
should be added to the current vector or replace the current elements. The last parameter can be used if the same communication pattern is used multiple times. This can be used to improve performance.
trilinos_vec
is not allowed to have ghost entries. void LinearAlgebra::ReadWriteVector< Number >::import | ( | const EpetraWrappers::Vector & | epetra_vec, |
VectorOperation::values | operation, | ||
std::shared_ptr< const CommunicationPatternBase > | communication_pattern = std::shared_ptr< const CommunicationPatternBase >() |
||
) |
Imports all the elements present in the vector's IndexSet from the input vector epetra_vec
. VectorOperation::values operation
is used to decide if the elements in V
should be added to the current vector or replace the current elements. The last parameter can be used if the same communication pattern is used multiple times. This can be used to improve performance.
void LinearAlgebra::ReadWriteVector< Number >::import | ( | const CUDAWrappers::Vector< Number > & | cuda_vec, |
VectorOperation::values | operation, | ||
std::shared_ptr< const CommunicationPatternBase > | communication_pattern = std::shared_ptr< const CommunicationPatternBase >() |
||
) |
Import all the elements present in the vector's IndexSet from the input vector cuda_vec
. VectorOperation::values operation
is used to decide if the elements in V
should be added to the current vector or replace the current elements. The last parameter is not used.
size_type LinearAlgebra::ReadWriteVector< Number >::size | ( | ) | const |
The value returned by this function denotes the dimension of the vector spaces that are modeled by objects of this kind. However, objects of the current class do not actually stores all elements of vectors of this space but may, in fact store only a subset. The number of elements stored is returned by n_elements() and is smaller or equal to the number returned by the current function.
size_type LinearAlgebra::ReadWriteVector< Number >::n_elements | ( | ) | const |
This function returns the number of elements stored. It is smaller or equal to the dimension of the vector space that is modeled by an object of this kind. This dimension is return by size().
const IndexSet& LinearAlgebra::ReadWriteVector< Number >::get_stored_elements | ( | ) | const |
Return the IndexSet that represents the indices of the elements stored.
iterator LinearAlgebra::ReadWriteVector< Number >::begin | ( | ) |
Make the ReadWriteVector
class a bit like the vector<>
class of the C++ standard library by returning iterators to the start and end of the locally stored elements of this vector.
const_iterator LinearAlgebra::ReadWriteVector< Number >::begin | ( | ) | const |
Return constant iterator to the start of the locally stored elements of the vector.
iterator LinearAlgebra::ReadWriteVector< Number >::end | ( | ) |
Return an iterator pointing to the element past the end of the array of locally stored entries.
const_iterator LinearAlgebra::ReadWriteVector< Number >::end | ( | ) | const |
Return a constant iterator pointing to the element past the end of the array of the locally stored entries.
Number LinearAlgebra::ReadWriteVector< Number >::operator() | ( | const size_type | global_index | ) | const |
Read access to the data in the position corresponding to global_index
. An exception is thrown if global_index
is not stored by the current object.
Number& LinearAlgebra::ReadWriteVector< Number >::operator() | ( | const size_type | global_index | ) |
Read and write access to the data in the position corresponding to global_index
. An exception is thrown if global_index
is not stored by the current object.
Number LinearAlgebra::ReadWriteVector< Number >::operator[] | ( | const size_type | global_index | ) | const |
Read access to the data in the position corresponding to global_index
. An exception is thrown if global_index
is not stored by the current object.
This function does the same thing as operator().
Number& LinearAlgebra::ReadWriteVector< Number >::operator[] | ( | const size_type | global_index | ) |
Read and write access to the data in the position corresponding to global_index
. An exception is thrown if global_index
is not stored by the current object.
This function does the same thing as operator().
void LinearAlgebra::ReadWriteVector< Number >::extract_subvector_to | ( | const std::vector< size_type > & | indices, |
std::vector< Number2 > & | values | ||
) | const |
Instead of getting individual elements of a vector via operator(), this function allows getting a whole set of elements at once. The indices of the elements to be read are stated in the first argument, the corresponding values are returned in the second.
If the current vector is called v
, then this function is the equivalent to the code
indices
and values
arrays must be identical. void LinearAlgebra::ReadWriteVector< Number >::extract_subvector_to | ( | ForwardIterator | indices_begin, |
const ForwardIterator | indices_end, | ||
OutputIterator | values_begin | ||
) | const |
Instead of getting individual elements of a vector via operator(), this function allows getting a whole set of elements at once. In contrast to the previous function, this function obtains the indices of the elements by dereferencing all elements of the iterator range provided by the first two arguments, and puts the vector values into memory locations obtained by dereferencing a range of iterators starting at the location pointed to by the third argument.
If the current vector is called v
, then this function is the equivalent to the code
values_begin
as there are iterators between indices_begin
and indices_end
. Number LinearAlgebra::ReadWriteVector< Number >::local_element | ( | const size_type | local_index | ) | const |
Read access to the data field specified by local_index
. When you access elements in the order in which they are stored, it is necessary that you know in which they are stored. In other words, you need to know the map between the global indices of the elements this class stores, and the local indices into the contiguous array of these global elements. For this, see the general documentation of this class.
Performance: Direct array access (fast).
Number& LinearAlgebra::ReadWriteVector< Number >::local_element | ( | const size_type | local_index | ) |
Read and write access to the data field specified by local_index
. When you access elements in the order in which they are stored, it is necessary that you know in which they are stored. In other words, you need to know the map between the global indices of the elements this class stores, and the local indices into the contiguous array of these global elements. For this, see the general documentation of this class.
Performance: Direct array access (fast).
void LinearAlgebra::ReadWriteVector< Number >::add | ( | const std::vector< size_type > & | indices, |
const std::vector< Number2 > & | values | ||
) |
This function adds a whole set of values stored in values
to the vector components specified by indices
.
void LinearAlgebra::ReadWriteVector< Number >::add | ( | const std::vector< size_type > & | indices, |
const ReadWriteVector< Number2 > & | values | ||
) |
This function is similar to the previous one but takes a ReadWriteVector of values.
void LinearAlgebra::ReadWriteVector< Number >::add | ( | const size_type | n_elements, |
const size_type * | indices, | ||
const Number2 * | values | ||
) |
Take an address where n_elements
are stored contiguously and add them into the vector. Handles all cases which are not covered by the other two add()
functions above.
void LinearAlgebra::ReadWriteVector< Number >::print | ( | std::ostream & | out, |
const unsigned int | precision = 3 , |
||
const bool | scientific = true |
||
) | const |
Prints the vector to the output stream out
.
std::size_t LinearAlgebra::ReadWriteVector< Number >::memory_consumption | ( | ) | const |
Return the memory consumption of this class in bytes.
|
protected |
Import all the elements present in the vector's IndexSet from the input vector multivector
. This is an helper function and it should not be used directly.
|
inlineprotected |
Return the local position of global_index
.
Definition at line 573 of file read_write_vector.h.
|
protected |
A helper function that is used to resize the val array.
|
protected |
Return a EpetraWrappers::Communication pattern and store it for future use.
Make all other ReadWriteVector types friends.
Definition at line 625 of file read_write_vector.h.
|
protected |
Indices of the elements stored.
Definition at line 598 of file read_write_vector.h.
|
protected |
IndexSet of the elements of the last imported vector;
Definition at line 603 of file read_write_vector.h.
|
protected |
CommunicationPattern for the communication between the source_stored_elements IndexSet and the current vector.
Definition at line 609 of file read_write_vector.h.
|
protected |
Pointer to the array of local elements of this vector.
Definition at line 614 of file read_write_vector.h.
|
mutableprotected |
For parallel loops with TBB, this member variable stores the affinity information of loops.
Definition at line 620 of file read_write_vector.h.