// and we cannot use the sadd function directly, so we have to create
// a completely distributed vector first and copy the local entries
// from the vector with ghost entries
- TrilinosWrappers::MPI::Vector u1_completely_distributed (u1_difference.vector_partitioner ());
-
+ TrilinosWrappers::MPI::Vector u1_completely_distributed;
+ u1_completely_distributed.reinit(u1_difference, true);
u1_completely_distributed = u1;
u1_difference.sadd(-1, u1_completely_distributed);
*/
BlockVector ();
- /**
- * Constructor. Generate a block vector with as many blocks as there are
- * entries in Input_Maps. For this non-distributed vector, the %parallel
- * partitioning is not used, just the global size of the partitioner.
- */
- explicit BlockVector (const std::vector<Epetra_Map> &partitioner) DEAL_II_DEPRECATED;
-
- /**
- * Constructor. Generate a block vector with as many blocks as there are
- * entries in Input_Maps. For this non-distributed vector, the %parallel
- * partitioning is not used, just the global size of the partitioner.
- */
- explicit BlockVector (const std::vector<IndexSet> &partitioner,
- const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
-
- /**
- * Copy-Constructor. Set all the properties of the non-%parallel vector to
- * those of the given %parallel vector and import the elements.
- */
- BlockVector (const MPI::BlockVector &V) DEAL_II_DEPRECATED;
-
- /**
- * Copy-Constructor. Set all the properties of the vector to those of the
- * given input vector and copy the elements.
- */
- BlockVector (const BlockVector &V) DEAL_II_DEPRECATED;
-
- /**
- * Creates a block vector consisting of <tt>num_blocks</tt> components,
- * but there is no content in the individual components and the user has
- * to fill appropriate data using a reinit of the blocks.
- */
- explicit BlockVector (const size_type num_blocks) DEAL_II_DEPRECATED;
-
- /**
- * Constructor. Set the number of blocks to <tt>n.size()</tt> and
- * initialize each block with <tt>n[i]</tt> zero elements.
- *
- * References BlockVector.reinit().
- */
- explicit BlockVector (const std::vector<size_type> &N) DEAL_II_DEPRECATED;
-
- /**
- * Constructor. Set the number of blocks to <tt>n.size()</tt>. Initialize
- * the vector with the elements pointed to by the range of iterators given
- * as second and third argument. Apart from the first argument, this
- * constructor is in complete analogy to the respective constructor of the
- * <tt>std::vector</tt> class, but the first argument is needed in order
- * to know how to subdivide the block vector into different blocks.
- */
- template <typename InputIterator>
- BlockVector (const std::vector<size_type> &n,
- const InputIterator first,
- const InputIterator end) DEAL_II_DEPRECATED;
-
/**
* Destructor. Clears memory
*/
inline
- BlockVector::BlockVector ()
- {}
+ BlockVector::BlockVector () = default;
inline
- BlockVector::BlockVector (const std::vector<Epetra_Map> &partitioning)
- {
- reinit (partitioning);
- }
-
-
+ BlockVector::~BlockVector() = default;
- inline
- BlockVector::BlockVector (const std::vector<IndexSet> &partitioning,
- const MPI_Comm &communicator)
- {
- reinit (partitioning, communicator);
- }
-
-
-
- inline
- BlockVector::BlockVector (const std::vector<size_type> &N)
- {
- reinit (N);
- }
-
-
-
- template <typename InputIterator>
- BlockVector::BlockVector (const std::vector<size_type> &n,
- const InputIterator first,
- const InputIterator end)
- {
- // first set sizes of blocks, but
- // don't initialize them as we will
- // copy elements soon
- (void)end;
- reinit (n, true);
- InputIterator start = first;
- for (size_type b=0; b<n.size(); ++b)
- {
- InputIterator end = start;
- std::advance (end, static_cast<size_type>(n[b]));
-
- for (size_type i=0; i<n[b]; ++i, ++start)
- this->block(b)(i) = *start;
- }
- Assert (start == end, ExcIteratorRangeDoesNotMatchVectorSize());
- }
-
-
-
- inline
- BlockVector::BlockVector (const size_type num_blocks)
- {
- reinit (num_blocks);
- }
-
-
-
- inline
- BlockVector::~BlockVector()
- {}
-
-
-
- inline
- BlockVector::BlockVector (const MPI::BlockVector &v)
- {
- reinit (v);
- }
-
-
-
- inline
- BlockVector::BlockVector (const BlockVector &v)
- :
- BlockVectorBase<Vector > ()
- {
- this->components.resize (v.n_blocks());
- this->block_indices = v.block_indices;
-
- for (size_type i=0; i<this->n_blocks(); ++i)
- this->components[i] = v.components[i];
- }
inline
*/
BlockVector ();
- /**
- * Constructor. Generate a block vector with as many blocks as there are
- * entries in @p partitioning. Each Epetra_Map contains the layout of
- * the distribution of data among the MPI processes.
- *
- * This function is deprecated.
- */
- explicit BlockVector (const std::vector<Epetra_Map> ¶llel_partitioning) DEAL_II_DEPRECATED;
-
/**
* Constructor. Generate a block vector with as many blocks as there are
* entries in @p partitioning. Each IndexSet together with the MPI
void import_nonlocal_data_for_fe (const TrilinosWrappers::BlockSparseMatrix &m,
const BlockVector &v);
- /**
- * Return the state of the vector, i.e., whether compress() needs to be
- * called after an operation requiring data exchange. Does only return
- * non-true values when used in <tt>debug</tt> mode, since it is quite
- * expensive to keep track of all operations that lead to the need for
- * compress().
- *
- * This function is deprecated.
- */
- bool is_compressed () const DEAL_II_DEPRECATED;
-
/**
* Return if this Vector contains ghost elements.
*
void import_nonlocal_data_for_fe
(const dealii::TrilinosWrappers::SparseMatrix &matrix,
const Vector &vector);
-//@}
- /**
- * @name Initialization with an Epetra_Map
- */
-//@{
- /**
- * This constructor takes an Epetra_Map that already knows how to
- * distribute the individual components among the MPI processors. Since
- * it also includes information about the size of the vector, this is
- * all we need to generate a parallel vector.
- *
- * Depending on whether the @p parallel_partitioning argument uniquely
- * subdivides elements among processors or not, the resulting vector may
- * or may not have ghost elements. See the general documentation of this
- * class for more information.
- *
- * This function is deprecated.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
- explicit Vector (const Epetra_Map ¶llel_partitioning) DEAL_II_DEPRECATED;
-
- /**
- * Copy constructor from the TrilinosWrappers vector class. Since a
- * vector of this class does not necessarily need to be distributed
- * among processes, the user needs to supply us with an Epetra_Map that
- * sets the partitioning details.
- *
- * Depending on whether the @p parallel_partitioning argument uniquely
- * subdivides elements among processors or not, the resulting vector may
- * or may not have ghost elements. See the general documentation of this
- * class for more information.
- *
- * This function is deprecated.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
- Vector (const Epetra_Map ¶llel_partitioning,
- const VectorBase &v) DEAL_II_DEPRECATED;
-
- /**
- * Reinitialize from a deal.II vector. The Epetra_Map specifies the
- * %parallel partitioning.
- *
- * Depending on whether the @p parallel_partitioning argument uniquely
- * subdivides elements among processors or not, the resulting vector may
- * or may not have ghost elements. See the general documentation of this
- * class for more information.
- *
- * This function is deprecated.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
- template <typename number>
- void reinit (const Epetra_Map ¶llel_partitioner,
- const dealii::Vector<number> &v) DEAL_II_DEPRECATED;
-
- /**
- * Reinit functionality. This function destroys the old vector content
- * and generates a new one based on the input map.
- *
- * Depending on whether the @p parallel_partitioning argument uniquely
- * subdivides elements among processors or not, the resulting vector may
- * or may not have ghost elements. See the general documentation of this
- * class for more information.
- *
- * This function is deprecated.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
- void reinit (const Epetra_Map ¶llel_partitioning,
- const bool omit_zeroing_entries = false) DEAL_II_DEPRECATED;
-
- /**
- * Copy-constructor from deal.II vectors. Sets the dimension to that of
- * the given vector, and copies all elements.
- *
- * Depending on whether the @p parallel_partitioning argument uniquely
- * subdivides elements among processors or not, the resulting vector may
- * or may not have ghost elements. See the general documentation of this
- * class for more information.
- *
- * This function is deprecated.
- *
- * @see
- * @ref GlossGhostedVector "vectors with ghost elements"
- */
- template <typename Number>
- Vector (const Epetra_Map ¶llel_partitioning,
- const dealii::Vector<Number> &v) DEAL_II_DEPRECATED;
//@}
/**
* @name Initialization with an IndexSet
#ifndef DOXYGEN
- template <typename number>
- Vector::Vector (const Epetra_Map &input_map,
- const dealii::Vector<number> &v)
- {
- reinit (input_map, v);
- }
-
-
-
template <typename number>
Vector::Vector (const IndexSet ¶llel_partitioner,
const dealii::Vector<number> &v,
-
- template <typename number>
- void Vector::reinit (const Epetra_Map ¶llel_partitioner,
- const dealii::Vector<number> &v)
- {
- if (vector.get() == nullptr || vector->Map().SameAs(parallel_partitioner) == false)
- reinit(parallel_partitioner);
-
- const int size = parallel_partitioner.NumMyElements();
-
- // Need to copy out values, since the deal.II might not use doubles, so
- // that a direct access is not possible.
- for (int i=0; i<size; ++i)
- (*vector)[0][i] = v(gid(parallel_partitioner,i));
-
- }
-
-
inline
Vector &
Vector::operator= (const TrilinosScalar s)
*/
typedef dealii::types::global_dof_index size_type;
- /**
- * Default constructor that generates an empty (zero size) vector. The
- * function <tt>reinit()</tt> will have to give the vector the correct
- * size.
- */
- Vector () DEAL_II_DEPRECATED;
-
- /**
- * This constructor takes as input the number of elements in the vector.
- */
- explicit Vector (const size_type n) DEAL_II_DEPRECATED;
-
- /**
- * This constructor takes as input the number of elements in the vector.
- * If the map is not localized, i.e., if there are some elements that are
- * not present on all processes, only the global size of the map will be
- * taken and a localized map will be generated internally. In other words,
- * which element of the @p partitioning argument are set is in fact
- * ignored, the only thing that matters is the size of the index space
- * described by this argument.
- */
- explicit Vector (const Epetra_Map &partitioning) DEAL_II_DEPRECATED;
-
- /**
- * This constructor takes as input the number of elements in the vector.
- * If the index set is not localized, i.e., if there are some elements
- * that are not present on all processes, only the global size of the
- * index set will be taken and a localized version will be generated
- * internally. In other words, which element of the @p partitioning
- * argument are set is in fact ignored, the only thing that matters is the
- * size of the index space described by this argument.
- */
- explicit Vector (const IndexSet &partitioning,
- const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
-
- /**
- * This constructor takes a (possibly parallel) Trilinos Vector and
- * generates a localized version of the whole content on each processor.
- */
- explicit Vector (const VectorBase &V) DEAL_II_DEPRECATED;
-
- /**
- * Copy-constructor from deal.II vectors. Sets the dimension to that of
- * the given vector, and copies all elements.
- */
- template <typename Number>
- explicit Vector (const dealii::Vector<Number> &v) DEAL_II_DEPRECATED;
-
/**
* Reinit function that resizes the vector to the size specified by
* <tt>n</tt>. The flag <tt>omit_zeroing_entries</tt> specifies whether
#ifndef DOXYGEN
- template <typename number>
- Vector::Vector (const dealii::Vector<number> &v)
- {
- Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
- vector.reset (new Epetra_FEVector(map));
- *this = v;
- }
-
-
-
inline
Vector &
Vector::operator= (const TrilinosScalar s)
*/
void compress (::dealii::VectorOperation::values operation);
- /**
- * Return the state of the vector, i.e., whether compress() has already
- * been called after an operation requiring data exchange.
- *
- * This function is deprecated.
- */
- bool is_compressed () const DEAL_II_DEPRECATED;
-
/**
* Set all components of the vector to the given number @p s. Simply pass
* this down to the Trilinos Epetra object, but we still need to declare
*/
TrilinosScalar mean_value () const;
- /**
- * Compute the minimal value of the elements of this vector.
- *
- * This function is deprecated use min() instead.
- */
- TrilinosScalar minimal_value () const DEAL_II_DEPRECATED;
-
/**
* Compute the minimal value of the elements of this vector.
*/
const ForwardIterator indices_end,
OutputIterator values_begin) const;
- /**
- * Return the value of the vector entry <i>i</i>. Note that this function
- * does only work properly when we request a data stored on the local
- * processor. In case the elements sits on another process, this function
- * returns 0 which might or might not be appropriate in a given situation.
- * If you rely on consistent results, use the access functions () or []
- * that throw an assertion in case a non-local element is used.
- *
- * This function is deprecated.
- */
- TrilinosScalar el (const size_type index) const DEAL_II_DEPRECATED;
-
/**
* Make the 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
const TrilinosScalar a,
const VectorBase &V);
- /**
- * Scaling and multiple addition.
- *
- * This function is deprecated.
- */
- void sadd (const TrilinosScalar s,
- const TrilinosScalar a,
- const VectorBase &V,
- const TrilinosScalar b,
- const VectorBase &W) DEAL_II_DEPRECATED;
-
- /**
- * Scaling and multiple addition. <tt>*this = s*(*this) + a*V + b*W +
- * c*X</tt>.
- *
- * This function is deprecated.
- */
- void sadd (const TrilinosScalar s,
- const TrilinosScalar a,
- const VectorBase &V,
- const TrilinosScalar b,
- const VectorBase &W,
- const TrilinosScalar c,
- const VectorBase &X) DEAL_II_DEPRECATED;
-
/**
* Scale each element of this vector by the corresponding element in the
* argument. This function is mostly meant to simulate multiplication (and
*/
void equ (const TrilinosScalar a,
const VectorBase &V);
-
- /**
- * Assignment <tt>*this = a*V + b*W</tt>.
- *
- * This function is deprecated.
- */
- void equ (const TrilinosScalar a,
- const VectorBase &V,
- const TrilinosScalar b,
- const VectorBase &W) DEAL_II_DEPRECATED;
-
- /**
- * Compute the elementwise ratio of the two given vectors, that is let
- * <tt>this[i] = a[i]/b[i]</tt>. This is useful for example if you want to
- * compute the cellwise ratio of true to estimated error.
- *
- * This vector is appropriately scaled to hold the result.
- *
- * If any of the <tt>b[i]</tt> is zero, the result is undefined. No
- * attempt is made to catch such situations.
- */
- void ratio (const VectorBase &a,
- const VectorBase &b) DEAL_II_DEPRECATED;
//@}
*/
const Epetra_Map &vector_partitioner () const;
- /**
- * Output of vector in user-defined format in analogy to the
- * dealii::Vector class.
- *
- * This function is deprecated.
- */
- void print (const char *format = nullptr) const DEAL_II_DEPRECATED;
-
/**
* Print to a stream. @p precision denotes the desired precision with
* which values shall be printed, @p scientific whether scientific
- inline
- bool
- VectorBase::is_compressed () const
- {
- return compressed;
- }
-
-
-
inline
bool
VectorBase::in_local_range (const size_type index) const
- inline
- TrilinosScalar
- VectorBase::minimal_value () const
- {
- return min();
- }
-
-
-
inline
TrilinosScalar
VectorBase::min () const
- inline
- void
- VectorBase::sadd (const TrilinosScalar s,
- const TrilinosScalar a,
- const VectorBase &v,
- const TrilinosScalar b,
- const VectorBase &w)
- {
- // if we have ghost values, do not allow
- // writing to this vector at all.
- Assert (!has_ghost_elements(), ExcGhostsPresent());
- Assert (size() == v.size(),
- ExcDimensionMismatch (size(), v.size()));
- Assert (size() == w.size(),
- ExcDimensionMismatch (size(), w.size()));
- AssertIsFinite(s);
- AssertIsFinite(a);
- AssertIsFinite(b);
-
- // We assume that the vectors have the same Map
- // if the local size is the same and if the vectors are not ghosted
- if (local_size() == v.local_size() && !v.has_ghost_elements() &&
- local_size() == w.local_size() && !w.has_ghost_elements())
- {
- Assert (this->vector->Map().SameAs(v.vector->Map())==true,
- ExcDifferentParallelPartitioning());
- Assert (this->vector->Map().SameAs(w.vector->Map())==true,
- ExcDifferentParallelPartitioning());
- const int ierr = vector->Update(a, *(v.vector), b, *(w.vector), s);
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
- }
- else
- {
- this->sadd( s, a, v);
- this->sadd(1., b, w);
- }
- }
-
-
-
- inline
- void
- VectorBase::sadd (const TrilinosScalar s,
- const TrilinosScalar a,
- const VectorBase &v,
- const TrilinosScalar b,
- const VectorBase &w,
- const TrilinosScalar c,
- const VectorBase &x)
- {
- // if we have ghost values, do not allow
- // writing to this vector at all.
- Assert (!has_ghost_elements(), ExcGhostsPresent());
- Assert (size() == v.size(),
- ExcDimensionMismatch (size(), v.size()));
- Assert (size() == w.size(),
- ExcDimensionMismatch (size(), w.size()));
- Assert (size() == x.size(),
- ExcDimensionMismatch (size(), x.size()));
- AssertIsFinite(s);
- AssertIsFinite(a);
- AssertIsFinite(b);
- AssertIsFinite(c);
-
- // We assume that the vectors have the same Map
- // if the local size is the same and if the vectors are not ghosted
- if (local_size() == v.local_size() && !v.has_ghost_elements() &&
- local_size() == w.local_size() && !w.has_ghost_elements() &&
- local_size() == x.local_size() && !x.has_ghost_elements())
- {
- Assert (this->vector->Map().SameAs(v.vector->Map())==true,
- ExcDifferentParallelPartitioning());
- Assert (this->vector->Map().SameAs(w.vector->Map())==true,
- ExcDifferentParallelPartitioning());
- Assert (this->vector->Map().SameAs(x.vector->Map())==true,
- ExcDifferentParallelPartitioning());
-
- // Update member can only
- // input two other vectors so
- // do it in two steps
- const int ierr = vector->Update(a, *(v.vector), b, *(w.vector), s);
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-
- const int jerr = vector->Update(c, *(x.vector), 1.);
- Assert (jerr == 0, ExcTrilinosError(jerr));
- (void)jerr; // removes -Wunused-parameter warning in optimized mode
- }
- else
- {
- this->sadd( s, a, v);
- this->sadd(1., b, w);
- this->sadd(1., c, x);
- }
- }
-
-
-
inline
void
VectorBase::scale (const VectorBase &factors)
- inline
- void
- VectorBase::ratio (const VectorBase &v,
- const VectorBase &w)
- {
- Assert (v.local_size() == w.local_size(),
- ExcDimensionMismatch (v.local_size(), w.local_size()));
- Assert (local_size() == w.local_size(),
- ExcDimensionMismatch (local_size(), w.local_size()));
-
- const int ierr = vector->ReciprocalMultiply(1.0, *(w.vector),
- *(v.vector), 0.0);
-
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
- }
-
-
-
inline
const Epetra_MultiVector &
VectorBase::trilinos_vector () const
// be a better solution than
// this, but it has not yet been
// found.
- TrilinosWrappers::Vector localized_vector (v);
+ TrilinosWrappers::Vector localized_vector;
+ localized_vector.reinit (v, false, true);
// get a representation of the vector
// and copy it
// of the Trilinos vectors and
// then call the other =
// operator.
- TrilinosWrappers::Vector localized_vector (v);
+ TrilinosWrappers::Vector localized_vector;
+ localized_vector.reinit(v, false, true);
*this = localized_vector;
return *this;
}
trilinos_block_vector.cc
trilinos_epetra_communication_pattern.cc
trilinos_epetra_vector.cc
- trilinos_parallel_block_vector.cc
trilinos_precondition.cc
trilinos_precondition_ml.cc
trilinos_precondition_muelu.cc
- void
- BlockVector::reinit (const std::vector<Epetra_Map> &input_maps,
- const bool omit_zeroing_entries)
- {
- const size_type no_blocks = input_maps.size();
- std::vector<size_type> block_sizes (no_blocks);
-
- for (size_type i=0; i<no_blocks; ++i)
- {
- block_sizes[i] = n_global_elements(input_maps[i]);
- }
-
- this->block_indices.reinit (block_sizes);
- if (components.size() != n_blocks())
- components.resize(n_blocks());
-
- for (size_type i=0; i<n_blocks(); ++i)
- components[i].reinit(input_maps[i], omit_zeroing_entries);
-
- collect_sizes();
- }
-
-
-
void
BlockVector::reinit (const std::vector<IndexSet> ¶llel_partitioning,
const MPI_Comm &communicator,
+++ /dev/null
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 2015 by the deal.II authors
-//
-// This file is part of the deal.II library.
-//
-// The deal.II library is free software; you can use it, redistribute
-// it, and/or modify it under the terms of the GNU Lesser General
-// Public License as published by the Free Software Foundation; either
-// version 2.1 of the License, or (at your option) any later version.
-// The full text of the license can be found in the file LICENSE at
-// the top level of the deal.II distribution.
-//
-// ---------------------------------------------------------------------
-
-
-#include <deal.II/lac/trilinos_parallel_block_vector.h>
-
-#ifdef DEAL_II_WITH_TRILINOS
-
-DEAL_II_NAMESPACE_OPEN
-
-namespace TrilinosWrappers
-{
- namespace MPI
- {
- BlockVector::BlockVector (const std::vector<Epetra_Map> ¶llel_partitioning)
- {
- reinit (parallel_partitioning, false);
- }
-
-
-
- bool
- BlockVector::is_compressed () const
- {
- bool compressed = true;
- for (unsigned int row=0; row<n_blocks(); ++row)
- if (block(row).is_compressed() == false)
- {
- compressed = false;
- break;
- }
-
- return compressed;
- }
- }
-}
-
-DEAL_II_NAMESPACE_CLOSE
-
-#endif // DEAL_II_WITH_TRILINOS
- Vector::Vector (const Epetra_Map ¶llel_partitioning)
- {
- reinit (parallel_partitioning);
- }
-
-
-
Vector::Vector (const IndexSet ¶llel_partitioning,
const MPI_Comm &communicator)
{
- Vector::Vector (const Epetra_Map &input_map,
- const VectorBase &v)
- :
- VectorBase()
- {
- AssertThrow (n_global_elements(input_map) == n_global_elements(v.vector->Map()),
- ExcDimensionMismatch (n_global_elements(input_map),
- n_global_elements(v.vector->Map())));
-
- last_action = Zero;
-
- if (input_map.SameAs(v.vector->Map()) == true)
- vector.reset (new Epetra_FEVector(*v.vector));
- else
- {
- vector.reset (new Epetra_FEVector(input_map));
- reinit (v, false, true);
- }
- }
-
-
-
Vector::Vector (const IndexSet ¶llel_partitioner,
const VectorBase &v,
const MPI_Comm &communicator)
- Vector::~Vector ()
- {}
-
-
-
- void
- Vector::reinit (const Epetra_Map &input_map,
- const bool /*omit_zeroing_entries*/)
- {
- nonlocal_vector.reset();
-
- vector.reset (new Epetra_FEVector(input_map));
-
- has_ghosts = vector->Map().UniqueGIDs()==false;
-
- // If the IndexSets are overlapping, we don't really know
- // which process owns what. So we decide that no process
- // owns anything in that case. In particular asking for
- // the locally owned elements is not allowed.
- owned_elements.clear();
- if (has_ghosts)
- owned_elements.set_size(0);
- else
- {
- owned_elements.set_size(size());
-
- // easy case: local range is contiguous
- if (vector->Map().LinearMap())
- {
- const std::pair<size_type, size_type> x = local_range();
- owned_elements.add_range (x.first, x.second);
- }
- else if (vector->Map().NumMyElements() > 0)
- {
- const size_type n_indices = vector->Map().NumMyElements();
-#ifndef DEAL_II_WITH_64BIT_INDICES
- unsigned int *vector_indices = (unsigned int *)vector->Map().MyGlobalElements();
-#else
- size_type *vector_indices = (size_type *)vector->Map().MyGlobalElements64();
-#endif
- owned_elements.add_indices(vector_indices, vector_indices+n_indices);
- owned_elements.compress();
- }
- }
-#if defined(DEBUG) && defined(DEAL_II_WITH_MPI)
- const Epetra_MpiComm *comm_ptr
- = dynamic_cast<const Epetra_MpiComm *>(&(vector->Comm()));
- Assert (comm_ptr != nullptr, ExcInternalError());
- const size_type n_elements_global
- = Utilities::MPI::sum (owned_elements.n_elements(), comm_ptr->Comm());
-
- Assert (has_ghosts || n_elements_global == size(), ExcInternalError());
-#endif
-
- last_action = Zero;
- }
+ Vector::~Vector () = default;
Assert (comm_ptr != nullptr, ExcInternalError());
const size_type n_elements_global
= Utilities::MPI::sum (owned_elements.n_elements(), comm_ptr->Comm());
-
Assert (has_ghosts || n_elements_global == size(), ExcInternalError());
#endif
}
- void
- Vector::reinit (const BlockVector &v,
- const bool import_data)
- {
- nonlocal_vector.reset();
- owned_elements.clear();
- owned_elements.set_size(v.size());
-
- // In case we do not allow to have different maps, this call means that
- // we have to reset the vector. So clear the vector, initialize our map
- // with the map in v, and generate the vector.
- if (v.n_blocks() == 0)
- return;
-
- // create a vector that holds all the elements contained in the block
- // vector. need to manually create an Epetra_Map.
- size_type n_elements = 0, added_elements = 0, block_offset = 0;
- for (size_type block=0; block<v.n_blocks(); ++block)
- n_elements += v.block(block).local_size();
- std::vector<TrilinosWrappers::types::int_type> global_ids (n_elements, -1);
- for (size_type block=0; block<v.n_blocks(); ++block)
- {
- TrilinosWrappers::types::int_type *glob_elements =
- my_global_elements(v.block(block).vector_partitioner());
- for (size_type i=0; i<v.block(block).local_size(); ++i)
- global_ids[added_elements++] = glob_elements[i] + block_offset;
- owned_elements.add_indices(v.block(block).owned_elements,
- block_offset);
- block_offset += v.block(block).size();
- }
-
- Assert (n_elements == added_elements, ExcInternalError());
- Epetra_Map new_map (v.size(), n_elements, &global_ids[0], 0,
- v.block(0).vector_partitioner().Comm());
-
- std::shared_ptr<Epetra_FEVector> actual_vec;
- if ( import_data == true )
- actual_vec.reset (new Epetra_FEVector (new_map));
- else
- {
- vector.reset (new Epetra_FEVector (new_map));
- actual_vec = vector;
- }
-
- TrilinosScalar *entries = (*actual_vec)[0];
- block_offset = 0;
- for (size_type block=0; block<v.n_blocks(); ++block)
- {
- v.block(block).trilinos_vector().ExtractCopy (entries, 0);
- entries += v.block(block).local_size();
- }
-
- if (import_data == true)
- {
- AssertThrow (static_cast<size_type>(global_length(*actual_vec))
- == v.size(),
- ExcDimensionMismatch (global_length(*actual_vec),
- v.size()));
-
- Epetra_Import data_exchange (vector->Map(), actual_vec->Map());
-
- const int ierr = vector->Import(*actual_vec, data_exchange, Insert);
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-
- last_action = Insert;
- }
-#if defined(DEBUG) && defined(DEAL_II_WITH_MPI)
- const Epetra_MpiComm *comm_ptr
- = dynamic_cast<const Epetra_MpiComm *>(&(vector->Comm()));
- Assert (comm_ptr != nullptr, ExcInternalError());
- const size_type n_elements_global
- = Utilities::MPI::sum (owned_elements.n_elements(), comm_ptr->Comm());
-
- Assert (has_ghosts || n_elements_global == size(), ExcInternalError());
-#endif
- }
-
-
void Vector::reinit(const IndexSet &locally_owned_entries,
const IndexSet &ghost_entries,
const MPI_Comm &communicator,
Epetra_Map map = parallel_partitioner.make_trilinos_map (communicator,
true);
vector.reset (new Epetra_FEVector(map));
- has_ghosts = vector->Map().UniqueGIDs()==false;
-
- last_action = Zero;
}
else
{
Assert (map.IsOneToOne(),
ExcMessage("A writable vector must not have ghost entries in "
"its parallel partitioning"));
- reinit (map);
+
+ if (vector->Map().SameAs(map)==false)
+ vector.reset (new Epetra_FEVector(map));
+ else
+ {
+ const int ierr = vector->PutScalar(0.);
+ (void)ierr;
+ Assert (ierr == 0, ExcTrilinosError(ierr));
+ }
IndexSet nonlocal_entries(ghost_entries);
nonlocal_entries.subtract_set(locally_owned_entries);
nonlocal_vector.reset(new Epetra_MultiVector(nonlocal_map, 1));
}
}
+
+ has_ghosts = vector->Map().UniqueGIDs()==false;
+
+ last_action = Zero;
+
#ifdef DEBUG
const size_type n_elements_global
= Utilities::MPI::sum (owned_elements.n_elements(), communicator);
-
- Vector::Vector ()
- {
- last_action = Zero;
- Epetra_LocalMap map (0, 0, Utilities::Trilinos::comm_self());
- vector.reset (new Epetra_FEVector(map));
- }
-
-
-
- Vector::Vector (const size_type n)
- {
- reinit(n);
- }
-
-
-
- Vector::Vector (const Epetra_Map &input_map)
- {
- reinit(input_map);
- }
-
-
-
- Vector::Vector (const IndexSet &partitioning,
- const MPI_Comm &communicator)
- {
- reinit (partitioning, communicator);
- }
-
-
-
- Vector::Vector (const VectorBase &v)
- {
- last_action = Zero;
- Epetra_LocalMap map (n_global_elements(v.vector->Map()),
- v.vector->Map().IndexBase(),
- v.vector->Map().Comm());
- vector.reset (new Epetra_FEVector(map));
-
- if (vector->Map().SameAs(v.vector->Map()) == true)
- {
- const int ierr = vector->Update(1.0, *v.vector, 0.0);
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
- }
- else
- reinit (v, false, true);
- }
-
-
-
void
Vector::reinit (const size_type n,
const bool /*omit_zeroing_entries*/)
- TrilinosScalar
- VectorBase::el (const size_type index) const
- {
- // Extract local indices in the vector.
- TrilinosWrappers::types::int_type trilinos_i =
- vector->Map().LID(static_cast<TrilinosWrappers::types::int_type>(index));
-
- // If the element is not present on the current processor, we can't
- // continue. Just print out 0 as opposed to the () method below.
- if (trilinos_i == -1)
- return 0.;
- else
- return (*vector)[0][trilinos_i];
- }
-
-
-
TrilinosScalar
VectorBase::operator () (const size_type index) const
{
- void
- VectorBase::equ (const TrilinosScalar a,
- const VectorBase &v,
- const TrilinosScalar b,
- const VectorBase &w)
- {
- // if we have ghost values, do not allow
- // writing to this vector at all.
- Assert (!has_ghost_elements(), ExcGhostsPresent());
- Assert (v.local_size() == w.local_size(),
- ExcDimensionMismatch (v.local_size(), w.local_size()));
-
- AssertIsFinite(a);
- AssertIsFinite(b);
-
- // If we don't have the same map, copy.
- if (vector->Map().SameAs(v.vector->Map())==false)
- {
- sadd(0., a, v, b, w);
- }
- else
- {
- // Otherwise, just update. verify
- // that *this does not only have
- // the same map as v (the
- // if-condition above) but also as
- // w
- Assert (vector->Map().SameAs(w.vector->Map()),
- ExcDifferentParallelPartitioning());
- int ierr = vector->Update(a, *v.vector, b, *w.vector, 0.0);
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-
- last_action = Zero;
- }
- }
-
-
-
- // TODO: up to now only local
- // data printed out! Find a
- // way to neatly output
- // distributed data...
- void
- VectorBase::print (const char *format) const
- {
- Assert (global_length(*vector)!=0, ExcEmptyObject());
- (void)global_length;
-
- for (size_type j=0; j<size(); ++j)
- {
- double t = (*vector)[0][j];
-
- if (format != nullptr)
- std::printf (format, t);
- else
- std::printf (" %5.2f", double(t));
- }
- std::printf ("\n");
- }
-
-
-
void
VectorBase::print (std::ostream &out,
const unsigned int precision,
projects / dealii.git / commitdiff