From 8e152fd47d7bc341481dac912bd1e5714f4cf4a5 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Thu, 4 May 2017 22:53:24 +0200 Subject: [PATCH] Remove deprecated functions in TrilinosWrappers vectors --- .../fe/fe_tools_interpolate.templates.h | 4 +- include/deal.II/lac/trilinos_block_vector.h | 139 +---------- .../lac/trilinos_parallel_block_vector.h | 20 -- include/deal.II/lac/trilinos_vector.h | 179 -------------- include/deal.II/lac/trilinos_vector_base.h | 216 ---------------- include/deal.II/lac/vector.templates.h | 6 +- source/lac/CMakeLists.txt | 1 - source/lac/trilinos_block_vector.cc | 24 -- source/lac/trilinos_parallel_block_vector.cc | 52 ---- source/lac/trilinos_vector.cc | 234 ++---------------- source/lac/trilinos_vector_base.cc | 79 ------ 11 files changed, 23 insertions(+), 931 deletions(-) delete mode 100644 source/lac/trilinos_parallel_block_vector.cc diff --git a/include/deal.II/fe/fe_tools_interpolate.templates.h b/include/deal.II/fe/fe_tools_interpolate.templates.h index 2eddceae25..4408d76724 100644 --- a/include/deal.II/fe/fe_tools_interpolate.templates.h +++ b/include/deal.II/fe/fe_tools_interpolate.templates.h @@ -579,8 +579,8 @@ namespace FETools // 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); diff --git a/include/deal.II/lac/trilinos_block_vector.h b/include/deal.II/lac/trilinos_block_vector.h index 6d011b3730..b10402aaf1 100644 --- a/include/deal.II/lac/trilinos_block_vector.h +++ b/include/deal.II/lac/trilinos_block_vector.h @@ -97,61 +97,6 @@ namespace TrilinosWrappers */ 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 &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 &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 num_blocks 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 n.size() and - * initialize each block with n[i] zero elements. - * - * References BlockVector.reinit(). - */ - explicit BlockVector (const std::vector &N) DEAL_II_DEPRECATED; - - /** - * Constructor. Set the number of blocks to n.size(). 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 - * std::vector class, but the first argument is needed in order - * to know how to subdivide the block vector into different blocks. - */ - template - BlockVector (const std::vector &n, - const InputIterator first, - const InputIterator end) DEAL_II_DEPRECATED; - /** * Destructor. Clears memory */ @@ -317,93 +262,13 @@ namespace TrilinosWrappers inline - BlockVector::BlockVector () - {} + BlockVector::BlockVector () = default; inline - BlockVector::BlockVector (const std::vector &partitioning) - { - reinit (partitioning); - } - - + BlockVector::~BlockVector() = default; - inline - BlockVector::BlockVector (const std::vector &partitioning, - const MPI_Comm &communicator) - { - reinit (partitioning, communicator); - } - - - - inline - BlockVector::BlockVector (const std::vector &N) - { - reinit (N); - } - - - - template - BlockVector::BlockVector (const std::vector &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[b])); - - for (size_type i=0; iblock(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 () - { - this->components.resize (v.n_blocks()); - this->block_indices = v.block_indices; - - for (size_type i=0; in_blocks(); ++i) - this->components[i] = v.components[i]; - } inline diff --git a/include/deal.II/lac/trilinos_parallel_block_vector.h b/include/deal.II/lac/trilinos_parallel_block_vector.h index 0597cce0a9..81b22369a8 100644 --- a/include/deal.II/lac/trilinos_parallel_block_vector.h +++ b/include/deal.II/lac/trilinos_parallel_block_vector.h @@ -98,15 +98,6 @@ namespace TrilinosWrappers */ 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 ¶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 @@ -266,17 +257,6 @@ namespace TrilinosWrappers 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 debug 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. * diff --git a/include/deal.II/lac/trilinos_vector.h b/include/deal.II/lac/trilinos_vector.h index a4fe2c3e6e..2d85c1e8cf 100644 --- a/include/deal.II/lac/trilinos_vector.h +++ b/include/deal.II/lac/trilinos_vector.h @@ -381,100 +381,6 @@ namespace TrilinosWrappers 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 - void reinit (const Epetra_Map ¶llel_partitioner, - const dealii::Vector &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 - Vector (const Epetra_Map ¶llel_partitioning, - const dealii::Vector &v) DEAL_II_DEPRECATED; //@} /** * @name Initialization with an IndexSet @@ -633,15 +539,6 @@ namespace TrilinosWrappers #ifndef DOXYGEN - template - Vector::Vector (const Epetra_Map &input_map, - const dealii::Vector &v) - { - reinit (input_map, v); - } - - - template Vector::Vector (const IndexSet ¶llel_partitioner, const dealii::Vector &v, @@ -654,24 +551,6 @@ namespace TrilinosWrappers - - template - void Vector::reinit (const Epetra_Map ¶llel_partitioner, - const dealii::Vector &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; ireinit() 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 - explicit Vector (const dealii::Vector &v) DEAL_II_DEPRECATED; - /** * Reinit function that resizes the vector to the size specified by * n. The flag omit_zeroing_entries specifies whether @@ -916,16 +747,6 @@ namespace TrilinosWrappers #ifndef DOXYGEN - template - Vector::Vector (const dealii::Vector &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) diff --git a/include/deal.II/lac/trilinos_vector_base.h b/include/deal.II/lac/trilinos_vector_base.h index 082c69a93f..3a8718b5b0 100644 --- a/include/deal.II/lac/trilinos_vector_base.h +++ b/include/deal.II/lac/trilinos_vector_base.h @@ -289,14 +289,6 @@ namespace TrilinosWrappers */ 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 @@ -436,13 +428,6 @@ namespace TrilinosWrappers */ 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. */ @@ -573,18 +558,6 @@ namespace TrilinosWrappers const ForwardIterator indices_end, OutputIterator values_begin) const; - /** - * Return the value of the vector entry 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 vector<> class of the C++ * standard library by returning iterators to the start and end of the @@ -740,31 +713,6 @@ namespace TrilinosWrappers 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. *this = s*(*this) + a*V + b*W + - * c*X. - * - * 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 @@ -777,29 +725,6 @@ namespace TrilinosWrappers */ void equ (const TrilinosScalar a, const VectorBase &V); - - /** - * Assignment *this = a*V + b*W. - * - * 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 - * this[i] = a[i]/b[i]. 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 b[i] 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; //@} @@ -826,14 +751,6 @@ namespace TrilinosWrappers */ 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 @@ -1068,15 +985,6 @@ namespace TrilinosWrappers - inline - bool - VectorBase::is_compressed () const - { - return compressed; - } - - - inline bool VectorBase::in_local_range (const size_type index) const @@ -1485,15 +1393,6 @@ namespace TrilinosWrappers - inline - TrilinosScalar - VectorBase::minimal_value () const - { - return min(); - } - - - inline TrilinosScalar VectorBase::min () const @@ -1798,103 +1697,6 @@ namespace TrilinosWrappers - 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) @@ -1939,24 +1741,6 @@ namespace TrilinosWrappers - 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 diff --git a/include/deal.II/lac/vector.templates.h b/include/deal.II/lac/vector.templates.h index df73ec9ed9..73ce84f6e1 100644 --- a/include/deal.II/lac/vector.templates.h +++ b/include/deal.II/lac/vector.templates.h @@ -179,7 +179,8 @@ Vector::Vector (const TrilinosWrappers::MPI::Vector &v) // 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 @@ -895,7 +896,8 @@ Vector::operator= (const TrilinosWrappers::MPI::Vector &v) // 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; } diff --git a/source/lac/CMakeLists.txt b/source/lac/CMakeLists.txt index 7bfd5c700b..5597406a10 100644 --- a/source/lac/CMakeLists.txt +++ b/source/lac/CMakeLists.txt @@ -115,7 +115,6 @@ IF(DEAL_II_WITH_TRILINOS) 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 diff --git a/source/lac/trilinos_block_vector.cc b/source/lac/trilinos_block_vector.cc index 8897cb05a9..3f0b2fc64a 100644 --- a/source/lac/trilinos_block_vector.cc +++ b/source/lac/trilinos_block_vector.cc @@ -104,30 +104,6 @@ namespace TrilinosWrappers - void - BlockVector::reinit (const std::vector &input_maps, - const bool omit_zeroing_entries) - { - const size_type no_blocks = input_maps.size(); - std::vector block_sizes (no_blocks); - - for (size_type i=0; iblock_indices.reinit (block_sizes); - if (components.size() != n_blocks()) - components.resize(n_blocks()); - - for (size_type i=0; i ¶llel_partitioning, const MPI_Comm &communicator, diff --git a/source/lac/trilinos_parallel_block_vector.cc b/source/lac/trilinos_parallel_block_vector.cc deleted file mode 100644 index 1017f44557..0000000000 --- a/source/lac/trilinos_parallel_block_vector.cc +++ /dev/null @@ -1,52 +0,0 @@ -// --------------------------------------------------------------------- -// -// 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 - -#ifdef DEAL_II_WITH_TRILINOS - -DEAL_II_NAMESPACE_OPEN - -namespace TrilinosWrappers -{ - namespace MPI - { - BlockVector::BlockVector (const std::vector ¶llel_partitioning) - { - reinit (parallel_partitioning, false); - } - - - - bool - BlockVector::is_compressed () const - { - bool compressed = true; - for (unsigned int row=0; rowMap()), - 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) @@ -191,62 +162,7 @@ namespace TrilinosWrappers - 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 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(&(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; @@ -362,91 +278,12 @@ namespace TrilinosWrappers 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 global_ids (n_elements, -1); - for (size_type block=0; block 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(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(&(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, @@ -461,9 +298,6 @@ namespace TrilinosWrappers 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 { @@ -472,7 +306,15 @@ namespace TrilinosWrappers 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); @@ -483,6 +325,11 @@ namespace TrilinosWrappers 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); @@ -621,57 +468,6 @@ namespace TrilinosWrappers - - 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*/) diff --git a/source/lac/trilinos_vector_base.cc b/source/lac/trilinos_vector_base.cc index ed7689b8f8..c4f0a7b979 100644 --- a/source/lac/trilinos_vector_base.cc +++ b/source/lac/trilinos_vector_base.cc @@ -244,23 +244,6 @@ namespace TrilinosWrappers - 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(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 { @@ -423,68 +406,6 @@ namespace TrilinosWrappers - 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