// ---------------------------------------------------------------------
//
-// Copyright (C) 2011 - 2019 by the deal.II authors
+// Copyright (C) 2020 by the deal.II authors
//
// This file is part of the deal.II library.
//
{
namespace SharedMPI
{
- /*! @addtogroup Vectors
- *@{
- */
-
- /**
- * Implementation of a parallel vector class. The design of this class is
- * similar to the standard ::dealii::Vector class in deal.II, with the
- * exception that storage is SharedMPI with MPI.
- *
- * The vector is designed for the following scheme of parallel
- * partitioning:
- * <ul>
- * <li> The indices held by individual processes (locally owned part) in
- * the MPI parallelization form a contiguous range
- * <code>[my_first_index,my_last_index)</code>.
- * <li> Ghost indices residing on arbitrary positions of other processors
- * are allowed. It is in general more efficient if ghost indices are
- * clustered, since they are stored as a set of intervals. The
- * communication pattern of the ghost indices is determined when calling
- * the function <code>reinit (locally_owned, ghost_indices,
- * communicator)</code>, and retained until the partitioning is changed.
- * This allows for efficient parallel communication of indices. In
- * particular, it stores the communication pattern, rather than having to
- * compute it again for every communication. For more information on ghost
- * vectors, see also the
- * @ref GlossGhostedVector "glossary entry on vectors with ghost elements".
- * <li> Besides the usual global access operator() it is also possible to
- * access vector entries in the local index space with the function @p
- * local_element(). Locally owned indices are placed first, [0,
- * local_size()), and then all ghost indices follow after them
- * contiguously, [local_size(), local_size()+n_ghost_entries()).
- * </ul>
- *
- * Functions related to parallel functionality:
- * <ul>
- * <li> The function <code>compress()</code> goes through the data
- * associated with ghost indices and communicates it to the owner process,
- * which can then add it to the correct position. This can be used e.g.
- * after having run an assembly routine involving ghosts that fill this
- * vector. Note that the @p insert mode of @p compress() does not set the
- * elements included in ghost entries but simply discards them, assuming
- * that the owning processor has set them to the desired value already
- * (See also the
- * @ref GlossCompress "glossary entry on compress").
- * <li> The <code>update_ghost_values()</code> function imports the data
- * from the owning processor to the ghost indices in order to provide read
- * access to the data associated with ghosts.
- * <li> It is possible to split the above functions into two phases, where
- * the first initiates the communication and the second one finishes it.
- * These functions can be used to overlap communication with computations
- * in other parts of the code.
- * <li> Of course, reduction operations (like norms) make use of
- * collective all-to-all MPI communications.
- * </ul>
- *
- * This vector can take two different states with respect to ghost
- * elements:
- * <ul>
- * <li> After creation and whenever zero_out_ghosts() is called (or
- * <code>operator= (0.)</code>), the vector does only allow writing into
- * ghost elements but not reading from ghost elements.
- * <li> After a call to update_ghost_values(), the vector does not allow
- * writing into ghost elements but only reading from them. This is to
- * avoid undesired ghost data artifacts when calling compress() after
- * modifying some vector entries. The current status of the ghost entries
- * (read mode or write mode) can be queried by the method
- * has_ghost_elements(), which returns <code>true</code> exactly when
- * ghost elements have been updated and <code>false</code> otherwise,
- * irrespective of the actual number of ghost entries in the vector layout
- * (for that information, use n_ghost_entries() instead).
- * </ul>
- *
- * This vector uses the facilities of the class dealii::Vector<Number> for
- * implementing the operations on the local range of the vector. In
- * particular, it also inherits thread parallelism that splits most
- * vector-vector operations into smaller chunks if the program uses
- * multiple threads. This may or may not be desired when working also with
- * MPI.
- *
- * <h4>Limitations regarding the vector size</h4>
- *
- * This vector class is based on two different number types for indexing.
- * The so-called global index type encodes the overall size of the vector.
- * Its type is types::global_dof_index. The largest possible value is
- * <code>2^32-1</code> or approximately 4 billion in case 64 bit integers
- * are disabled at configuration of deal.II (default case) or
- * <code>2^64-1</code> or approximately <code>10^19</code> if 64 bit
- * integers are enabled (see the glossary entry on
- * @ref GlobalDoFIndex
- * for further information).
- *
- * The second relevant index type is the local index used within one MPI
- * rank. As opposed to the global index, the implementation assumes 32-bit
- * unsigned integers unconditionally. In other words, to actually use a
- * vector with more than four billion entries, you need to use MPI with
- * more than one rank (which in general is a safe assumption since four
- * billion entries consume at least 16 GB of memory for floats or 32 GB of
- * memory for doubles) and enable 64-bit indices. If more than 4 billion
- * local elements are present, the implementation tries to detect that,
- * which triggers an exception and aborts the code. Note, however, that
- * the detection of overflow is tricky and the detection mechanism might
- * fail in some circumstances. Therefore, it is strongly recommended to
- * not rely on this class to automatically detect the unsupported case.
- *
- * <h4>CUDA support</h4>
- *
- * This vector class supports two different memory spaces: Host and CUDA. By
- * default, the memory space is Host and all the data are allocated on the
- * CPU. When the memory space is CUDA, all the data is allocated on the GPU.
- * The operations on the vector are performed on the chosen memory space. *
- * From the host, there are two methods to access the elements of the Vector
- * when using the CUDA memory space:
- * <ul>
- * <li> use get_values():
- * @code
- * Vector<double, MemorySpace::CUDA> vector(local_range, comm);
- * double* vector_dev = vector.get_values();
- * std::vector<double> vector_host(local_range.n_elements(), 1.);
- * Utilities::CUDA::copy_to_dev(vector_host, vector_dev);
- * @endcode
- * <li> use import():
- * @code
- * Vector<double, MemorySpace::CUDA> vector(local_range, comm);
- * ReadWriteVector<double> rw_vector(local_range);
- * for (auto & val : rw_vector)
- * val = 1.;
- * vector.import(rw_vector, VectorOperations::insert);
- * @endcode
- * </ul>
- * The import method is a lot safer and will perform an MPI communication if
- * necessary. Since an MPI communication may be performed, import needs to
- * be called on all the processors.
- *
- * @note By default, all the ranks will try to access the device 0. This is
- * fine is if you have one rank per node and one gpu per node. If you
- * have multiple GPUs on one node, we need each process to access a
- * different GPU. If each node has the same number of GPUs, this can be done
- * as follows:
- * <code> int n_devices = 0; cudaGetDeviceCount(&n_devices); int
- * device_id = my_rank % n_devices;
- * cudaSetDevice(device_id);
- * </code>
- * @see CUDAWrappers
- *
- * @author Katharina Kormann, Martin Kronbichler, Bruno Turcksin 2010, 2011,
- * 2016, 2018
- */
template <typename Number, typename MemorySpace = MemorySpace::Host>
class Vector : public ::dealii::LinearAlgebra::VectorSpaceVector<Number>,
public Subscriptor
std::is_same<MemorySpace, ::dealii::MemorySpace::CUDA>::value,
"MemorySpace should be Host or CUDA");
- /**
- * @name 1: Basic Object-handling
- */
- //@{
- /**
- * Empty constructor.
- */
Vector();
- /**
- * Copy constructor. Uses the parallel partitioning of @p in_vector.
- * It should be noted that this constructor automatically sets ghost
- * values to zero. Call @p update_ghost_values() directly following
- * construction if a ghosted vector is required.
- */
Vector(const Vector<Number, MemorySpace> &in_vector);
- /**
- * Construct a parallel vector of the given global size without any
- * actual parallel distribution.
- */
Vector(const size_type size);
- /**
- * Construct a parallel vector. The local range is specified by @p
- * locally_owned_set (note that this must be a contiguous interval,
- * multiple intervals are not possible). The IndexSet @p ghost_indices
- * specifies ghost indices, i.e., indices which one might need to read
- * data from or accumulate data from. It is allowed that the set of
- * ghost indices also contains the local range, but it does not need to.
- *
- * This function involves global communication, so it should only be
- * called once for a given layout. Use the constructor with
- * Vector<Number> argument to create additional vectors with the same
- * parallel layout.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
Vector(const IndexSet &local_range,
const IndexSet &ghost_indices,
const MPI_Comm communicator);
- /**
- * Same constructor as above but without any ghost indices.
- */
Vector(const IndexSet &local_range, const MPI_Comm communicator);
- /**
- * Create the vector based on the parallel partitioning described in @p
- * partitioner. The input argument is a shared pointer, which store the
- * partitioner data only once and share it between several vectors with
- * the same layout.
- */
Vector(
const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
- /**
- * Destructor.
- */
virtual ~Vector() override;
- /**
- * Set the global size of the vector to @p size without any actual
- * parallel distribution.
- */
void
reinit(const size_type size, const bool omit_zeroing_entries = false);
- /**
- * Uses the parallel layout of the input vector @p in_vector and
- * allocates memory for this vector. Recommended initialization function
- * when several vectors with the same layout should be created.
- *
- * If the flag @p 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).
- */
template <typename Number2>
void
reinit(const Vector<Number2, MemorySpace> &in_vector,
const bool omit_zeroing_entries = false);
- /**
- * Initialize the vector. The local range is specified by @p
- * locally_owned_set (note that this must be a contiguous interval,
- * multiple intervals are not possible). The IndexSet @p ghost_indices
- * specifies ghost indices, i.e., indices which one might need to read
- * data from or accumulate data from. It is allowed that the set of
- * ghost indices also contains the local range, but it does not need to.
- *
- * This function involves global communication, so it should only be
- * called once for a given layout. Use the @p reinit function with
- * Vector<Number> argument to create additional vectors with the same
- * parallel layout.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
void
reinit(const IndexSet &local_range,
const IndexSet &ghost_indices,
const MPI_Comm communicator);
- /**
- * Same as above, but without ghost entries.
- */
void
reinit(const IndexSet &local_range, const MPI_Comm communicator);
- /**
- * Initialize the vector given to the parallel partitioning described in
- * @p partitioner. The input argument is a shared pointer, which store
- * the partitioner data only once and share it between several vectors
- * with the same layout.
- */
void
reinit(
const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
- /**
- * Swap the contents of this vector and the other vector @p 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 @p swap function of all C++
- * standard containers. Also, there is a global function
- * <tt>swap(u,v)</tt> that simply calls <tt>u.swap(v)</tt>, again in
- * analogy to standard functions.
- *
- * This function is virtual in order to allow for derived classes to
- * handle memory separately.
- */
void
swap(Vector<Number, MemorySpace> &v);
- /**
- * Assigns the vector to the parallel partitioning of the input vector
- * @p in_vector, and copies all the data.
- *
- * If one of the input vector or the calling vector (to the left of the
- * assignment operator) had ghost elements set before this operation,
- * the calling vector will have ghost values set. Otherwise, it will be
- * in write mode. If the input vector does not have any ghost elements
- * at all, the vector will also update its ghost values in analogy to
- * the respective setting the Trilinos and PETSc vectors.
- */
Vector<Number, MemorySpace> &
operator=(const Vector<Number, MemorySpace> &in_vector);
- /**
- * Assigns the vector to the parallel partitioning of the input vector
- * @p in_vector, and copies all the data.
- *
- * If one of the input vector or the calling vector (to the left of the
- * assignment operator) had ghost elements set before this operation,
- * the calling vector will have ghost values set. Otherwise, it will be
- * in write mode. If the input vector does not have any ghost elements
- * at all, the vector will also update its ghost values in analogy to
- * the respective setting the Trilinos and PETSc vectors.
- */
template <typename Number2>
Vector<Number, MemorySpace> &
operator=(const Vector<Number2, MemorySpace> &in_vector);
-#ifdef DEAL_II_WITH_PETSC
- /**
- * Copy the content of a PETSc vector into the calling vector. This
- * function assumes that the vectors layouts have already been
- * initialized to match.
- *
- * This operator is only available if deal.II was configured with PETSc.
- *
- * This function is deprecated. Use the interface through
- * ReadWriteVector instead.
- */
- DEAL_II_DEPRECATED
- Vector<Number, MemorySpace> &
- operator=(const PETScWrappers::MPI::Vector &petsc_vec);
-#endif
-
-#ifdef DEAL_II_WITH_TRILINOS
- /**
- * Copy the content of a Trilinos vector into the calling vector. This
- * function assumes that the vectors layouts have already been
- * initialized to match.
- *
- * This operator is only available if deal.II was configured with
- * Trilinos.
- *
- * This function is deprecated. Use the interface through
- * ReadWriteVector instead.
- */
- DEAL_II_DEPRECATED
- Vector<Number, MemorySpace> &
- operator=(const TrilinosWrappers::MPI::Vector &trilinos_vec);
-#endif
- //@}
-
- /**
- * @name 2: Parallel data exchange
- */
- //@{
- /**
- * This function copies the data that has accumulated in the data buffer
- * for ghost indices to the owning processor. For the meaning of the
- * argument @p operation, see the entry on
- * @ref GlossCompress "Compressing SharedMPI vectors and matrices"
- * in the glossary.
- *
- * There are four variants for this function. If called with argument @p
- * VectorOperation::add adds all the data accumulated in ghost elements
- * to the respective elements on the owning processor and clears the
- * ghost array afterwards. If called with argument @p
- * VectorOperation::insert, a set operation is performed. Since setting
- * elements in a vector with ghost elements is ambiguous (as one can set
- * both the element on the ghost site as well as the owning site), this
- * operation makes the assumption that all data is set correctly on the
- * owning processor. Upon call of compress(VectorOperation::insert), all
- * ghost entries are thus simply zeroed out (using zero_ghost_values()).
- * In debug mode, a check is performed for whether the data set is
- * actually consistent between processors, i.e., whenever a non-zero
- * ghost element is found, it is compared to the value on the owning
- * processor and an exception is thrown if these elements do not agree.
- * If called with VectorOperation::min or VectorOperation::max, the
- * minimum or maximum on all elements across the processors is set.
- * @note This vector class has a fixed set of ghost entries attached to
- * the local representation. As a consequence, all ghost entries are
- * assumed to be valid and will be exchanged unconditionally according
- * to the given VectorOperation. Make sure to initialize all ghost
- * entries with the neutral element of the given VectorOperation or
- * touch all ghost entries. The neutral element is zero for
- * VectorOperation::add and VectorOperation::insert, `+inf` for
- * VectorOperation::min, and `-inf` for VectorOperation::max. If all
- * values are initialized with values below zero and compress is called
- * with VectorOperation::max two times subsequently, the maximal value
- * after the second calculation will be zero.
- */
virtual void
compress(::dealii::VectorOperation::values operation) override;
- /**
- * Fills the data field for ghost indices with the values stored in the
- * respective positions of the owning processor. This function is needed
- * before reading from ghosts. The function is @p const even though
- * ghost data is changed. This is needed to allow functions with a @p
- * const vector to perform the data exchange without creating
- * temporaries.
- *
- * After calling this method, write access to ghost elements of the
- * vector is forbidden and an exception is thrown. Only read access to
- * ghost elements is allowed in this state. Note that all subsequent
- * operations on this vector, like global vector addition, etc., will
- * also update the ghost values by a call to this method after the
- * operation. However, global reduction operations like norms or the
- * inner product will always ignore ghost elements in order to avoid
- * counting the ghost data more than once. To allow writing to ghost
- * elements again, call zero_out_ghosts().
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
void
update_ghost_values() const;
- /**
- * Initiates communication for the @p compress() function with non-
- * blocking communication. This function does not wait for the transfer
- * to finish, in order to allow for other computations during the time
- * it takes until all data arrives.
- *
- * Before the data is actually exchanged, the function must be followed
- * by a call to @p compress_finish().
- *
- * In case this function is called for more than one vector before @p
- * compress_finish() is invoked, it is mandatory to specify a unique
- * communication channel to each such call, in order to avoid several
- * messages with the same ID that will corrupt this operation. Any
- * communication channel less than 100 is a valid value (in particular,
- * the range $[100, 200)$ is reserved for
- * LinearAlgebra::SharedMPI::BlockVector).
- */
void
compress_start(
const unsigned int communication_channel = 0,
::dealii::VectorOperation::values operation = VectorOperation::add);
- /**
- * For all requests that have been initiated in compress_start, wait for
- * the communication to finish. Once it is finished, add or set the data
- * (depending on the flag operation) to the respective positions in the
- * owning processor, and clear the contents in the ghost data fields.
- * The meaning of this argument is the same as in compress().
- *
- * This function should be called exactly once per vector after calling
- * compress_start, otherwise the result is undefined. In particular, it
- * is not well-defined to call compress_start on the same vector again
- * before compress_finished has been called. However, there is no
- * warning to prevent this situation.
- *
- * Must follow a call to the @p compress_start function.
- */
void
compress_finish(::dealii::VectorOperation::values operation);
- /**
- * Initiates communication for the @p update_ghost_values() function
- * with non-blocking communication. This function does not wait for the
- * transfer to finish, in order to allow for other computations during
- * the time it takes until all data arrives.
- *
- * Before the data is actually exchanged, the function must be followed
- * by a call to @p update_ghost_values_finish().
- *
- * In case this function is called for more than one vector before @p
- * update_ghost_values_finish() is invoked, it is mandatory to specify a
- * unique communication channel to each such call, in order to avoid
- * several messages with the same ID that will corrupt this operation.
- * Any communication channel less than 100 is a valid value (in
- * particular, the range $[100, 200)$ is reserved for
- * LinearAlgebra::SharedMPI::BlockVector).
- */
void
update_ghost_values_start(
const unsigned int communication_channel = 0) const;
-
- /**
- * For all requests that have been started in update_ghost_values_start,
- * wait for the communication to finish.
- *
- * Must follow a call to the @p update_ghost_values_start function
- * before reading data from ghost indices.
- */
void
update_ghost_values_finish() const;
- /**
- * This method zeros the entries on ghost dofs, but does not touch
- * locally owned DoFs.
- *
- * After calling this method, read access to ghost elements of the
- * vector is forbidden and an exception is thrown. Only write access to
- * ghost elements is allowed in this state.
- */
void
zero_out_ghosts() const;
- /**
- * Return whether the vector currently is in a state where ghost values
- * can be read or not. This is the same functionality as other parallel
- * vectors have. If this method returns false, this only means that
- * read-access to ghost elements is prohibited whereas write access is
- * still possible (to those entries specified as ghosts during
- * initialization), not that there are no ghost elements at all.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
bool
has_ghost_elements() const;
- /**
- * This method copies the data in the locally owned range from another
- * SharedMPI vector @p src into the calling vector. As opposed to
- * operator= that also includes ghost entries, this operation ignores
- * the ghost range. The only prerequisite is that the local range on the
- * calling vector and the given vector @p src are the same on all
- * processors. It is explicitly allowed that the two vectors have
- * different ghost elements that might or might not be related to each
- * other.
- *
- * Since no data exchange is performed, make sure that neither @p src
- * nor the calling vector have pending communications in order to obtain
- * correct results.
- */
template <typename Number2>
void
copy_locally_owned_data_from(const Vector<Number2, MemorySpace> &src);
- /**
- * Import all the elements present in the SharedMPI vector @p src.
- * VectorOperation::values @p operation is used to decide if the elements
- * in @p V should be added to the current vector or replace the current
- * elements. The main purpose of this function is to get data from one
- * memory space, e.g. CUDA, to the other, e.g. the Host.
- *
- * @note The partitioners of the two SharedMPI vectors need to be the
- * same as no MPI communication is performed.
- */
template <typename MemorySpace2>
void
import(const Vector<Number, MemorySpace2> &src,
VectorOperation::values operation);
- //@}
-
- /**
- * @name 3: Implementation of VectorSpaceVector
- */
- //@{
-
- /**
- * Change the dimension to that of the vector V. The elements of V are not
- * copied.
- */
virtual void
reinit(const VectorSpaceVector<Number> &V,
const bool omit_zeroing_entries = false) override;
- /**
- * Multiply the entire vector by a fixed factor.
- */
virtual Vector<Number, MemorySpace> &
operator*=(const Number factor) override;
- /**
- * Divide the entire vector by a fixed factor.
- */
virtual Vector<Number, MemorySpace> &
operator/=(const Number factor) override;
- /**
- * Add the vector @p V to the present one.
- */
virtual Vector<Number, MemorySpace> &
operator+=(const VectorSpaceVector<Number> &V) override;
- /**
- * Subtract the vector @p V from the present one.
- */
virtual Vector<Number, MemorySpace> &
operator-=(const VectorSpaceVector<Number> &V) override;
- /**
- * Import all the elements present in the vector's IndexSet from the input
- * vector @p V. VectorOperation::values @p operation is used to decide if
- * the elements in @p 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.
- *
- * @note If the MemorySpace is CUDA, the data in the ReadWriteVector will
- * be moved to the device.
- */
virtual void
import(
const LinearAlgebra::ReadWriteVector<Number> & V,
std::shared_ptr<const CommunicationPatternBase> communication_pattern =
std::shared_ptr<const CommunicationPatternBase>()) override;
- /**
- * Return the scalar product of two vectors.
- */
virtual Number
operator*(const VectorSpaceVector<Number> &V) const override;
- /**
- * Add @p a to all components. Note that @p a is a scalar not a vector.
- */
virtual void
add(const Number a) override;
- /**
- * Simple addition of a multiple of a vector, i.e. <tt>*this += a*V</tt>.
- */
virtual void
add(const Number a, const VectorSpaceVector<Number> &V) override;
- /**
- * Multiple addition of scaled vectors, i.e. <tt>*this += a*V+b*W</tt>.
- */
virtual void
add(const Number a,
const VectorSpaceVector<Number> &V,
const Number b,
const VectorSpaceVector<Number> &W) override;
- /**
- * A collective add operation: This function adds a whole set of values
- * stored in @p values to the vector components specified by @p indices.
- */
virtual void
add(const std::vector<size_type> &indices,
const std::vector<Number> & values);
- /**
- * Scaling and simple addition of a multiple of a vector, i.e. <tt>*this =
- * s*(*this)+a*V</tt>.
- */
virtual void
sadd(const Number s,
const Number a,
const VectorSpaceVector<Number> &V) override;
- /**
- * Scale each element of this vector by the corresponding element in the
- * argument. This function is mostly meant to simulate multiplication (and
- * immediate re-assignment) by a diagonal scaling matrix.
- */
virtual void
scale(const VectorSpaceVector<Number> &scaling_factors) override;
- /**
- * Assignment <tt>*this = a*V</tt>.
- */
virtual void
equ(const Number a, const VectorSpaceVector<Number> &V) override;
- /**
- * Return the l<sub>1</sub> norm of the vector (i.e., the sum of the
- * absolute values of all entries among all processors).
- */
virtual real_type
l1_norm() const override;
- /**
- * Return the $l_2$ norm of the vector (i.e., the square root of
- * the sum of the square of all entries among all processors).
- */
virtual real_type
l2_norm() const override;
- /**
- * Return the square of the $l_2$ norm of the vector.
- */
real_type
norm_sqr() const;
- /**
- * Return the maximum norm of the vector (i.e., the maximum absolute value
- * among all entries and among all processors).
- */
virtual real_type
linfty_norm() const override;
- /**
- * Perform a combined operation of a vector addition and a subsequent
- * inner product, returning the value of the inner product. In other
- * words, the result of this function is the same as if the user called
- * @code
- * this->add(a, V);
- * return_value = *this * W;
- * @endcode
- *
- * The reason this function exists is that this operation involves less
- * memory transfer than calling the two functions separately. This method
- * only needs to load three vectors, @p this, @p V, @p W, whereas calling
- * separate methods means to load the calling vector @p this twice. Since
- * most vector operations are memory transfer limited, this reduces the
- * time by 25\% (or 50\% if @p W equals @p this).
- *
- * For complex-valued vectors, the scalar product in the second step is
- * implemented as
- * $\left<v,w\right>=\sum_i v_i \bar{w_i}$.
- */
virtual Number
add_and_dot(const Number a,
const VectorSpaceVector<Number> &V,
const VectorSpaceVector<Number> &W) override;
- /**
- * Return the global size of the vector, equal to the sum of the number of
- * locally owned indices among all processors.
- */
virtual size_type
size() const override;
- /**
- * Return an index set that describes which elements of this vector are
- * owned by the current processor. As a consequence, the index sets
- * returned on different processors if this is a SharedMPI vector will
- * form disjoint sets that add up to the complete index set. Obviously, if
- * a vector is created on only one processor, then the result would
- * satisfy
- * @code
- * vec.locally_owned_elements() == complete_index_set(vec.size())
- * @endcode
- */
virtual dealii::IndexSet
locally_owned_elements() const override;
- /**
- * Print the vector to the output stream @p out.
- */
virtual void
print(std::ostream & out,
const unsigned int precision = 3,
const bool scientific = true,
const bool across = true) const override;
- /**
- * Return the memory consumption of this class in bytes.
- */
virtual std::size_t
memory_consumption() const override;
- //@}
-
- /**
- * @name 4: Other vector operations not included in VectorSpaceVector
- */
- //@{
-
- /**
- * Sets all elements of the vector to the scalar @p s. If the scalar is
- * zero, also ghost elements are set to zero, otherwise they remain
- * unchanged.
- */
+
virtual Vector<Number, MemorySpace> &
operator=(const Number s) override;
- /**
- * This is a collective add operation that adds a whole set of values
- * stored in @p values to the vector components specified by @p indices.
- */
template <typename OtherNumber>
void
add(const std::vector<size_type> & indices,
const ::dealii::Vector<OtherNumber> &values);
- /**
- * Take an address where n_elements are stored contiguously and add them
- * into the vector.
- */
template <typename OtherNumber>
void
add(const size_type n_elements,
const size_type * indices,
const OtherNumber *values);
- /**
- * Scaling and simple vector addition, i.e. <tt>*this =
- * s*(*this)+V</tt>.
- */
void
sadd(const Number s, const Vector<Number, MemorySpace> &V);
- /**
- * Scaling and multiple addition.
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
void
sadd(const Number s,
const Number b,
const Vector<Number, MemorySpace> &W);
- /**
- * Assignment <tt>*this = a*u + b*v</tt>.
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
void
equ(const Number a,
const Number b,
const Vector<Number, MemorySpace> &v);
- //@}
-
-
- /**
- * @name 5: Entry access and local data representation
- */
- //@{
-
- /**
- * Return the local size of the vector, i.e., the number of indices
- * owned locally.
- */
size_type
local_size() const;
- /**
- * Return the half-open interval that specifies the locally owned range
- * of the vector. Note that <code>local_size() == local_range().second -
- * local_range().first</code>.
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
std::pair<size_type, size_type>
local_range() const;
- /**
- * Return true if the given global index is in the local range of this
- * processor.
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
bool
in_local_range(const size_type global_index) const;
- /**
- * Return the number of ghost elements present on the vector.
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
size_type
n_ghost_entries() const;
- /**
- * Return an index set that describes which elements of this vector are
- * not owned by the current processor but can be written into or read
- * from locally (ghost elements).
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
const IndexSet &
ghost_elements() const;
- /**
- * Return whether the given global index is a ghost index on the
- * present processor. Returns false for indices that are owned locally
- * and for indices not present at all.
- *
- * This function is deprecated.
- */
DEAL_II_DEPRECATED
bool
is_ghost_entry(const types::global_dof_index global_index) const;
- /**
- * Make the @p Vector class a bit like the <tt>vector<></tt> class of
- * the C++ standard library by returning iterators to the start and end
- * of the <i>locally owned</i> elements of this vector.
- *
- * It holds that end() - begin() == local_size().
- *
- * @note For the CUDA memory space, the iterator points to memory on the
- * device.
- */
iterator
begin();
- /**
- * Return constant iterator to the start of the locally owned elements
- * of the vector.
- *
- * @note For the CUDA memory space, the iterator points to memory on the
- * device.
- */
const_iterator
begin() const;
- /**
- * Return an iterator pointing to the element past the end of the array
- * of locally owned entries.
- *
- * @note For the CUDA memory space, the iterator points to memory on the
- * device.
- */
iterator
end();
- /**
- * Return a constant iterator pointing to the element past the end of
- * the array of the locally owned entries.
- *
- * @note For the CUDA memory space, the iterator points to memory on the
- * device.
- */
const_iterator
end() const;
- /**
- * Read access to the data in the position corresponding to @p
- * global_index. The index must be either in the local range of the
- * vector or be specified as a ghost index at construction.
- *
- * Performance: <tt>O(1)</tt> for locally owned elements that represent
- * a contiguous range and <tt>O(log(n<sub>ranges</sub>))</tt> for ghost
- * elements (quite fast, but slower than local_element()).
- */
Number
operator()(const size_type global_index) const;
- /**
- * Read and write access to the data in the position corresponding to @p
- * global_index. The index must be either in the local range of the
- * vector or be specified as a ghost index at construction.
- *
- * Performance: <tt>O(1)</tt> for locally owned elements that represent
- * a contiguous range and <tt>O(log(n<sub>ranges</sub>))</tt> for ghost
- * elements (quite fast, but slower than local_element()).
- */
Number &
operator()(const size_type global_index);
- /**
- * Read access to the data in the position corresponding to @p
- * global_index. The index must be either in the local range of the
- * vector or be specified as a ghost index at construction.
- *
- * This function does the same thing as operator().
- */
Number operator[](const size_type global_index) const;
- /**
- * Read and write access to the data in the position corresponding to @p
- * global_index. The index must be either in the local range of the
- * vector or be specified as a ghost index at construction.
- *
- * This function does the same thing as operator().
- */
+
Number &operator[](const size_type global_index);
- /**
- * Read access to the data field specified by @p local_index. Locally
- * owned indices can be accessed with indices
- * <code>[0,local_size)</code>, and ghost indices with indices
- * <code>[local_size,local_size+ n_ghost_entries]</code>.
- *
- * Performance: Direct array access (fast).
- */
Number
local_element(const size_type local_index) const;
- /**
- * Read and write access to the data field specified by @p local_index.
- * Locally owned indices can be accessed with indices
- * <code>[0,local_size)</code>, and ghost indices with indices
- * <code>[local_size,local_size+n_ghosts]</code>.
- *
- * Performance: Direct array access (fast).
- */
Number &
local_element(const size_type local_index);
- /**
- * Return the pointer to the underlying raw array.
- *
- * @note For the CUDA memory space, the pointer points to memory on the
- * device.
- */
Number *
get_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 @p v, then this function is the equivalent
- * to the code
- * @code
- * for (unsigned int i=0; i<indices.size(); ++i)
- * values[i] = v[indices[i]];
- * @endcode
- *
- * @pre The sizes of the @p indices and @p values arrays must be identical.
- *
- * @note This function is not implemented for CUDA memory space.
- */
template <typename OtherNumber>
void
extract_subvector_to(const std::vector<size_type> &indices,
std::vector<OtherNumber> & values) 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 @p v, then this function is the equivalent
- * to the code
- * @code
- * ForwardIterator indices_p = indices_begin;
- * OutputIterator values_p = values_begin;
- * while (indices_p != indices_end)
- * {
- * *values_p = v[*indices_p];
- * ++indices_p;
- * ++values_p;
- * }
- * @endcode
- *
- * @pre It must be possible to write into as many memory locations
- * starting at @p values_begin as there are iterators between
- * @p indices_begin and @p indices_end.
- */
template <typename ForwardIterator, typename OutputIterator>
void
extract_subvector_to(ForwardIterator indices_begin,
const ForwardIterator indices_end,
OutputIterator values_begin) const;
- /**
- * Return whether the vector contains only elements with value zero.
- * This is a collective operation. This function is expensive, because
- * potentially all elements have to be checked.
- */
+
virtual bool
all_zero() const override;
- /**
- * Compute the mean value of all the entries in the vector.
- */
virtual Number
mean_value() const override;
- /**
- * $l_p$-norm of the vector. The pth root of the sum of the pth powers
- * of the absolute values of the elements.
- */
real_type
lp_norm(const real_type p) const;
- //@}
-
- /**
- * @name 6: Mixed stuff
- */
- //@{
- /**
- * Return a reference to the MPI communicator object in use with this
- * vector.
- */
const MPI_Comm &
get_mpi_communicator() const;
- /**
- * Return the MPI partitioner that describes the parallel layout of the
- * vector. This object can be used to initialize another vector with the
- * respective reinit() call, for additional queries regarding the
- * parallel communication, or the compatibility of partitioners.
- */
const std::shared_ptr<const Utilities::MPI::Partitioner> &
get_partitioner() const;
- /**
- * Check whether the given partitioner is compatible with the
- * partitioner used for this vector. Two partitioners are compatible if
- * they have the same local size and the same ghost indices. They do not
- * necessarily need to be the same data field of the shared pointer.
- * This is a local operation only, i.e., if only some processors decide
- * that the partitioning is not compatible, only these processors will
- * return @p false, whereas the other processors will return @p true.
- */
bool
partitioners_are_compatible(
const Utilities::MPI::Partitioner &part) const;
- /**
- * Check whether the given partitioner is compatible with the
- * partitioner used for this vector. Two partitioners are compatible if
- * they have the same local size and the same ghost indices. They do not
- * necessarily need to be the same data field. As opposed to
- * partitioners_are_compatible(), this method checks for compatibility
- * among all processors and the method only returns @p true if the
- * partitioner is the same on all processors.
- *
- * This method performs global communication, so make sure to use it
- * only in a context where all processors call it the same number of
- * times.
- */
bool
partitioners_are_globally_compatible(
const Utilities::MPI::Partitioner &part) const;
- /**
- * Change the ghost state of this vector to @p ghosted.
- */
void
set_ghost_state(const bool ghosted) const;
- //@}
-
- /**
- * Attempt to perform an operation between two incompatible vector types.
- *
- * @ingroup Exceptions
- */
DeclException0(ExcVectorTypeNotCompatible);
- /**
- * Attempt to perform an operation not implemented on the device.
- *
- * @ingroup Exceptions
- */
DeclException0(ExcNotAllowedForCuda);
- /**
- * Exception
- */
DeclException3(ExcNonMatchingElements,
Number,
Number,
<< std::setprecision(16) << arg2
<< " on the owner processor " << arg3);
- /**
- * Exception
- */
DeclException4(ExcAccessToNonLocalElement,
size_type,
size_type,
<< "that this vector can access.");
private:
- /**
- * Simple addition of a multiple of a vector, i.e. <tt>*this += a*V</tt>
- * without MPI communication.
- */
void
add_local(const Number a, const VectorSpaceVector<Number> &V);
- /**
- * Scaling and simple addition of a multiple of a vector, i.e. <tt>*this =
- * s*(*this)+a*V</tt> without MPI communication.
- */
void
sadd_local(const Number s,
const Number a,
const VectorSpaceVector<Number> &V);
- /**
- * Local part of the inner product of two vectors.
- */
template <typename Number2>
Number
inner_product_local(const Vector<Number2, MemorySpace> &V) const;
- /**
- * Local part of norm_sqr().
- */
real_type
norm_sqr_local() const;
- /**
- * Local part of mean_value().
- */
Number
mean_value_local() const;
- /**
- * Local part of l1_norm().
- */
real_type
l1_norm_local() const;
- /**
- * Local part of lp_norm().
- */
real_type
lp_norm_local(const real_type p) const;
- /**
- * Local part of linfty_norm().
- */
real_type
linfty_norm_local() const;
- /**
- * Local part of the addition followed by an inner product of two
- * vectors. The same applies for complex-valued vectors as for
- * the add_and_dot() function.
- */
Number
add_and_dot_local(const Number a,
const Vector<Number, MemorySpace> &V,
const Vector<Number, MemorySpace> &W);
- /**
- * Shared pointer to store the parallel partitioning information. This
- * information can be shared between several vectors that have the same
- * partitioning.
- */
std::shared_ptr<const Utilities::MPI::Partitioner> partitioner;
- /**
- * The size that is currently allocated in the val array.
- */
size_type allocated_size;
- /**
- * Underlying data structure storing the local elements of this vector.
- */
mutable ::dealii::MemorySpace::MemorySpaceData<Number, MemorySpace> data;
- /**
- * For parallel loops with TBB, this member variable stores the affinity
- * information of loops.
- */
mutable std::shared_ptr<::dealii::parallel::internal::TBBPartitioner>
thread_loop_partitioner;
- /**
- * Temporary storage that holds the data that is sent to this processor
- * in @p compress() or sent from this processor in
- * @p update_ghost_values.
- */
mutable ::dealii::MemorySpace::MemorySpaceData<Number, MemorySpace>
import_data;
- /**
- * Stores whether the vector currently allows for reading ghost elements
- * or not. Note that this is to ensure consistent ghost data and does
- * not indicate whether the vector actually can store ghost elements. In
- * particular, when assembling a vector we do not allow reading
- * elements, only writing them.
- */
mutable bool vector_is_ghosted;
-#ifdef DEAL_II_WITH_MPI
- /**
- * A vector that collects all requests from @p compress() operations.
- * This class uses persistent MPI communicators, i.e., the communication
- * channels are stored during successive calls to a given function. This
- * reduces the overhead involved with setting up the MPI machinery, but
- * it does not remove the need for a receive operation to be posted
- * before the data can actually be sent.
- */
- std::vector<MPI_Request> compress_requests;
-
- /**
- * A vector that collects all requests from @p update_ghost_values()
- * operations. This class uses persistent MPI communicators.
- */
- mutable std::vector<MPI_Request> update_ghost_values_requests;
-#endif
-
- /**
- * A lock that makes sure that the @p compress and @p
- * update_ghost_values functions give reasonable results also when used
- * with several threads.
- */
mutable std::mutex mutex;
- /**
- * A helper function that clears the compress_requests and
- * update_ghost_values_requests field. Used in reinit functions.
- */
void
clear_mpi_requests();
- /**
- * A helper function that is used to resize the val array.
- */
void
resize_val(const size_type new_allocated_size);
- // Make all other vector types friends.
template <typename Number2, typename MemorySpace2>
friend class Vector;
- // Make BlockVector type friends.
template <typename Number2>
friend class BlockVector;
};
- /*@}*/
/*-------------------- Inline functions ---------------------------------*/
static inline typename Vector<Number, MemorySpace>::iterator
begin(::dealii::MemorySpace::MemorySpaceData<Number, MemorySpace> &)
{
+ Assert(false, ExcNotImplemented());
return nullptr;
}
begin(
const ::dealii::MemorySpace::MemorySpaceData<Number, MemorySpace> &)
{
+ Assert(false, ExcNotImplemented());
return nullptr;
}
get_values(
::dealii::MemorySpace::MemorySpaceData<Number, MemorySpace> &)
{
+ Assert(false, ExcNotImplemented());
return nullptr;
}
};
begin(::dealii::MemorySpace::
MemorySpaceData<Number, ::dealii::MemorySpace::Host> &data)
{
+ Assert(false, ExcNotImplemented());
return data.values.get();
}
begin(const ::dealii::MemorySpace::
MemorySpaceData<Number, ::dealii::MemorySpace::Host> &data)
{
+ Assert(false, ExcNotImplemented());
return data.values.get();
}
get_values(::dealii::MemorySpace::
MemorySpaceData<Number, ::dealii::MemorySpace::Host> &data)
{
+ Assert(false, ExcNotImplemented());
return data.values.get();
}
};
begin(::dealii::MemorySpace::
MemorySpaceData<Number, ::dealii::MemorySpace::CUDA> &data)
{
+ Assert(false, ExcNotImplemented());
return data.values_dev.get();
}
begin(const ::dealii::MemorySpace::
MemorySpaceData<Number, ::dealii::MemorySpace::CUDA> &data)
{
+ Assert(false, ExcNotImplemented());
return data.values_dev.get();
}
get_values(::dealii::MemorySpace::
MemorySpaceData<Number, ::dealii::MemorySpace::CUDA> &data)
{
+ Assert(false, ExcNotImplemented());
return data.values_dev.get();
}
};
inline bool
Vector<Number, MemorySpace>::has_ghost_elements() const
{
- return vector_is_ghosted;
+ Assert(false, ExcNotImplemented());
+ return false;
}
inline typename Vector<Number, MemorySpace>::size_type
Vector<Number, MemorySpace>::size() const
{
- return partitioner->size();
+ Assert(false, ExcNotImplemented());
+ return 0;
}
inline typename Vector<Number, MemorySpace>::size_type
Vector<Number, MemorySpace>::local_size() const
{
- return partitioner->local_size();
+ Assert(false, ExcNotImplemented());
+ return 0;
}
typename Vector<Number, MemorySpace>::size_type>
Vector<Number, MemorySpace>::local_range() const
{
- return partitioner->local_range();
+ Assert(false, ExcNotImplemented());
+ return {0, 0};
}
Vector<Number, MemorySpace>::in_local_range(
const size_type global_index) const
{
- return partitioner->in_local_range(global_index);
+ Assert(false, ExcNotImplemented());
+ (void)global_index;
+ return true;
}
inline IndexSet
Vector<Number, MemorySpace>::locally_owned_elements() const
{
- IndexSet is(size());
-
- is.add_range(partitioner->local_range().first,
- partitioner->local_range().second);
+ Assert(false, ExcNotImplemented());
+ IndexSet is;
return is;
}
inline typename Vector<Number, MemorySpace>::size_type
Vector<Number, MemorySpace>::n_ghost_entries() const
{
- return partitioner->n_ghost_indices();
+ Assert(false, ExcNotImplemented());
+ return 0;
}
inline const IndexSet &
Vector<Number, MemorySpace>::ghost_elements() const
{
+ Assert(false, ExcNotImplemented());
return partitioner->ghost_indices();
}
Vector<Number, MemorySpace>::is_ghost_entry(
const size_type global_index) const
{
- return partitioner->is_ghost_entry(global_index);
+ Assert(false, ExcNotImplemented());
+ (void)global_index;
+ return false;
}
inline typename Vector<Number, MemorySpace>::iterator
Vector<Number, MemorySpace>::begin()
{
+ Assert(false, ExcNotImplemented());
return internal::Policy<Number, MemorySpace>::begin(data);
}
inline typename Vector<Number, MemorySpace>::const_iterator
Vector<Number, MemorySpace>::begin() const
{
+ Assert(false, ExcNotImplemented());
return internal::Policy<Number, MemorySpace>::begin(data);
}
inline typename Vector<Number, MemorySpace>::iterator
Vector<Number, MemorySpace>::end()
{
+ Assert(false, ExcNotImplemented());
return internal::Policy<Number, MemorySpace>::begin(data) +
partitioner->local_size();
}
inline typename Vector<Number, MemorySpace>::const_iterator
Vector<Number, MemorySpace>::end() const
{
+ Assert(false, ExcNotImplemented());
return internal::Policy<Number, MemorySpace>::begin(data) +
partitioner->local_size();
}
inline Number
Vector<Number, MemorySpace>::operator()(const size_type global_index) const
{
- Assert((std::is_same<MemorySpace, ::dealii::MemorySpace::Host>::value),
- ExcMessage(
- "This function is only implemented for the Host memory space"));
- Assert(
- partitioner->in_local_range(global_index) ||
- partitioner->ghost_indices().is_element(global_index),
- ExcAccessToNonLocalElement(global_index,
- partitioner->local_range().first,
- partitioner->local_range().second,
- partitioner->ghost_indices().n_elements()));
- // do not allow reading a vector which is not in ghost mode
- Assert(partitioner->in_local_range(global_index) ||
- vector_is_ghosted == true,
- ExcMessage("You tried to read a ghost element of this vector, "
- "but it has not imported its ghost values."));
+ Assert(false, ExcNotImplemented());
return data.values[partitioner->global_to_local(global_index)];
}
inline Number &
Vector<Number, MemorySpace>::operator()(const size_type global_index)
{
- Assert((std::is_same<MemorySpace, ::dealii::MemorySpace::Host>::value),
- ExcMessage(
- "This function is only implemented for the Host memory space"));
- Assert(
- partitioner->in_local_range(global_index) ||
- partitioner->ghost_indices().is_element(global_index),
- ExcAccessToNonLocalElement(global_index,
- partitioner->local_range().first,
- partitioner->local_range().second,
- partitioner->ghost_indices().n_elements()));
- // we would like to prevent reading ghosts from a vector that does not
- // have them imported, but this is not possible because we might be in a
- // part of the code where the vector has enabled ghosts but is non-const
- // (then, the compiler picks this method according to the C++ rule book
- // even if a human would pick the const method when this subsequent use
- // is just a read)
+ Assert(false, ExcNotImplemented());
return data.values[partitioner->global_to_local(global_index)];
}
inline Number Vector<Number, MemorySpace>::
operator[](const size_type global_index) const
{
+ Assert(false, ExcNotImplemented());
return operator()(global_index);
}
inline Number &Vector<Number, MemorySpace>::
operator[](const size_type global_index)
{
+ Assert(false, ExcNotImplemented());
return operator()(global_index);
}
Vector<Number, MemorySpace>::local_element(
const size_type local_index) const
{
- Assert((std::is_same<MemorySpace, ::dealii::MemorySpace::Host>::value),
- ExcMessage(
- "This function is only implemented for the Host memory space"));
- AssertIndexRange(local_index,
- partitioner->local_size() +
- partitioner->n_ghost_indices());
- // do not allow reading a vector which is not in ghost mode
- Assert(local_index < local_size() || vector_is_ghosted == true,
- ExcMessage("You tried to read a ghost element of this vector, "
- "but it has not imported its ghost values."));
+ Assert(false, ExcNotImplemented());
return data.values[local_index];
}
inline Number &
Vector<Number, MemorySpace>::local_element(const size_type local_index)
{
- Assert((std::is_same<MemorySpace, ::dealii::MemorySpace::Host>::value),
- ExcMessage(
- "This function is only implemented for the Host memory space"));
-
- AssertIndexRange(local_index,
- partitioner->local_size() +
- partitioner->n_ghost_indices());
+ Assert(false, ExcNotImplemented());
return data.values[local_index];
}
inline Number *
Vector<Number, MemorySpace>::get_values() const
{
+ Assert(false, ExcNotImplemented());
+
return internal::Policy<Number, MemorySpace>::get_values(data);
}
const std::vector<size_type> &indices,
std::vector<OtherNumber> & values) const
{
- for (size_type i = 0; i < indices.size(); ++i)
- values[i] = operator()(indices[i]);
+ Assert(false, ExcNotImplemented());
+ (void)indices;
+ (void)values;
}
const ForwardIterator indices_end,
OutputIterator values_begin) const
{
- while (indices_begin != indices_end)
- {
- *values_begin = operator()(*indices_begin);
- indices_begin++;
- values_begin++;
- }
+ Assert(false, ExcNotImplemented());
+ (void)indices_begin;
+ (void)indices_end;
+ (void)values_begin;
}
const std::vector<size_type> & indices,
const ::dealii::Vector<OtherNumber> &values)
{
- AssertDimension(indices.size(), values.size());
- for (size_type i = 0; i < indices.size(); ++i)
- {
- Assert(
- numbers::is_finite(values[i]),
- ExcMessage(
- "The given value is not finite but either infinite or Not A Number (NaN)"));
- this->operator()(indices[i]) += values(i);
- }
+ Assert(false, ExcNotImplemented());
+ (void)values;
}
const size_type * indices,
const OtherNumber *values)
{
- for (size_type i = 0; i < n_elements; ++i, ++indices, ++values)
- {
- Assert(
- numbers::is_finite(*values),
- ExcMessage(
- "The given value is not finite but either infinite or Not A Number (NaN)"));
- this->operator()(*indices) += *values;
- }
+ Assert(false, ExcNotImplemented());
+ (void)n_elements;
+ (void)indices;
+ (void)values;
}
inline const MPI_Comm &
Vector<Number, MemorySpace>::get_mpi_communicator() const
{
+ Assert(false, ExcNotImplemented());
return partitioner->get_mpi_communicator();
}
inline const std::shared_ptr<const Utilities::MPI::Partitioner> &
Vector<Number, MemorySpace>::get_partitioner() const
{
+ Assert(false, ExcNotImplemented());
return partitioner;
}
inline void
Vector<Number, MemorySpace>::set_ghost_state(const bool ghosted) const
{
- vector_is_ghosted = ghosted;
+ Assert(false, ExcNotImplemented());
+ (void)ghosted;
}
#endif
} // namespace LinearAlgebra
-/**
- * Global function @p swap which overloads the default implementation of the
- * C++ standard library which uses a temporary object. The function simply
- * exchanges the data of the two vectors.
- *
- * @relatesalso Vector
- * @author Katharina Kormann, Martin Kronbichler, 2011
- */
template <typename Number, typename MemorySpace>
inline void
swap(LinearAlgebra::SharedMPI::Vector<Number, MemorySpace> &u,
LinearAlgebra::SharedMPI::Vector<Number, MemorySpace> &v)
{
- u.swap(v);
+ Assert(false, ExcNotImplemented());
+ (void)u;
+ (void)v;
}
-/**
- * Declare dealii::LinearAlgebra::Vector< Number > as SharedMPI vector.
- *
- * @author Uwe Koecher, 2017
- */
template <typename Number, typename MemorySpace>
struct is_serial_vector<LinearAlgebra::SharedMPI::Vector<Number, MemorySpace>>
: std::false_type
template <typename>
class ReinitHelper;
- /**
- * A helper class used internally in linear_operator.h. Specialization for
- * LinearAlgebra::SharedMPI::Vector<Number>.
- */
template <typename Number>
class ReinitHelper<LinearAlgebra::SharedMPI::Vector<Number>>
{
LinearAlgebra::SharedMPI::Vector<Number> &v,
bool omit_zeroing_entries)
{
- matrix.initialize_dof_vector(v);
- if (!omit_zeroing_entries)
- v = Number();
+ Assert(false, ExcNotImplemented());
+ (void)matrix;
+ (void)v;
+ (void)omit_zeroing_entries;
}
template <typename Matrix>
LinearAlgebra::SharedMPI::Vector<Number> &v,
bool omit_zeroing_entries)
{
- matrix.initialize_dof_vector(v);
- if (!omit_zeroing_entries)
- v = Number();
+ Assert(false, ExcNotImplemented());
+ (void)matrix;
+ (void)v;
+ (void)omit_zeroing_entries;
}
};
// ---------------------------------------------------------------------
//
-// Copyright (C) 2011 - 2019 by the deal.II authors
+// Copyright (C) 2020 by the deal.II authors
//
// This file is part of the deal.II library.
//
{
namespace SharedMPI
{
- namespace internal
- {
- // In the import_from_ghosted_array_finish we might need to calculate the
- // maximal and minimal value for the given number type, which is not
- // straightforward for complex numbers. Therefore, comparison of complex
- // numbers is prohibited and throws an exception.
- template <typename Number>
- Number
- get_min(const Number a, const Number b)
- {
- return std::min(a, b);
- }
-
- template <typename Number>
- std::complex<Number>
- get_min(const std::complex<Number> a, const std::complex<Number>)
- {
- AssertThrow(false,
- ExcMessage("VectorOperation::min not "
- "implemented for complex numbers"));
- return a;
- }
-
- template <typename Number>
- Number
- get_max(const Number a, const Number b)
- {
- return std::max(a, b);
- }
-
- template <typename Number>
- std::complex<Number>
- get_max(const std::complex<Number> a, const std::complex<Number>)
- {
- AssertThrow(false,
- ExcMessage("VectorOperation::max not "
- "implemented for complex numbers"));
- return a;
- }
-
-
-
- // Resize the underlying array on the host or on the device
- template <typename Number, typename MemorySpaceType>
- struct la_parallel_vector_templates_functions
- {
- static_assert(std::is_same<MemorySpaceType, MemorySpace::Host>::value ||
- std::is_same<MemorySpaceType, MemorySpace::CUDA>::value,
- "MemorySpace should be Host or CUDA");
-
- static void
- resize_val(
- const types::global_dof_index /*new_alloc_size*/,
- types::global_dof_index & /*allocated_size*/,
- ::dealii::MemorySpace::MemorySpaceData<Number, MemorySpaceType>
- & /*data*/)
- {}
-
- static void
- import(
- const ::dealii::LinearAlgebra::ReadWriteVector<Number> & /*V*/,
- ::dealii::VectorOperation::values /*operation*/,
- const std::shared_ptr<const ::dealii::Utilities::MPI::Partitioner> &
- /*communication_pattern*/,
- const IndexSet & /*locally_owned_elem*/,
- ::dealii::MemorySpace::MemorySpaceData<Number, MemorySpaceType>
- & /*data*/)
- {}
-
- template <typename RealType>
- static void
- linfty_norm_local(
- const ::dealii::MemorySpace::MemorySpaceData<Number, MemorySpaceType>
- & /*data*/,
- const unsigned int /*size*/,
- RealType & /*max*/)
- {}
- };
-
- template <typename Number>
- struct la_parallel_vector_templates_functions<Number,
- ::dealii::MemorySpace::Host>
- {
- using size_type = types::global_dof_index;
-
- static void
- resize_val(const types::global_dof_index new_alloc_size,
- types::global_dof_index & allocated_size,
- ::dealii::MemorySpace::
- MemorySpaceData<Number, ::dealii::MemorySpace::Host> &data)
- {
- if (new_alloc_size > allocated_size)
- {
- Assert(((allocated_size > 0 && data.values != nullptr) ||
- data.values == nullptr),
- ExcInternalError());
-
- Number *new_val;
- Utilities::System::posix_memalign(
- reinterpret_cast<void **>(&new_val),
- 64,
- sizeof(Number) * new_alloc_size);
- data.values.reset(new_val);
-
- allocated_size = new_alloc_size;
- }
- else if (new_alloc_size == 0)
- {
- data.values.reset();
- allocated_size = 0;
- }
- }
-
- static void
- import(
- const ::dealii::LinearAlgebra::ReadWriteVector<Number> &V,
- ::dealii::VectorOperation::values operation,
- const std::shared_ptr<const ::dealii::Utilities::MPI::Partitioner>
- & communication_pattern,
- const IndexSet &locally_owned_elem,
- ::dealii::MemorySpace::MemorySpaceData<Number,
- ::dealii::MemorySpace::Host>
- &data)
- {
- Assert(
- (operation == ::dealii::VectorOperation::add) ||
- (operation == ::dealii::VectorOperation::insert),
- ExcMessage(
- "Only VectorOperation::add and VectorOperation::insert are allowed"));
-
- ::dealii::LinearAlgebra::SharedMPI::
- Vector<Number, ::dealii::MemorySpace::Host>
- tmp_vector(communication_pattern);
-
- // fill entries from ReadWriteVector into the SharedMPI vector,
- // including ghost entries. this is not really efficient right now
- // because indices are translated twice, once by nth_index_in_set(i)
- // and once for operator() of tmp_vector
- const IndexSet &v_stored = V.get_stored_elements();
- for (size_type i = 0; i < v_stored.n_elements(); ++i)
- tmp_vector(v_stored.nth_index_in_set(i)) = V.local_element(i);
-
- tmp_vector.compress(operation);
-
- // Copy the local elements of tmp_vector to the right place in val
- IndexSet tmp_index_set = tmp_vector.locally_owned_elements();
- if (operation == VectorOperation::add)
- {
- for (size_type i = 0; i < tmp_index_set.n_elements(); ++i)
- {
- data.values[locally_owned_elem.index_within_set(
- tmp_index_set.nth_index_in_set(i))] +=
- tmp_vector.local_element(i);
- }
- }
- else
- {
- for (size_type i = 0; i < tmp_index_set.n_elements(); ++i)
- {
- data.values[locally_owned_elem.index_within_set(
- tmp_index_set.nth_index_in_set(i))] =
- tmp_vector.local_element(i);
- }
- }
- }
-
- template <typename RealType>
- static void
- linfty_norm_local(const ::dealii::MemorySpace::MemorySpaceData<
- Number,
- ::dealii::MemorySpace::Host> &data,
- const unsigned int size,
- RealType & max)
- {
- for (size_type i = 0; i < size; ++i)
- max =
- std::max(numbers::NumberTraits<Number>::abs(data.values[i]), max);
- }
- };
-
-#ifdef DEAL_II_COMPILER_CUDA_AWARE
- template <typename Number>
- struct la_parallel_vector_templates_functions<Number,
- ::dealii::MemorySpace::CUDA>
- {
- using size_type = types::global_dof_index;
-
- static void
- resize_val(const types::global_dof_index new_alloc_size,
- types::global_dof_index & allocated_size,
- ::dealii::MemorySpace::
- MemorySpaceData<Number, ::dealii::MemorySpace::CUDA> &data)
- {
- static_assert(
- std::is_same<Number, float>::value ||
- std::is_same<Number, double>::value,
- "Number should be float or double for CUDA memory space");
-
- if (new_alloc_size > allocated_size)
- {
- Assert(((allocated_size > 0 && data.values_dev != nullptr) ||
- data.values_dev == nullptr),
- ExcInternalError());
-
- Number *new_val_dev;
- Utilities::CUDA::malloc(new_val_dev, new_alloc_size);
- data.values_dev.reset(new_val_dev);
-
- allocated_size = new_alloc_size;
- }
- else if (new_alloc_size == 0)
- {
- data.values_dev.reset();
- allocated_size = 0;
- }
- }
-
- static void
- import(const ReadWriteVector<Number> &V,
- VectorOperation::values operation,
- std::shared_ptr<const Utilities::MPI::Partitioner>
- communication_pattern,
- const IndexSet &locally_owned_elem,
- ::dealii::MemorySpace::
- MemorySpaceData<Number, ::dealii::MemorySpace::CUDA> &data)
- {
- Assert(
- (operation == ::dealii::VectorOperation::add) ||
- (operation == ::dealii::VectorOperation::insert),
- ExcMessage(
- "Only VectorOperation::add and VectorOperation::insert are allowed"));
-
- ::dealii::LinearAlgebra::SharedMPI::
- Vector<Number, ::dealii::MemorySpace::CUDA>
- tmp_vector(communication_pattern);
-
- // fill entries from ReadWriteVector into the SharedMPI vector,
- // including ghost entries. this is not really efficient right now
- // because indices are translated twice, once by nth_index_in_set(i)
- // and once for operator() of tmp_vector
- const IndexSet & v_stored = V.get_stored_elements();
- const size_type n_elements = v_stored.n_elements();
- std::vector<size_type> indices(n_elements);
- for (size_type i = 0; i < n_elements; ++i)
- indices[i] = communication_pattern->global_to_local(
- v_stored.nth_index_in_set(i));
- // Move the indices to the device
- size_type *indices_dev;
- ::dealii::Utilities::CUDA::malloc(indices_dev, n_elements);
- ::dealii::Utilities::CUDA::copy_to_dev(indices, indices_dev);
- // Move the data to the device
- Number *V_dev;
- ::dealii::Utilities::CUDA::malloc(V_dev, n_elements);
- cudaError_t cuda_error_code = cudaMemcpy(V_dev,
- V.begin(),
- n_elements * sizeof(Number),
- cudaMemcpyHostToDevice);
- AssertCuda(cuda_error_code);
-
- // Set the values in tmp_vector
- const int n_blocks =
- 1 + n_elements / (::dealii::CUDAWrappers::chunk_size *
- ::dealii::CUDAWrappers::block_size);
- ::dealii::LinearAlgebra::CUDAWrappers::kernel::set_permutated<Number>
- <<<n_blocks, ::dealii::CUDAWrappers::block_size>>>(
- indices_dev, tmp_vector.begin(), V_dev, n_elements);
-
- tmp_vector.compress(operation);
-
- // Copy the local elements of tmp_vector to the right place in val
- IndexSet tmp_index_set = tmp_vector.locally_owned_elements();
- const size_type tmp_n_elements = tmp_index_set.n_elements();
- indices.resize(tmp_n_elements);
- for (size_type i = 0; i < tmp_n_elements; ++i)
- indices[i] = locally_owned_elem.index_within_set(
- tmp_index_set.nth_index_in_set(i));
- ::dealii::Utilities::CUDA::free(indices_dev);
- ::dealii::Utilities::CUDA::malloc(indices_dev, tmp_n_elements);
- ::dealii::Utilities::CUDA::copy_to_dev(indices, indices_dev);
-
- if (operation == VectorOperation::add)
- ::dealii::LinearAlgebra::CUDAWrappers::kernel::add_permutated<
- Number><<<n_blocks, ::dealii::CUDAWrappers::block_size>>>(
- indices_dev,
- data.values_dev.get(),
- tmp_vector.begin(),
- tmp_n_elements);
- else
- ::dealii::LinearAlgebra::CUDAWrappers::kernel::set_permutated<
- Number><<<n_blocks, ::dealii::CUDAWrappers::block_size>>>(
- indices_dev,
- data.values_dev.get(),
- tmp_vector.begin(),
- tmp_n_elements);
-
- ::dealii::Utilities::CUDA::free(indices_dev);
- ::dealii::Utilities::CUDA::free(V_dev);
- }
-
- template <typename RealType>
- static void
- linfty_norm_local(const ::dealii::MemorySpace::MemorySpaceData<
- Number,
- ::dealii::MemorySpace::CUDA> &data,
- const unsigned int size,
- RealType & result)
- {
- static_assert(std::is_same<Number, RealType>::value,
- "RealType should be the same type as Number");
-
- Number * result_device;
- cudaError_t error_code = cudaMalloc(&result_device, sizeof(Number));
- AssertCuda(error_code);
- error_code = cudaMemset(result_device, 0, sizeof(Number));
-
- const int n_blocks = 1 + size / (::dealii::CUDAWrappers::chunk_size *
- ::dealii::CUDAWrappers::block_size);
- ::dealii::LinearAlgebra::CUDAWrappers::kernel::reduction<
- Number,
- ::dealii::LinearAlgebra::CUDAWrappers::kernel::LInfty<Number>>
- <<<dim3(n_blocks, 1), dim3(::dealii::CUDAWrappers::block_size)>>>(
- result_device, data.values_dev.get(), size);
-
- // Copy the result back to the host
- error_code = cudaMemcpy(&result,
- result_device,
- sizeof(Number),
- cudaMemcpyDeviceToHost);
- AssertCuda(error_code);
- // Free the memory on the device
- error_code = cudaFree(result_device);
- AssertCuda(error_code);
- }
- };
-#endif
- } // namespace internal
-
-
template <typename Number, typename MemorySpaceType>
void
Vector<Number, MemorySpaceType>::clear_mpi_requests()
{
-#ifdef DEAL_II_WITH_MPI
- for (size_type j = 0; j < compress_requests.size(); j++)
- {
- const int ierr = MPI_Request_free(&compress_requests[j]);
- AssertThrowMPI(ierr);
- }
- compress_requests.clear();
- for (size_type j = 0; j < update_ghost_values_requests.size(); j++)
- {
- const int ierr = MPI_Request_free(&update_ghost_values_requests[j]);
- AssertThrowMPI(ierr);
- }
- update_ghost_values_requests.clear();
-#endif
+ Assert(false, ExcNotImplemented());
}
void
Vector<Number, MemorySpaceType>::resize_val(const size_type new_alloc_size)
{
- internal::la_parallel_vector_templates_functions<
- Number,
- MemorySpaceType>::resize_val(new_alloc_size, allocated_size, data);
-
- thread_loop_partitioner =
- std::make_shared<::dealii::parallel::internal::TBBPartitioner>();
+ Assert(false, ExcNotImplemented());
+ (void)new_alloc_size;
}
Vector<Number, MemorySpaceType>::reinit(const size_type size,
const bool omit_zeroing_entries)
{
- clear_mpi_requests();
-
- // check whether we need to reallocate
- resize_val(size);
-
- // delete previous content in import data
- import_data.values.reset();
- import_data.values_dev.reset();
-
- // set partitioner to serial version
- partitioner = std::make_shared<Utilities::MPI::Partitioner>(size);
-
- // set entries to zero if so requested
- if (omit_zeroing_entries == false)
- this->operator=(Number());
- else
- zero_out_ghosts();
+ Assert(false, ExcNotImplemented());
+ (void)size;
+ (void)omit_zeroing_entries;
}
const Vector<Number2, MemorySpaceType> &v,
const bool omit_zeroing_entries)
{
- clear_mpi_requests();
- Assert(v.partitioner.get() != nullptr, ExcNotInitialized());
-
- // check whether the partitioners are
- // different (check only if the are allocated
- // differently, not if the actual data is
- // different)
- if (partitioner.get() != v.partitioner.get())
- {
- partitioner = v.partitioner;
- const size_type new_allocated_size =
- partitioner->local_size() + partitioner->n_ghost_indices();
- resize_val(new_allocated_size);
- }
-
- if (omit_zeroing_entries == false)
- this->operator=(Number());
- else
- zero_out_ghosts();
-
- // do not reallocate import_data directly, but only upon request. It
- // is only used as temporary storage for compress() and
- // update_ghost_values, and we might have vectors where we never
- // call these methods and hence do not need to have the storage.
- import_data.values.reset();
- import_data.values_dev.reset();
-
- thread_loop_partitioner = v.thread_loop_partitioner;
+ Assert(false, ExcNotImplemented());
+ (void)v;
+ (void)omit_zeroing_entries;
}
const IndexSet &ghost_indices,
const MPI_Comm communicator)
{
- // set up parallel partitioner with index sets and communicator
- std::shared_ptr<const Utilities::MPI::Partitioner> new_partitioner(
- new Utilities::MPI::Partitioner(locally_owned_indices,
- ghost_indices,
- communicator));
- reinit(new_partitioner);
+ Assert(false, ExcNotImplemented());
+ (void)locally_owned_indices;
+ (void)ghost_indices;
+ (void)communicator;
}
const IndexSet &locally_owned_indices,
const MPI_Comm communicator)
{
- // set up parallel partitioner with index sets and communicator
- std::shared_ptr<const Utilities::MPI::Partitioner> new_partitioner(
- new Utilities::MPI::Partitioner(locally_owned_indices, communicator));
- reinit(new_partitioner);
+ Assert(false, ExcNotImplemented());
+ (void)locally_owned_indices;
+ (void)communicator;
}
Vector<Number, MemorySpaceType>::reinit(
const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner_in)
{
- clear_mpi_requests();
- partitioner = partitioner_in;
-
- // set vector size and allocate memory
- const size_type new_allocated_size =
- partitioner->local_size() + partitioner->n_ghost_indices();
- resize_val(new_allocated_size);
-
- // initialize to zero
- this->operator=(Number());
-
-
- // do not reallocate import_data directly, but only upon request. It
- // is only used as temporary storage for compress() and
- // update_ghost_values, and we might have vectors where we never
- // call these methods and hence do not need to have the storage.
- import_data.values.reset();
- import_data.values_dev.reset();
-
- vector_is_ghosted = false;
+ Assert(false, ExcNotImplemented());
+ (void)partitioner_in;
}
: partitioner(new Utilities::MPI::Partitioner())
, allocated_size(0)
{
- reinit(0);
+ Assert(false, ExcNotImplemented());
}
, allocated_size(0)
, vector_is_ghosted(false)
{
- reinit(v, true);
-
- thread_loop_partitioner = v.thread_loop_partitioner;
-
- const size_type this_size = local_size();
- if (this_size > 0)
- {
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::copy(
- thread_loop_partitioner, partitioner->local_size(), v.data, data);
- }
+ Assert(false, ExcNotImplemented());
+ (void)v;
}
: allocated_size(0)
, vector_is_ghosted(false)
{
- reinit(local_range, ghost_indices, communicator);
+ Assert(false, ExcNotImplemented());
+ (void)local_range;
+ (void)ghost_indices;
+ (void)communicator;
}
: allocated_size(0)
, vector_is_ghosted(false)
{
- reinit(local_range, communicator);
+ Assert(false, ExcNotImplemented());
+ (void)local_range;
+ (void)communicator;
}
: allocated_size(0)
, vector_is_ghosted(false)
{
- reinit(size, false);
+ Assert(false, ExcNotImplemented());
+ (void)size;
}
: allocated_size(0)
, vector_is_ghosted(false)
{
- reinit(partitioner);
+ Assert(false, ExcNotImplemented());
+ (void)partitioner;
}
template <typename Number, typename MemorySpaceType>
inline Vector<Number, MemorySpaceType>::~Vector()
{
+ Assert(false, ExcNotImplemented());
try
{
clear_mpi_requests();
Vector<Number, MemorySpaceType>::
operator=(const Vector<Number, MemorySpaceType> &c)
{
-#ifdef _MSC_VER
- return this->operator=<Number>(c);
-#else
- return this->template operator=<Number>(c);
-#endif
+ Assert(false, ExcNotImplemented());
+ (void)c;
+ return *this;
}
Vector<Number, MemorySpaceType>::
operator=(const Vector<Number2, MemorySpaceType> &c)
{
- Assert(c.partitioner.get() != nullptr, ExcNotInitialized());
-
- // we update ghost values whenever one of the input or output vector
- // already held ghost values or when we import data from a vector with
- // the same local range but different ghost layout
- bool must_update_ghost_values = c.vector_is_ghosted;
-
- // check whether the two vectors use the same parallel partitioner. if
- // not, check if all local ranges are the same (that way, we can
- // exchange data between different parallel layouts). One variant which
- // is included here and necessary for compatibility with the other
- // SharedMPI vector classes (Trilinos, PETSc) is the case when vector
- // c does not have any ghosts (constructed without ghost elements given)
- // but the current vector does: In that case, we need to exchange data
- // also when none of the two vector had updated its ghost values before.
- if (partitioner.get() == nullptr)
- reinit(c, true);
- else if (partitioner.get() != c.partitioner.get())
- {
- // local ranges are also the same if both partitioners are empty
- // (even if they happen to define the empty range as [0,0) or [c,c)
- // for some c!=0 in a different way).
- int local_ranges_are_identical =
- (partitioner->local_range() == c.partitioner->local_range() ||
- (partitioner->local_range().second ==
- partitioner->local_range().first &&
- c.partitioner->local_range().second ==
- c.partitioner->local_range().first));
- if ((c.partitioner->n_mpi_processes() > 1 &&
- Utilities::MPI::min(local_ranges_are_identical,
- c.partitioner->get_mpi_communicator()) ==
- 0) ||
- !local_ranges_are_identical)
- reinit(c, true);
- else
- must_update_ghost_values |= vector_is_ghosted;
-
- must_update_ghost_values |=
- (c.partitioner->ghost_indices_initialized() == false &&
- partitioner->ghost_indices_initialized() == true);
- }
- else
- must_update_ghost_values |= vector_is_ghosted;
-
- thread_loop_partitioner = c.thread_loop_partitioner;
-
- const size_type this_size = partitioner->local_size();
- if (this_size > 0)
- {
- dealii::internal::VectorOperations::
- functions<Number, Number2, MemorySpaceType>::copy(
- thread_loop_partitioner, this_size, c.data, data);
- }
-
- if (must_update_ghost_values)
- update_ghost_values();
- else
- zero_out_ghosts();
+ Assert(false, ExcNotImplemented());
+ (void)c;
return *this;
}
Vector<Number, MemorySpaceType>::copy_locally_owned_data_from(
const Vector<Number2, MemorySpaceType> &src)
{
- AssertDimension(partitioner->local_size(), src.partitioner->local_size());
- if (partitioner->local_size() > 0)
- {
- dealii::internal::VectorOperations::
- functions<Number, Number2, MemorySpaceType>::copy(
- thread_loop_partitioner,
- partitioner->local_size(),
- src.data,
- data);
- }
+ Assert(false, ExcNotImplemented());
+ (void)src;
}
const Vector<Number, MemorySpaceType2> &src,
VectorOperation::values operation)
{
- Assert(src.partitioner.get() != nullptr, ExcNotInitialized());
- Assert(partitioner->locally_owned_range() ==
- src.partitioner->locally_owned_range(),
- ExcMessage("Locally owned indices should be identical."));
- Assert(partitioner->ghost_indices() == src.partitioner->ghost_indices(),
- ExcMessage("Ghost indices should be identical."));
- ::dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::import(
- thread_loop_partitioner, allocated_size, operation, src.data, data);
- }
-
-
-
-#ifdef DEAL_II_WITH_PETSC
-
- namespace petsc_helpers
- {
- template <typename PETSC_Number, typename Number>
- void
- copy_petsc_vector(const PETSC_Number *petsc_start_ptr,
- const PETSC_Number *petsc_end_ptr,
- Number * ptr)
- {
- std::copy(petsc_start_ptr, petsc_end_ptr, ptr);
- }
-
- template <typename PETSC_Number, typename Number>
- void
- copy_petsc_vector(const std::complex<PETSC_Number> *petsc_start_ptr,
- const std::complex<PETSC_Number> *petsc_end_ptr,
- std::complex<Number> * ptr)
- {
- std::copy(petsc_start_ptr, petsc_end_ptr, ptr);
- }
-
- template <typename PETSC_Number, typename Number>
- void
- copy_petsc_vector(const std::complex<PETSC_Number> * /*petsc_start_ptr*/,
- const std::complex<PETSC_Number> * /*petsc_end_ptr*/,
- Number * /*ptr*/)
- {
- AssertThrow(false, ExcMessage("Tried to copy complex -> real"));
- }
- } // namespace petsc_helpers
-
- template <typename Number, typename MemorySpaceType>
- Vector<Number, MemorySpaceType> &
- Vector<Number, MemorySpaceType>::
- operator=(const PETScWrappers::MPI::Vector &petsc_vec)
- {
- // TODO: We would like to use the same compact infrastructure as for the
- // Trilinos vector below, but the interface through ReadWriteVector does
- // not support overlapping (ghosted) PETSc vectors, which we need for
- // backward compatibility.
-
- Assert(petsc_vec.locally_owned_elements() == locally_owned_elements(),
- StandardExceptions::ExcInvalidState());
-
- // get a representation of the vector and copy it
- PetscScalar * start_ptr;
- PetscErrorCode ierr =
- VecGetArray(static_cast<const Vec &>(petsc_vec), &start_ptr);
- AssertThrow(ierr == 0, ExcPETScError(ierr));
-
- const size_type vec_size = local_size();
- petsc_helpers::copy_petsc_vector(start_ptr,
- start_ptr + vec_size,
- begin());
-
- // restore the representation of the vector
- ierr = VecRestoreArray(static_cast<const Vec &>(petsc_vec), &start_ptr);
- AssertThrow(ierr == 0, ExcPETScError(ierr));
-
- // spread ghost values between processes?
- if (vector_is_ghosted || petsc_vec.has_ghost_elements())
- update_ghost_values();
-
- // return a reference to this object per normal c++ operator overloading
- // semantics
- return *this;
- }
-
-#endif
-
-
-
-#ifdef DEAL_II_WITH_TRILINOS
-
- template <typename Number, typename MemorySpaceType>
- Vector<Number, MemorySpaceType> &
- Vector<Number, MemorySpaceType>::
- operator=(const TrilinosWrappers::MPI::Vector &trilinos_vec)
- {
-# ifdef DEAL_II_WITH_MPI
- IndexSet combined_set = partitioner->locally_owned_range();
- combined_set.add_indices(partitioner->ghost_indices());
- ReadWriteVector<Number> rw_vector(combined_set);
- rw_vector.import(trilinos_vec, VectorOperation::insert);
- import(rw_vector, VectorOperation::insert);
-
- if (vector_is_ghosted || trilinos_vec.has_ghost_elements())
- update_ghost_values();
-# else
- AssertThrow(false, ExcNotImplemented());
-# endif
-
- return *this;
+ Assert(false, ExcNotImplemented());
+ (void)src;
+ (void)operation;
}
-#endif
-
template <typename Number, typename MemorySpaceType>
Vector<Number, MemorySpaceType>::compress(
::dealii::VectorOperation::values operation)
{
+ Assert(false, ExcNotImplemented());
compress_start(0, operation);
compress_finish(operation);
}
void
Vector<Number, MemorySpaceType>::update_ghost_values() const
{
+ Assert(false, ExcNotImplemented());
update_ghost_values_start();
update_ghost_values_finish();
}
void
Vector<Number, MemorySpaceType>::zero_out_ghosts() const
{
- if (data.values != nullptr)
- std::fill_n(data.values.get() + partitioner->local_size(),
- partitioner->n_ghost_indices(),
- Number());
-#ifdef DEAL_II_COMPILER_CUDA_AWARE
- if (data.values_dev != nullptr)
- {
- const cudaError_t cuda_error_code =
- cudaMemset(data.values_dev.get() + partitioner->local_size(),
- 0,
- partitioner->n_ghost_indices() * sizeof(Number));
- AssertCuda(cuda_error_code);
- }
-#endif
-
- vector_is_ghosted = false;
+ Assert(false, ExcNotImplemented());
}
const unsigned int communication_channel,
::dealii::VectorOperation::values operation)
{
- AssertIndexRange(communication_channel, 200);
- Assert(vector_is_ghosted == false,
- ExcMessage("Cannot call compress() on a ghosted vector"));
-
-#ifdef DEAL_II_WITH_MPI
- // make this function thread safe
- std::lock_guard<std::mutex> lock(mutex);
-
- // allocate import_data in case it is not set up yet
- if (partitioner->n_import_indices() > 0)
- {
-# if defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- if (std::is_same<MemorySpaceType, dealii::MemorySpace::CUDA>::value)
- {
- if (import_data.values_dev == nullptr)
- import_data.values_dev.reset(
- Utilities::CUDA::allocate_device_data<Number>(
- partitioner->n_import_indices()));
- }
- else
-# endif
- {
-# if !defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- static_assert(
- std::is_same<MemorySpaceType, dealii::MemorySpace::Host>::value,
- "This code path should only be compiled for CUDA-aware-MPI for MemorySpace::Host!");
-# endif
- if (import_data.values == nullptr)
- {
- Number *new_val;
- Utilities::System::posix_memalign(
- reinterpret_cast<void **>(&new_val),
- 64,
- sizeof(Number) * partitioner->n_import_indices());
- import_data.values.reset(new_val);
- }
- }
- }
-
-# if defined DEAL_II_COMPILER_CUDA_AWARE && \
- !defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- if (std::is_same<MemorySpaceType, dealii::MemorySpace::CUDA>::value)
- {
- // Move the data to the host and then move it back to the
- // device. We use values to store the elements because the function
- // uses a view of the array and thus we need the data on the host to
- // outlive the scope of the function.
- Number *new_val;
- Utilities::System::posix_memalign(reinterpret_cast<void **>(&new_val),
- 64,
- sizeof(Number) * allocated_size);
-
- data.values.reset(new_val);
-
- cudaError_t cuda_error_code =
- cudaMemcpy(data.values.get(),
- data.values_dev.get(),
- allocated_size * sizeof(Number),
- cudaMemcpyDeviceToHost);
- AssertCuda(cuda_error_code);
- }
-# endif
-
-# if defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- if (std::is_same<MemorySpaceType, dealii::MemorySpace::CUDA>::value)
- {
- partitioner->import_from_ghosted_array_start(
- operation,
- communication_channel,
- ArrayView<Number, MemorySpace::CUDA>(
- data.values_dev.get() + partitioner->local_size(),
- partitioner->n_ghost_indices()),
- ArrayView<Number, MemorySpace::CUDA>(
- import_data.values_dev.get(), partitioner->n_import_indices()),
- compress_requests);
- }
- else
-# endif
- {
- partitioner->import_from_ghosted_array_start(
- operation,
- communication_channel,
- ArrayView<Number, MemorySpace::Host>(
- data.values.get() + partitioner->local_size(),
- partitioner->n_ghost_indices()),
- ArrayView<Number, MemorySpace::Host>(
- import_data.values.get(), partitioner->n_import_indices()),
- compress_requests);
- }
-#else
+ Assert(false, ExcNotImplemented());
(void)communication_channel;
(void)operation;
-#endif
}
Vector<Number, MemorySpaceType>::compress_finish(
::dealii::VectorOperation::values operation)
{
-#ifdef DEAL_II_WITH_MPI
- vector_is_ghosted = false;
-
- // in order to zero ghost part of the vector, we need to call
- // import_from_ghosted_array_finish() regardless of
- // compress_requests.size() == 0
-
- // make this function thread safe
- std::lock_guard<std::mutex> lock(mutex);
-# if defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- if (std::is_same<MemorySpaceType, MemorySpace::CUDA>::value)
- {
- Assert(partitioner->n_import_indices() == 0 ||
- import_data.values_dev != nullptr,
- ExcNotInitialized());
- partitioner
- ->import_from_ghosted_array_finish<Number, MemorySpace::CUDA>(
- operation,
- ArrayView<const Number, MemorySpace::CUDA>(
- import_data.values_dev.get(), partitioner->n_import_indices()),
- ArrayView<Number, MemorySpace::CUDA>(data.values_dev.get(),
- partitioner->local_size()),
- ArrayView<Number, MemorySpace::CUDA>(
- data.values_dev.get() + partitioner->local_size(),
- partitioner->n_ghost_indices()),
- compress_requests);
- }
- else
-# endif
- {
- Assert(partitioner->n_import_indices() == 0 ||
- import_data.values != nullptr,
- ExcNotInitialized());
- partitioner
- ->import_from_ghosted_array_finish<Number, MemorySpace::Host>(
- operation,
- ArrayView<const Number, MemorySpace::Host>(
- import_data.values.get(), partitioner->n_import_indices()),
- ArrayView<Number, MemorySpace::Host>(data.values.get(),
- partitioner->local_size()),
- ArrayView<Number, MemorySpace::Host>(
- data.values.get() + partitioner->local_size(),
- partitioner->n_ghost_indices()),
- compress_requests);
- }
-
-# if defined DEAL_II_COMPILER_CUDA_AWARE && \
- !defined DEAL_II_MPI_WITH_CUDA_SUPPORT
- // The communication is done on the host, so we need to
- // move the data back to the device.
- if (std::is_same<MemorySpaceType, MemorySpace::CUDA>::value)
- {
- cudaError_t cuda_error_code =
- cudaMemcpy(data.values_dev.get(),
- data.values.get(),
- allocated_size * sizeof(Number),
- cudaMemcpyHostToDevice);
- AssertCuda(cuda_error_code);
-
- data.values.reset();
- }
-# endif
-#else
+ Assert(false, ExcNotImplemented());
(void)operation;
-#endif
}
Vector<Number, MemorySpaceType>::update_ghost_values_start(
const unsigned int communication_channel) const
{
- AssertIndexRange(communication_channel, 200);
-#ifdef DEAL_II_WITH_MPI
- // nothing to do when we neither have import nor ghost indices.
- if (partitioner->n_ghost_indices() == 0 &&
- partitioner->n_import_indices() == 0)
- return;
-
- // make this function thread safe
- std::lock_guard<std::mutex> lock(mutex);
-
- // allocate import_data in case it is not set up yet
- if (partitioner->n_import_indices() > 0)
- {
-# if defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- Assert(
- (std::is_same<MemorySpaceType, dealii::MemorySpace::CUDA>::value),
- ExcMessage(
- "Using MemorySpace::CUDA only allowed if the code is compiled with a CUDA compiler!"));
- if (import_data.values_dev == nullptr)
- import_data.values_dev.reset(
- Utilities::CUDA::allocate_device_data<Number>(
- partitioner->n_import_indices()));
-# else
-# ifdef DEAL_II_MPI_WITH_CUDA_SUPPORT
- static_assert(
- std::is_same<MemorySpaceType, dealii::MemorySpace::Host>::value,
- "This code path should only be compiled for CUDA-aware-MPI for MemorySpace::Host!");
-# endif
- if (import_data.values == nullptr)
- {
- Number *new_val;
- Utilities::System::posix_memalign(
- reinterpret_cast<void **>(&new_val),
- 64,
- sizeof(Number) * partitioner->n_import_indices());
- import_data.values.reset(new_val);
- }
-# endif
- }
-
-# if defined DEAL_II_COMPILER_CUDA_AWARE && \
- !defined(DEAL_II_MPI_WITH_CUDA_SUPPORT)
- // Move the data to the host and then move it back to the
- // device. We use values to store the elements because the function
- // uses a view of the array and thus we need the data on the host to
- // outlive the scope of the function.
- Number *new_val;
- Utilities::System::posix_memalign(reinterpret_cast<void **>(&new_val),
- 64,
- sizeof(Number) * allocated_size);
-
- data.values.reset(new_val);
-
- cudaError_t cuda_error_code = cudaMemcpy(data.values.get(),
- data.values_dev.get(),
- allocated_size * sizeof(Number),
- cudaMemcpyDeviceToHost);
- AssertCuda(cuda_error_code);
-# endif
-
-# if !(defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT))
- partitioner->export_to_ghosted_array_start<Number, MemorySpace::Host>(
- communication_channel,
- ArrayView<const Number, MemorySpace::Host>(data.values.get(),
- partitioner->local_size()),
- ArrayView<Number, MemorySpace::Host>(import_data.values.get(),
- partitioner->n_import_indices()),
- ArrayView<Number, MemorySpace::Host>(data.values.get() +
- partitioner->local_size(),
- partitioner->n_ghost_indices()),
- update_ghost_values_requests);
-# else
- partitioner->export_to_ghosted_array_start<Number, MemorySpace::CUDA>(
- communication_channel,
- ArrayView<const Number, MemorySpace::CUDA>(data.values_dev.get(),
- partitioner->local_size()),
- ArrayView<Number, MemorySpace::CUDA>(import_data.values_dev.get(),
- partitioner->n_import_indices()),
- ArrayView<Number, MemorySpace::CUDA>(data.values_dev.get() +
- partitioner->local_size(),
- partitioner->n_ghost_indices()),
- update_ghost_values_requests);
-# endif
-
-#else
+ Assert(false, ExcNotImplemented());
(void)communication_channel;
-#endif
}
void
Vector<Number, MemorySpaceType>::update_ghost_values_finish() const
{
-#ifdef DEAL_II_WITH_MPI
- // wait for both sends and receives to complete, even though only
- // receives are really necessary. this gives (much) better performance
- AssertDimension(partitioner->ghost_targets().size() +
- partitioner->import_targets().size(),
- update_ghost_values_requests.size());
- if (update_ghost_values_requests.size() > 0)
- {
- // make this function thread safe
- std::lock_guard<std::mutex> lock(mutex);
-
-# if !(defined(DEAL_II_COMPILER_CUDA_AWARE) && \
- defined(DEAL_II_MPI_WITH_CUDA_SUPPORT))
- partitioner->export_to_ghosted_array_finish(
- ArrayView<Number, MemorySpace::Host>(
- data.values.get() + partitioner->local_size(),
- partitioner->n_ghost_indices()),
- update_ghost_values_requests);
-# else
- partitioner->export_to_ghosted_array_finish(
- ArrayView<Number, MemorySpace::CUDA>(
- data.values_dev.get() + partitioner->local_size(),
- partitioner->n_ghost_indices()),
- update_ghost_values_requests);
-# endif
- }
-
-# if defined DEAL_II_COMPILER_CUDA_AWARE && \
- !defined DEAL_II_MPI_WITH_CUDA_SUPPORT
- // The communication is done on the host, so we need to
- // move the data back to the device.
- if (std::is_same<MemorySpaceType, MemorySpace::CUDA>::value)
- {
- cudaError_t cuda_error_code =
- cudaMemcpy(data.values_dev.get() + partitioner->local_size(),
- data.values.get() + partitioner->local_size(),
- partitioner->n_ghost_indices() * sizeof(Number),
- cudaMemcpyHostToDevice);
- AssertCuda(cuda_error_code);
-
- data.values.reset();
- }
-# endif
-
-#endif
- vector_is_ghosted = true;
+ Assert(false, ExcNotImplemented());
}
VectorOperation::values operation,
std::shared_ptr<const CommunicationPatternBase> communication_pattern)
{
- // If no communication pattern is given, create one. Otherwise, use the
- // given one.
- std::shared_ptr<const Utilities::MPI::Partitioner> comm_pattern;
- if (communication_pattern.get() == nullptr)
- {
- // Split the IndexSet of V in locally owned elements and ghost indices
- // then create the communication pattern
- IndexSet locally_owned_elem = locally_owned_elements();
- IndexSet ghost_indices = V.get_stored_elements();
- ghost_indices.subtract_set(locally_owned_elem);
- comm_pattern = std::make_shared<Utilities::MPI::Partitioner>(
- locally_owned_elem, ghost_indices, get_mpi_communicator());
- }
- else
- {
- comm_pattern =
- std::dynamic_pointer_cast<const Utilities::MPI::Partitioner>(
- communication_pattern);
- AssertThrow(comm_pattern != nullptr,
- ExcMessage("The communication pattern is not of type "
- "Utilities::MPI::Partitioner."));
- }
- Vector<Number, ::dealii::MemorySpace::Host> tmp_vector(comm_pattern);
-
- data.copy_to(tmp_vector.begin(), local_size());
-
- // fill entries from ReadWriteVector into the SharedMPI vector,
- // including ghost entries. this is not really efficient right now
- // because indices are translated twice, once by nth_index_in_set(i) and
- // once for operator() of tmp_vector
- const IndexSet &v_stored = V.get_stored_elements();
- const size_type v_n_elements = v_stored.n_elements();
- switch (operation)
- {
- case VectorOperation::insert:
- {
- for (size_type i = 0; i < v_n_elements; ++i)
- tmp_vector(v_stored.nth_index_in_set(i)) = V.local_element(i);
-
- break;
- }
- case VectorOperation::add:
- {
- for (size_type i = 0; i < v_n_elements; ++i)
- tmp_vector(v_stored.nth_index_in_set(i)) += V.local_element(i);
-
- break;
- }
- case VectorOperation::min:
- {
- for (size_type i = 0; i < v_n_elements; ++i)
- tmp_vector(v_stored.nth_index_in_set(i)) =
- internal::get_min(tmp_vector(v_stored.nth_index_in_set(i)),
- V.local_element(i));
-
- break;
- }
- case VectorOperation::max:
- {
- for (size_type i = 0; i < v_n_elements; ++i)
- tmp_vector(v_stored.nth_index_in_set(i)) =
- internal::get_max(tmp_vector(v_stored.nth_index_in_set(i)),
- V.local_element(i));
-
- break;
- }
- default:
- {
- Assert(false, ExcMessage("This operation is not supported."));
- }
- }
- tmp_vector.compress(operation);
-
- data.copy_from(tmp_vector.begin(), local_size());
+ Assert(false, ExcNotImplemented());
+ (void)V;
+ (void)operation;
+ (void)communication_pattern;
}
template <typename Number, typename MemorySpaceType>
void
Vector<Number, MemorySpaceType>::swap(Vector<Number, MemorySpaceType> &v)
{
-#ifdef DEAL_II_WITH_MPI
-
-# ifdef DEBUG
- if (Utilities::MPI::job_supports_mpi())
- {
- // make sure that there are not outstanding requests from updating
- // ghost values or compress
- int flag = 1;
- if (update_ghost_values_requests.size() > 0)
- {
- const int ierr = MPI_Testall(update_ghost_values_requests.size(),
- update_ghost_values_requests.data(),
- &flag,
- MPI_STATUSES_IGNORE);
- AssertThrowMPI(ierr);
- Assert(flag == 1,
- ExcMessage(
- "MPI found unfinished update_ghost_values() requests "
- "when calling swap, which is not allowed."));
- }
- if (compress_requests.size() > 0)
- {
- const int ierr = MPI_Testall(compress_requests.size(),
- compress_requests.data(),
- &flag,
- MPI_STATUSES_IGNORE);
- AssertThrowMPI(ierr);
- Assert(flag == 1,
- ExcMessage("MPI found unfinished compress() requests "
- "when calling swap, which is not allowed."));
- }
- }
-# endif
-
- std::swap(compress_requests, v.compress_requests);
- std::swap(update_ghost_values_requests, v.update_ghost_values_requests);
-#endif
-
- std::swap(partitioner, v.partitioner);
- std::swap(thread_loop_partitioner, v.thread_loop_partitioner);
- std::swap(allocated_size, v.allocated_size);
- std::swap(data, v.data);
- std::swap(import_data, v.import_data);
- std::swap(vector_is_ghosted, v.vector_is_ghosted);
+ Assert(false, ExcNotImplemented());
+ (void)v;
}
Vector<Number, MemorySpaceType> &
Vector<Number, MemorySpaceType>::operator=(const Number s)
{
- const size_type this_size = local_size();
- if (this_size > 0)
- {
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::set(
- thread_loop_partitioner, this_size, s, data);
- }
-
- // if we call Vector::operator=0, we want to zero out all the entries
- // plus ghosts.
- if (s == Number())
- zero_out_ghosts();
+ Assert(false, ExcNotImplemented());
+ (void)s;
return *this;
}
Vector<Number, MemorySpaceType>::reinit(const VectorSpaceVector<Number> &V,
const bool omit_zeroing_entries)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&V) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &down_V = dynamic_cast<const VectorType &>(V);
-
- reinit(down_V, omit_zeroing_entries);
+ Assert(false, ExcNotImplemented());
+ (void)V;
+ (void)omit_zeroing_entries;
}
Vector<Number, MemorySpaceType>::
operator+=(const VectorSpaceVector<Number> &vv)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&vv) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- AssertDimension(local_size(), v.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::add_vector(
- thread_loop_partitioner, partitioner->local_size(), v.data, data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)vv;
return *this;
}
Vector<Number, MemorySpaceType>::
operator-=(const VectorSpaceVector<Number> &vv)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&vv) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- AssertDimension(local_size(), v.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::subtract_vector(
- thread_loop_partitioner, partitioner->local_size(), v.data, data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)vv;
return *this;
}
void
Vector<Number, MemorySpaceType>::add(const Number a)
{
- AssertIsFinite(a);
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::add_factor(
- thread_loop_partitioner, partitioner->local_size(), a, data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)a;
}
const Number a,
const VectorSpaceVector<Number> &vv)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&vv) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- AssertIsFinite(a);
- AssertDimension(local_size(), v.local_size());
-
- // nothing to do if a is zero
- if (a == Number(0.))
- return;
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::add_av(
- thread_loop_partitioner, partitioner->local_size(), a, v.data, data);
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)vv;
}
Vector<Number, MemorySpaceType>::add(const Number a,
const VectorSpaceVector<Number> &vv)
{
- add_local(a, vv);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)vv;
}
const Number b,
const VectorSpaceVector<Number> &ww)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&vv) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
- Assert(dynamic_cast<const VectorType *>(&ww) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &w = dynamic_cast<const VectorType &>(ww);
-
- AssertIsFinite(a);
- AssertIsFinite(b);
-
- AssertDimension(local_size(), v.local_size());
- AssertDimension(local_size(), w.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::add_avpbw(
- thread_loop_partitioner,
- partitioner->local_size(),
- a,
- b,
- v.data,
- w.data,
- data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)vv;
+ (void)b;
+ (void)ww;
}
Vector<Number, MemorySpaceType>::add(const std::vector<size_type> &indices,
const std::vector<Number> & values)
{
- for (std::size_t i = 0; i < indices.size(); ++i)
- {
- this->operator()(indices[i]) += values[i];
- }
+ Assert(false, ExcNotImplemented());
+ (void)indices;
+ (void)values;
}
const Number x,
const Vector<Number, MemorySpaceType> &v)
{
- AssertIsFinite(x);
- AssertDimension(local_size(), v.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::sadd_xv(
- thread_loop_partitioner, partitioner->local_size(), x, v.data, data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)x;
+ (void)v;
}
const Number a,
const VectorSpaceVector<Number> &vv)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert((dynamic_cast<const VectorType *>(&vv) != nullptr),
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- AssertIsFinite(x);
- AssertIsFinite(a);
- AssertDimension(local_size(), v.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::sadd_xav(
- thread_loop_partitioner,
- partitioner->local_size(),
- x,
- a,
- v.data,
- data);
+ Assert(false, ExcNotImplemented());
+ (void)x;
+ (void)a;
+ (void)vv;
}
const Number a,
const VectorSpaceVector<Number> &vv)
{
- sadd_local(x, a, vv);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)x;
+ (void)a;
+ (void)vv;
}
const Number b,
const Vector<Number, MemorySpaceType> &w)
{
- AssertIsFinite(x);
- AssertIsFinite(a);
- AssertIsFinite(b);
-
- AssertDimension(local_size(), v.local_size());
- AssertDimension(local_size(), w.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::sadd_xavbw(
- thread_loop_partitioner,
- partitioner->local_size(),
- x,
- a,
- b,
- v.data,
- w.data,
- data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)x;
+ (void)a;
+ (void)v;
+ (void)b;
+ (void)w;
}
Vector<Number, MemorySpaceType> &
Vector<Number, MemorySpaceType>::operator*=(const Number factor)
{
- AssertIsFinite(factor);
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::multiply_factor(
- thread_loop_partitioner, partitioner->local_size(), factor, data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)factor;
return *this;
}
Vector<Number, MemorySpaceType> &
Vector<Number, MemorySpaceType>::operator/=(const Number factor)
{
- operator*=(static_cast<Number>(1.) / factor);
+ Assert(false, ExcNotImplemented());
+ (void)factor;
+
return *this;
}
void
Vector<Number, MemorySpaceType>::scale(const VectorSpaceVector<Number> &vv)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&vv) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- AssertDimension(local_size(), v.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::scale(
- thread_loop_partitioner, local_size(), v.data, data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)vv;
}
Vector<Number, MemorySpaceType>::equ(const Number a,
const VectorSpaceVector<Number> &vv)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert(dynamic_cast<const VectorType *>(&vv) != nullptr,
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- AssertIsFinite(a);
- AssertDimension(local_size(), v.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::equ_au(
- thread_loop_partitioner, partitioner->local_size(), a, v.data, data);
-
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)vv;
}
const Number b,
const Vector<Number, MemorySpaceType> &w)
{
- AssertIsFinite(a);
- AssertIsFinite(b);
-
- AssertDimension(local_size(), v.local_size());
- AssertDimension(local_size(), w.local_size());
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::equ_aubv(
- thread_loop_partitioner,
- partitioner->local_size(),
- a,
- b,
- v.data,
- w.data,
- data);
-
- if (vector_is_ghosted)
- update_ghost_values();
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)v;
+ (void)b;
+ (void)w;
}
bool
Vector<Number, MemorySpaceType>::all_zero() const
{
- return (linfty_norm() == 0) ? true : false;
+ Assert(false, ExcNotImplemented());
+ return false;
}
Vector<Number, MemorySpaceType>::inner_product_local(
const Vector<Number2, MemorySpaceType> &v) const
{
- if (PointerComparison::equal(this, &v))
- return norm_sqr_local();
-
- AssertDimension(partitioner->local_size(), v.partitioner->local_size());
-
- return dealii::internal::VectorOperations::
- functions<Number, Number2, MemorySpaceType>::dot(
- thread_loop_partitioner, partitioner->local_size(), v.data, data);
+ Assert(false, ExcNotImplemented());
+ (void)v;
}
Number Vector<Number, MemorySpaceType>::
operator*(const VectorSpaceVector<Number> &vv) const
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert((dynamic_cast<const VectorType *>(&vv) != nullptr),
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
-
- Number local_result = inner_product_local(v);
- if (partitioner->n_mpi_processes() > 1)
- return Utilities::MPI::sum(local_result,
- partitioner->get_mpi_communicator());
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ (void)vv;
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::norm_sqr_local() const
{
- real_type sum;
-
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::norm_2(
- thread_loop_partitioner, partitioner->local_size(), sum, data);
-
- AssertIsFinite(sum);
-
- return sum;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
Number
Vector<Number, MemorySpaceType>::mean_value_local() const
{
- Assert(size() != 0, ExcEmptyObject());
-
- if (partitioner->local_size() == 0)
- return Number();
-
- Number sum = ::dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::mean_value(
- thread_loop_partitioner, partitioner->local_size(), data);
-
- return sum / real_type(partitioner->local_size());
+ Assert(false, ExcNotImplemented());
+ return 0;
}
Number
Vector<Number, MemorySpaceType>::mean_value() const
{
- Number local_result = mean_value_local();
- if (partitioner->n_mpi_processes() > 1)
- return Utilities::MPI::sum(local_result * static_cast<real_type>(
- partitioner->local_size()),
- partitioner->get_mpi_communicator()) /
- static_cast<real_type>(partitioner->size());
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::l1_norm_local() const
{
- real_type sum;
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::norm_1(
- thread_loop_partitioner, partitioner->local_size(), sum, data);
-
- return sum;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::l1_norm() const
{
- real_type local_result = l1_norm_local();
- if (partitioner->n_mpi_processes() > 1)
- return Utilities::MPI::sum(local_result,
- partitioner->get_mpi_communicator());
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::norm_sqr() const
{
- real_type local_result = norm_sqr_local();
- if (partitioner->n_mpi_processes() > 1)
- return Utilities::MPI::sum(local_result,
- partitioner->get_mpi_communicator());
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::l2_norm() const
{
- return std::sqrt(norm_sqr());
+ Assert(false, ExcNotImplemented());
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::lp_norm_local(const real_type p) const
{
- real_type sum = 0.;
-
- dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::norm_p(
- thread_loop_partitioner, partitioner->local_size(), sum, p, data);
-
- return std::pow(sum, 1. / p);
+ Assert(false, ExcNotImplemented());
+ (void)p;
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::lp_norm(const real_type p) const
{
- const real_type local_result = lp_norm_local(p);
- if (partitioner->n_mpi_processes() > 1)
- return std::pow(
- Utilities::MPI::sum(std::pow(local_result, p),
- partitioner->get_mpi_communicator()),
- static_cast<real_type>(1.0 / p));
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ (void)p;
+ return 0;
}
typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::linfty_norm_local() const
{
- real_type max = 0.;
-
- const size_type local_size = partitioner->local_size();
- internal::la_parallel_vector_templates_functions<
- Number,
- MemorySpaceType>::linfty_norm_local(data, local_size, max);
-
- return max;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
inline typename Vector<Number, MemorySpaceType>::real_type
Vector<Number, MemorySpaceType>::linfty_norm() const
{
- const real_type local_result = linfty_norm_local();
- if (partitioner->n_mpi_processes() > 1)
- return Utilities::MPI::max(local_result,
- partitioner->get_mpi_communicator());
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
const Vector<Number, MemorySpaceType> &v,
const Vector<Number, MemorySpaceType> &w)
{
- const size_type vec_size = partitioner->local_size();
- AssertDimension(vec_size, v.local_size());
- AssertDimension(vec_size, w.local_size());
-
- Number sum = dealii::internal::VectorOperations::
- functions<Number, Number, MemorySpaceType>::add_and_dot(
- thread_loop_partitioner, vec_size, a, v.data, w.data, data);
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)v;
+ (void)w;
- AssertIsFinite(sum);
-
- return sum;
+ return 0;
}
const VectorSpaceVector<Number> &vv,
const VectorSpaceVector<Number> &ww)
{
- // Downcast. Throws an exception if invalid.
- using VectorType = Vector<Number, MemorySpaceType>;
- Assert((dynamic_cast<const VectorType *>(&vv) != nullptr),
- ExcVectorTypeNotCompatible());
- const VectorType &v = dynamic_cast<const VectorType &>(vv);
- Assert((dynamic_cast<const VectorType *>(&ww) != nullptr),
- ExcVectorTypeNotCompatible());
- const VectorType &w = dynamic_cast<const VectorType &>(ww);
-
- Number local_result = add_and_dot_local(a, v, w);
- if (partitioner->n_mpi_processes() > 1)
- return Utilities::MPI::sum(local_result,
- partitioner->get_mpi_communicator());
- else
- return local_result;
+ Assert(false, ExcNotImplemented());
+ (void)a;
+ (void)vv;
+ (void)ww;
+ return 0;
}
Vector<Number, MemorySpaceType>::partitioners_are_compatible(
const Utilities::MPI::Partitioner &part) const
{
+ Assert(false, ExcNotImplemented());
return partitioner->is_compatible(part);
}
Vector<Number, MemorySpaceType>::partitioners_are_globally_compatible(
const Utilities::MPI::Partitioner &part) const
{
+ Assert(false, ExcNotImplemented());
return partitioner->is_globally_compatible(part);
}
std::size_t
Vector<Number, MemorySpaceType>::memory_consumption() const
{
- std::size_t memory = sizeof(*this);
- memory += sizeof(Number) * static_cast<std::size_t>(allocated_size);
-
- // if the partitioner is shared between more processors, just count a
- // fraction of that memory, since we're not actually using more memory
- // for it.
- if (partitioner.use_count() > 0)
- memory +=
- partitioner->memory_consumption() / partitioner.use_count() + 1;
- if (import_data.values != nullptr || import_data.values_dev != nullptr)
- memory += (static_cast<std::size_t>(partitioner->n_import_indices()) *
- sizeof(Number));
- return memory;
+ Assert(false, ExcNotImplemented());
+ return 0;
}
const bool scientific,
const bool across) const
{
- Assert(partitioner.get() != nullptr, ExcInternalError());
- AssertThrow(out, ExcIO());
- std::ios::fmtflags old_flags = out.flags();
- unsigned int old_precision = out.precision(precision);
-
- out.precision(precision);
- if (scientific)
- out.setf(std::ios::scientific, std::ios::floatfield);
- else
- out.setf(std::ios::fixed, std::ios::floatfield);
-
- // to make the vector write out all the information in order, use as
- // many barriers as there are processors and start writing when it's our
- // turn
-#ifdef DEAL_II_WITH_MPI
- if (partitioner->n_mpi_processes() > 1)
- for (unsigned int i = 0; i < partitioner->this_mpi_process(); i++)
- {
- const int ierr = MPI_Barrier(partitioner->get_mpi_communicator());
- AssertThrowMPI(ierr);
- }
-#endif
-
- std::vector<Number> stored_elements(allocated_size);
- data.copy_to(stored_elements.data(), allocated_size);
-
- out << "Process #" << partitioner->this_mpi_process() << std::endl
- << "Local range: [" << partitioner->local_range().first << ", "
- << partitioner->local_range().second
- << "), global size: " << partitioner->size() << std::endl
- << "Vector data:" << std::endl;
- if (across)
- for (size_type i = 0; i < partitioner->local_size(); ++i)
- out << stored_elements[i] << ' ';
- else
- for (size_type i = 0; i < partitioner->local_size(); ++i)
- out << stored_elements[i] << std::endl;
- out << std::endl;
-
- if (vector_is_ghosted)
- {
- out << "Ghost entries (global index / value):" << std::endl;
- if (across)
- for (size_type i = 0; i < partitioner->n_ghost_indices(); ++i)
- out << '(' << partitioner->ghost_indices().nth_index_in_set(i)
- << '/' << stored_elements[partitioner->local_size() + i]
- << ") ";
- else
- for (size_type i = 0; i < partitioner->n_ghost_indices(); ++i)
- out << '(' << partitioner->ghost_indices().nth_index_in_set(i)
- << '/' << stored_elements[partitioner->local_size() + i]
- << ")" << std::endl;
- out << std::endl;
- }
- out << std::flush;
-
-#ifdef DEAL_II_WITH_MPI
- if (partitioner->n_mpi_processes() > 1)
- {
- int ierr = MPI_Barrier(partitioner->get_mpi_communicator());
- AssertThrowMPI(ierr);
-
- for (unsigned int i = partitioner->this_mpi_process() + 1;
- i < partitioner->n_mpi_processes();
- i++)
- {
- ierr = MPI_Barrier(partitioner->get_mpi_communicator());
- AssertThrowMPI(ierr);
- }
- }
-#endif
-
- AssertThrow(out, ExcIO());
- // reset output format
- out.flags(old_flags);
- out.precision(old_precision);
+ Assert(false, ExcNotImplemented());
+ (void)out;
+ (void)precision;
+ (void)scientific;
+ (void)across;
}
} // end of namespace SharedMPI
projects / dealii.git / commitdiff