From 9e8ace47198b190f3ab772a9f7ac29324ab5db6a Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Sun, 14 Feb 2016 08:16:10 -0500 Subject: [PATCH] remove experimental vector classes --- include/deal.II/lac/la_vector.h | 439 ---------- include/deal.II/lac/la_vector.templates.h | 407 --------- include/deal.II/lac/read_write_vector.h | 774 ------------------ .../deal.II/lac/read_write_vector.templates.h | 287 ------- include/deal.II/lac/vector_space_vector.h | 191 ----- source/lac/CMakeLists.txt | 4 - source/lac/la_vector.cc | 25 - source/lac/la_vector.inst.in | 29 - source/lac/read_write_vector.cc | 23 - source/lac/read_write_vector.inst.in | 52 -- tests/lac/la_vector_accumulation_01.cc | 66 -- tests/lac/la_vector_accumulation_01.output | 2 - tests/lac/la_vector_add_and_dot.cc | 91 -- tests/lac/la_vector_add_and_dot.output | 25 - tests/lac/la_vector_large_numbers.cc | 80 -- tests/lac/la_vector_large_numbers.output | 2 - tests/lac/la_vector_norms.cc | 118 --- tests/lac/la_vector_norms.output | 2 - tests/lac/la_vector_output.cc | 87 -- tests/lac/la_vector_output.output | 2 - tests/lac/la_vector_vector.cc | 158 ---- tests/lac/la_vector_vector.output | 76 -- tests/lac/readwritevector_add.cc | 83 -- tests/lac/readwritevector_add.output | 49 -- tests/lac/readwritevector_assignement.cc | 67 -- tests/lac/readwritevector_assignement.output | 55 -- tests/lac/readwritevector_assignment.cc | 67 -- tests/lac/readwritevector_assignment.output | 55 -- tests/trilinos/readwritevector.cc | 92 --- .../trilinos/readwritevector.mpirun=2.output | 2 - 30 files changed, 3410 deletions(-) delete mode 100644 include/deal.II/lac/la_vector.h delete mode 100644 include/deal.II/lac/la_vector.templates.h delete mode 100644 include/deal.II/lac/read_write_vector.h delete mode 100644 include/deal.II/lac/read_write_vector.templates.h delete mode 100644 include/deal.II/lac/vector_space_vector.h delete mode 100644 source/lac/la_vector.cc delete mode 100644 source/lac/la_vector.inst.in delete mode 100644 source/lac/read_write_vector.cc delete mode 100644 source/lac/read_write_vector.inst.in delete mode 100644 tests/lac/la_vector_accumulation_01.cc delete mode 100644 tests/lac/la_vector_accumulation_01.output delete mode 100644 tests/lac/la_vector_add_and_dot.cc delete mode 100644 tests/lac/la_vector_add_and_dot.output delete mode 100644 tests/lac/la_vector_large_numbers.cc delete mode 100644 tests/lac/la_vector_large_numbers.output delete mode 100644 tests/lac/la_vector_norms.cc delete mode 100644 tests/lac/la_vector_norms.output delete mode 100644 tests/lac/la_vector_output.cc delete mode 100644 tests/lac/la_vector_output.output delete mode 100644 tests/lac/la_vector_vector.cc delete mode 100644 tests/lac/la_vector_vector.output delete mode 100644 tests/lac/readwritevector_add.cc delete mode 100644 tests/lac/readwritevector_add.output delete mode 100644 tests/lac/readwritevector_assignement.cc delete mode 100644 tests/lac/readwritevector_assignement.output delete mode 100644 tests/lac/readwritevector_assignment.cc delete mode 100644 tests/lac/readwritevector_assignment.output delete mode 100644 tests/trilinos/readwritevector.cc delete mode 100644 tests/trilinos/readwritevector.mpirun=2.output diff --git a/include/deal.II/lac/la_vector.h b/include/deal.II/lac/la_vector.h deleted file mode 100644 index cfb00525fb..0000000000 --- a/include/deal.II/lac/la_vector.h +++ /dev/null @@ -1,439 +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. -// -// --------------------------------------------------------------------- - -#ifndef dealii__la_vector_h -#define dealii__la_vector_h - - -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include -#include -#include - - -DEAL_II_NAMESPACE_OPEN - -namespace LinearAlgebra -{ - /*! @addtogroup Vectors - *@{ - */ - - /** - * Numerical vector of data. This class derives from both - * ::dealii::LinearAlgebra::ReadWriteVector and - * ::dealii::LinearAlgebra::VectorSpaceVector. As opposed to the array of - * the C++ standard library, this class implements an element of a vector - * space suitable for numerical computations. - * - * @author Bruno Turcksin, 2015. - */ - template - class Vector : public ReadWriteVector, public VectorSpaceVector - { - public: - typedef types::global_dof_index size_type; - - /** - * Constructor. Create a vector of dimension zero. - */ - Vector(); - - /** - * Copy constructor. Sets the dimension to that of the given vector and - * copies all elements. - */ - Vector(const Vector &V); - - /** - * Constructor. Set dimension to @p n and initialize all elements with - * zero. - * - * The constructor is made explicit to avoid accident like this: - * v=0;. Presumably, the user wants to set every element of the - * vector to zero, but instead, what happens is this call: - * v=Vector@(0);, i.e. the vector is replaced by one of - * length zero. - */ - explicit Vector(const size_type n); - - /** - * Initialize the vector with a given range of values pointed to by the - * iterators. This function exists in analogy to the @p std::vector class. - */ - template - Vector(const InputIterator first, const InputIterator last); - - /** - * Destructor, deallocates memory. - */ - virtual ~Vector(); - - /** - * Copies the data of the input vector @p in_vector. - */ - Vector &operator= (const Vector &in_vector); - - /** - * Copies the data of the input vector @p in_vector. - */ - template - Vector &operator= (const Vector &in_vector); - - /** - * Sets all elements of the vector to the scalar @p s. This operation is - * only allowed if @p s is equal to zero. - */ - Vector &operator= (const Number s); - - /** - * Multiply the entire vector by a fixed factor. - */ - virtual VectorSpaceVector &operator*= (const Number factor); - - /** - * Divide the entire vector by a fixed factor. - */ - virtual VectorSpaceVector &operator/= (const Number factor); - - /** - * Add the vector @p V to the present one. - */ - virtual VectorSpaceVector &operator+= (const VectorSpaceVector &V); - - /** - * Substract the vector @p V from the present one. - */ - virtual VectorSpaceVector &operator-= (const VectorSpaceVector &V); - - /** - * Return the scalar product of two vectors. - */ - virtual Number operator* (const VectorSpaceVector &V); - - /** - * Add @p a to all components. Note that @p a is a scalar not a vector. - */ - virtual void add(const Number a); - - /** - * Simple addition of a multiple of a vector, i.e. *this += a*V. - */ - virtual void add(const Number a, const VectorSpaceVector &V); - - /** - * Multiple addition of a multiple of a vector, i.e. *this += - * a*V+b*W. - */ - virtual void add(const Number a, const VectorSpaceVector &V, - const Number b, const VectorSpaceVector &W); - - /** - * Scaling and simple addition of a multiple of a vector, i.e. *this = - * s*(*this)+a*V. - */ - virtual void sadd(const Number s, const Number a, - const VectorSpaceVector &V); - - /** - * Scale each element of this vector by the corresponding element in the - * argument. This function is mostly meant to simulate multiplication (and - * immediate re-assignement) by a diagonal scaling matrix. - */ - virtual void scale(const VectorSpaceVector &scaling_factors); - - /** - * Assignement *this = a*V. - */ - virtual void equ(const Number a, const VectorSpaceVector &V); - - /** - * Return the l1 norm of the vector (i.e., the sum of the - * absolute values of all entries). - */ - virtual typename VectorSpaceVector::real_type l1_norm(); - - /** - * Return the l2 norm of the vector (i.e., the square root of - * the sum of the square of all entries among all processors). - */ - virtual typename VectorSpaceVector::real_type l2_norm(); - - /** - * Return the maximum norm of the vector (i.e., the maximum absolute value - * among all entries and among all processors). - */ - virtual typename VectorSpaceVector::real_type linfty_norm(); - - /** - * 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 - */ - virtual Number add_and_dot(const Number a, - const VectorSpaceVector &V, - const VectorSpaceVector &W); - - /** - * Return the global size of the vector, equal to the sum of the number of - * locally owned indices among all processors. - */ - size_type size() const; - - /** - * 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 procesors if this is a distributed 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 - */ - const dealii::IndexSet &locally_owned_elements() const; - - /** - * Prints 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; - /** - * Write the vector en bloc to a file. This is done in a binary mode, so - * the output is neither readable by humans nor (probably) by other - * computers using a different operating system or number format. - */ - void block_write (std::ostream &out) const; - - /** - * Read a vector en block from a file. This is done using the inverse - * operations to the above function, so it is reasonably fast because the - * bitstream is not interpreted. - * - * The vector is resized if necessary. - * - * A primitive form of error checking is performed which will recognize - * the bluntest attempts to interpret some data as a vector stored bitwise - * to a file, but not more. - */ - void block_read (std::istream &in); - - /** - * Returns the memory consumption of this class in bytes. - */ - std::size_t memory_consumption() const; - - /** - * Attempt to perform an operation between two incompatible vector types. - * - * @ingroup Exceptions - */ - DeclException0(ExcVectorTypeNotCompatible); - - private: - /** - * Compute the L1 norm in a recursive way by dividing the vector on - * smaller and smaller intervals. This reduces the numerical error on - * large vector. - */ - typename VectorSpaceVector::real_type l1_norm_recursive(unsigned int i, - unsigned int j); - - /** - * Compute the squared L2 norm in a recursive way by dividing the vector - * on smaller and smaller intervals. This reduces the numerical error on - * large vector. - */ - typename VectorSpaceVector::real_type l2_norm_squared_recursive( - unsigned int i, - unsigned int j); - - /** - * Serialize the data of this object using boost. This function is - * necessary to use boost::archive::text_iarchive and - * boost::archive::text_oarchive. - */ - template - void serialize(Archive &ar, const unsigned int version); - - friend class boost::serialization::access; - }; - - /*@}*/ - /*----------------------- Inline functions ----------------------------------*/ - - template - inline - Vector::Vector() {} - - - - template - inline - Vector::Vector(const Vector &V) - : - ReadWriteVector(V) - {} - - - - template - inline - Vector::Vector(const size_type n) - : - ReadWriteVector(n) - {} - - - - template - template - inline - Vector::Vector(const InputIterator first, const InputIterator last) - { - this->reinit(complete_index_set(std::distance (first, last)), true); - std::copy(first, last, this->begin()); - } - - - - template - inline - Vector::~Vector() {} - - - - template - inline - Vector &Vector::operator= (const Vector &in_vector) - { - this->reinit(in_vector.size(),true); - std::copy(in_vector.begin(), in_vector.end(), this->begin()); - - return *this; - } - - - - template - template - inline - Vector &Vector::operator= (const Vector &in_vector) - { - this->reinit(in_vector.size(in_vector.size()),true); - std::copy(in_vector.begin(), in_vector.end(), this->begin()); - - return *this; - } - - - - template - inline - Vector &Vector::operator= (const Number s) - { - Assert(s==static_cast(0), ExcMessage("Only 0 can be assigned to a vector.")); - (void) s; - - std::fill(this->begin(),this->end(),Number()); - - return *this; - } - - - - template - inline - void Vector::add(const Number a) - { - AssertIsFinite(a); - for (unsigned int i=0; isize(); ++i) - this->val[i] += a; - } - - - - template - inline - typename Vector::size_type Vector::size() const - { - return ReadWriteVector::size(); - } - - - - template - inline - const dealii::IndexSet &Vector::locally_owned_elements() const - { - return ReadWriteVector::get_stored_elements(); - } - - - - template - inline - void Vector::print(std::ostream &out, - const unsigned int precision, - const bool scientific, - const bool) const - { - ReadWriteVector::print(out, precision, scientific); - } - - - - template - template - inline - void Vector::serialize(Archive &ar, const unsigned int) - { - unsigned int current_size = this->size(); - ar &static_cast(*this); - ar &this->stored_elements; - // If necessary, resize the vector during a read operation - if (this->size()!=current_size) - this->reinit(this->size()); - ar &boost::serialization::make_array(this->val,this->size()); - } - - - - template - inline - std::size_t Vector::memory_consumption() const - { - return ReadWriteVector::memory_consumption(); - } -} // end of namespace LinearAlgebra - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/lac/la_vector.templates.h b/include/deal.II/lac/la_vector.templates.h deleted file mode 100644 index 9c69b9c78d..0000000000 --- a/include/deal.II/lac/la_vector.templates.h +++ /dev/null @@ -1,407 +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. -// -// --------------------------------------------------------------------- - -#ifndef dealii__la_vector_templates_h -#define dealii__la_vector_templates_h - -#include -#include -#include - -DEAL_II_NAMESPACE_OPEN - -namespace LinearAlgebra -{ - template - VectorSpaceVector &Vector::operator*= (const Number factor) - { - AssertIsFinite(factor); - for (unsigned int i=0; isize(); ++i) - this->val[i] *= factor; - - return *this; - } - - - - template - VectorSpaceVector &Vector::operator/= (const Number factor) - { - AssertIsFinite(factor); - Assert(factor!=Number(0.), ExcZero()); - this->operator *= (Number(1.)/factor); - - return *this; - } - - - - template - VectorSpaceVector &Vector::operator+= (const VectorSpaceVector &V) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - Assert(down_V.size()==this->size(), - ExcMessage("Cannot add two vectors with different numbers of elements")); - for (unsigned int i=0; isize(); ++i) - { - AssertIsFinite(down_V[i]); - this->val[i] += down_V[i]; - } - - return *this; - } - - - - template - VectorSpaceVector &Vector::operator-= (const VectorSpaceVector &V) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - Assert(down_V.size()==this->size(), - ExcMessage("Cannot subtract two vectors with different numbers of elements")); - for (unsigned int i=0; isize(); ++i) - { - AssertIsFinite(down_V[i]); - this->val[i] -= down_V[i]; - } - - return *this; - } - - - - template - Number Vector::operator* (const VectorSpaceVector &V) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - Assert(down_V.size()==this->size(), - ExcMessage("Cannot compute the scalar product " - "of two vectors with different numbers of elements")); - Number value = 0.; - for (unsigned int i=0; isize(); ++i) - { - AssertIsFinite(down_V[i]); - value += this->val[i]*down_V[i]; - } - - return value; - } - - - - template - void Vector::add(const Number a, const VectorSpaceVector &V) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - AssertIsFinite(a); - Assert(down_V.size()==this->size(), - ExcMessage("Cannot add two vectors with different numbers of elements")); - for (unsigned int i=0; isize(); ++i) - { - AssertIsFinite(down_V[i]); - this->val[i] += a*down_V[i]; - } - } - - - - template - void Vector::add(const Number a, const VectorSpaceVector &V, - const Number b, const VectorSpaceVector &W) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, - ExcVectorTypeNotCompatible()); - // Check that casting will work. - Assert(dynamic_cast*>(&W)!=NULL, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - // Downcast W. If fails, throws an exception. - const Vector &down_W = dynamic_cast&>(W); - AssertIsFinite(a); - Assert(down_V.size()==this->size(), - ExcMessage("Cannot add two vectors with different numbers of elements")); - AssertIsFinite(b); - Assert(down_W.size()==this->size(), - ExcMessage("Cannot add two vectors with different numbers of elements")); - for (unsigned int i=0; isize(); ++i) - { - AssertIsFinite(down_V[i]); - AssertIsFinite(down_W[i]); - this->val[i] += a*down_V[i]+b*down_W[i]; - } - } - - - - template - void Vector::sadd(const Number s, const Number a, - const VectorSpaceVector &V) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, - ExcVectorTypeNotCompatible()); - - *this *= s; - // Downcast V. It fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - Vector tmp(down_V); - tmp *= a; - *this += tmp; - } - - - - template - void Vector::scale(const VectorSpaceVector &scaling_factors) - { - // Check that casting will work. - Assert(dynamic_cast*>(&scaling_factors)!=NULL, - ExcVectorTypeNotCompatible()); - - // Downcast scaling_factors. If fails, throws an exception. - const Vector &down_scaling_factors = - dynamic_cast&>(scaling_factors); - Assert(down_scaling_factors.size()==this->size(), - ExcMessage("Cannot add two vectors with different numbers of elements")); - for (unsigned int i=0; isize(); ++i) - { - AssertIsFinite(down_scaling_factors[i]); - this->val[i] *= down_scaling_factors[i]; - } - } - - - - template - void Vector::equ(const Number a, const VectorSpaceVector &V) - { - // Check that casting will work. - Assert(dynamic_cast*>(&V)!=NULL, - ExcVectorTypeNotCompatible()); - - // Downcast V. If fails, throws an exception. - const Vector &down_V = dynamic_cast&>(V); - *this = down_V; - *this *= a; - } - - - - template - typename VectorSpaceVector::real_type Vector::l1_norm() - { - Assert (this->size(), ExcEmptyObject()); - - typename ReadWriteVector::real_type norm = 0.; - norm = l1_norm_recursive(0,this->size()-1); - - return norm; - } - - - - - template - typename VectorSpaceVector::real_type Vector::l2_norm() - { - Assert (this->size(), ExcEmptyObject()); - - // if l2_norm()^2 is finite and non-zero, the answer is computed as - // std::sqrt(norm_sqr()). If norm_sqrt() is infinite or zero, the l2 norm - // might still be finite. In that case, recompute ig (this is a rare case, - // so working on the vector twice is uncritical and paid off by the extended - // precision) using the BLAS approach with a weight, see e.g. dnrm2.f. - typedef typename ReadWriteVector::real_type real_type; - real_type norm_square = 0.; - norm_square = l2_norm_squared_recursive(0,this->size()-1); - if (numbers::is_finite(norm_square) && - norm_square>=std::numeric_limits::min()) - return std::sqrt(norm_square); - else - { - real_type scale = 0.; - real_type sum = 1.; - for (unsigned int i=0; isize(); ++i) - { - if (this->val[i] != Number()) - { - const real_type abs_x = - numbers::NumberTraits::abs(this->val[i]); - if (scale < abs_x) - { - sum = 1. + sum * (scale/abs_x) * (scale/abs_x); - scale = abs_x; - } - else - sum += (abs_x/scale) * (abs_x/scale); - } - } - AssertIsFinite(scale*std::sqrt(sum)); - return scale*std::sqrt(sum); - } - } - - - - template - typename VectorSpaceVector::real_type Vector::linfty_norm() - { - typename ReadWriteVector::real_type norm = 0.; - for (unsigned int i=0; isize(); ++i) - norm = std::max(std::abs(this->val[i]),norm); - - return norm; - } - - - - template - Number Vector::add_and_dot(const Number a, - const VectorSpaceVector &V, - const VectorSpaceVector &W) - { - this->add(a,V); - - return *this*W; - } - - - - template - typename VectorSpaceVector::real_type Vector::l1_norm_recursive( - unsigned int i, - unsigned int j) - { - Assert(j>=i, ExcInternalError()); - typename ReadWriteVector::real_type norm = 0.; - - if ((j-i)!=0) - { - norm += l1_norm_recursive(i,(i+j)/2); - norm += l1_norm_recursive((i+j)/2+1,j); - } - else - norm += std::abs(this->val[i]); - - return norm; - } - - - - template - typename VectorSpaceVector::real_type Vector::l2_norm_squared_recursive( - unsigned int i, - unsigned int j) - { - Assert(j>=i, ExcInternalError()); - - typename ReadWriteVector::real_type norm = 0.; - if ((j-i)!=0) - { - norm += l2_norm_squared_recursive(i,(i+j)/2); - norm += l2_norm_squared_recursive((i+j)/2+1,j); - } - else - norm += std::pow(std::abs(this->val[i]),2); - - return norm; - } - - - - template - void Vector::block_write(std::ostream &out) const - { - AssertThrow(out, ExcIO()); - - // Other version of the following - // out << size() << std::endl << '['; - // Reason: operator<< seems to use some resources that lead to problems in - // a multithreaded environment. - const size_type sz = this->size(); - char buf[16]; -#ifdef DEAL_II_WITH_64BIT_INDICES - std::sprintf(buf, "%llu", sz); -#else - std::sprintf(buf, "%u", sz); -#endif - std::strcat(buf, "\n["); - - out.write(buf, std::strlen(buf)); - out.write(reinterpret_cast(this->begin()), - reinterpret_cast(this->end()) - - reinterpret_cast(this->begin())); - - // out << ']'; - const char outro = ']'; - out.write(&outro, 1); - - AssertThrow(out, ExcIO()); - } - - - - template - void Vector::block_read(std::istream &in) - { - AssertThrow(in, ExcIO()); - - size_type sz; - - char buf[16]; - - in.getline(buf,16,'\n'); - sz = std::atoi(buf); - this->reinit(sz,true); - - char c; - // in >> c; - in.read(&c,1); - AssertThrow(c=='[', ExcIO()); - - in.read(reinterpret_cast(this->begin()), - reinterpret_cast(this->end()) - - reinterpret_cast(this->begin())); - - // in >> c; - in.read(&c,1); - AssertThrow(c=']', ExcIO()); - } -} - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/lac/read_write_vector.h b/include/deal.II/lac/read_write_vector.h deleted file mode 100644 index 207ccb289c..0000000000 --- a/include/deal.II/lac/read_write_vector.h +++ /dev/null @@ -1,774 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2015 - 2016 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. -// -// --------------------------------------------------------------------- - -#ifndef dealii__read_write_vector_h -#define dealii__read_write_vector_h - -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include - - -DEAL_II_NAMESPACE_OPEN - -#ifdef DEAL_II_WITH_PETSC -namespace PETScWrappers -{ - namespace MPI - { - class Vector; - } -} -#endif - -#ifdef DEAL_II_WITH_TRILINOS -namespace TrilinosWrappers -{ - namespace MPI - { - class Vector; - } -} -#endif - -namespace LinearAlgebra -{ - /*! @addtogroup Vectors - *@{ - */ - - /** - * ReadWriteVector is intended to represent vectors in ${\mathbb R}^N$ for - * which it stores all or a subset of elements. The latter case in important - * in parallel computations, where $N$ may be so large that no processor can - * actually all elements of a solution vector, but where this is also not - * necessary: one typically only has to store the values of degrees of - * freedom that live on cells that are locally owned plus potentially those - * degrees of freedom that live on ghost cells. - * - * This class allows to access individual elements to be read or written. - * However, it does not allow global operations such as taking the norm. - * ReadWriteVector can be used to read and write elements in vectors derived - * from VectorSpaceVector such as TrilinosWrappers::MPI::Vector and - * PETScWrappers::MPI::Vector. - * - *

Storing elements

Most of the time, one will simply read from or - * write into a vector of the current class using the global numbers of - * these degrees of freedom. This is done using operator() or operator[] - * which call global_to_local() to transform the global index into a - * local one. In such cases, it is clear that one can only access - * elements of the vector that the current object indeed stores. - * - * However, it is also possible to access elements in the order in which - * they are stored by the current object. In other words, one is not - * interested in accessing elements with their global indices, but - * instead using an enumeration that only takes into account the elements - * that are actually stored. This is facilitated by the local_element() - * function. To this end, it is necessary to know in which order the - * current class stores its element. The elements of all the consecutive - * ranges are stored in ascending order of the first index of each range. - * The function largest_range_starting_index() of IndexSet can be used to - * get the first index of the largest range. - * - * @author Bruno Turcksin, 2015. - */ - template - class ReadWriteVector : public Subscriptor - { - public: - /** - * Declare standard types used in all containers. These types parallel - * those in the C++ standard libraries vector<...> - * class. - */ - typedef Number value_type; - typedef value_type *pointer; - typedef const value_type *const_pointer; - typedef value_type *iterator; - typedef const value_type *const_iterator; - typedef value_type &reference; - typedef const value_type &const_reference; - typedef types::global_dof_index size_type; - typedef typename numbers::NumberTraits::real_type real_type; - - /** - * @name 1: Basic Object-handling - */ - //@{ - /** - * Empty constructor. - */ - ReadWriteVector (); - - /** - * Copy constructor. - */ - ReadWriteVector (const ReadWriteVector &in_vector); - - /** - * Constructs a vector given the size, the stored elements have their - * index in [0,size). - */ - explicit ReadWriteVector (const size_type size); - - /** - * Constructs a vector whose stored elements indices are given by the - * IndexSet @p locally_stored_indices. - */ - explicit ReadWriteVector (const IndexSet &locally_stored_indices); - - /** - * Destructor. - */ - ~ReadWriteVector (); - - /** - * Sets the global size of the vector to @p size. The stored elements have - * their index in [0,size). - * - * 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). - */ - void reinit (const size_type size, - const bool omit_zeroing_entries = false); - - /** - * Uses the same IndexSet as the one of the input vector @p in_vector and - * allocates memory for this vector. - * - * 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 - void reinit(const ReadWriteVector &in_vector, - const bool omit_zeroing_entries = false); - - /** - * Initializes the vector. The indices are specified by @p - * locally_stored_indices. - * - * 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). - * locally_stored_indices. - */ - void reinit (const IndexSet &locally_stored_indices, - const bool omit_zeroing_entries = false); - - /** - * 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 analog to the the @p swap function of all C++ standard - * containers. Also, there is a global function swap(u,v) that - * simply calls u.swap(v), again in analogy to standard - * functions. - */ - void swap (ReadWriteVector &v); - - /** - * Copies the data and the IndexSet of the input vector @p in_vector. - */ - ReadWriteVector & - operator= (const ReadWriteVector &in_vector); - - /** - * Copies the data and the IndexSet of the input vector @p in_vector. - */ - template - ReadWriteVector & - operator= (const ReadWriteVector &in_vector); - -#ifdef DEAL_II_WITH_PETSC - /** - * Imports all the elements present in the vector's IndexSet from the - * input vector @p petsc_vec. - */ - ReadWriteVector & - operator= (const PETScWrappers::MPI::Vector &petsc_vec); -#endif - -#ifdef DEAL_II_WITH_TRILINOS - /** - * Imports all the elements present in the vector's IndexSet from the - * input vector @p trilinos_vec. - */ - ReadWriteVector & - operator= (const TrilinosWrappers::MPI::Vector &trilinos_vec); -#endif - - /** - * Sets all elements of the vector to the scalar @p s. This operation is - * only allowed if @p s is equal to zero. - */ - ReadWriteVector &operator = (const Number s); - - /** - * The value returned by this function denotes the dimension of the vector - * spaces that are modeled by objects of this kind. However, objects of - * the current class do not actually stores all elements of vectors of - * this space but may, in fact store only a subset. The number of elements - * stored is returned by n_elements() and is smaller or equal to the - * number returned by the current function. - */ - size_type size() const; - - /** - * This function returns the number of elements stored. It is smaller or - * equal to the dimension of the vector space that is modeled by an object - * of this kind. This dimension is return by size(). - */ - size_type n_elements() const; - - /** - * Return the IndexSet that represents the indices of the elements stored. - */ - const IndexSet &get_stored_elements () const; - - /** - * Make the @p ReadWriteVector class a bit like the vector<> - * class of the C++ standard library by returning iterators to the start - * and end of the locally stored elements of this vector. - */ - iterator begin (); - - /** - * Returns constant iterator to the start of the locally stored elements - * of the vector. - */ - const_iterator begin () const; - - /** - * Returns an iterator pointing to the element past the end of the array - * of locally stored entries. - */ - iterator end (); - - /** - * Returns a constant iterator pointing to the element past the end of the - * array of the locally stored entries. - */ - const_iterator end () const; - //@} - - - /** - * @name 2: Data-Access - */ - //@{ - - /** - * Read access to the data in the position corresponding to @p - * global_index. An exception is thrown if @p global_index is not stored - * by the current object. - */ - Number operator () (const size_type global_index) const; - - /** - * Read and write access to the data in the position corresponding to @p - * global_index. An exception is thrown if @p global_index is not stored - * by the current object. - */ - Number &operator () (const size_type global_index); - - /** - * Read access to the data in the position corresponding to @p - * global_index. An exception is thrown if @p global_index is not stored - * by the current object. - * - * 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. An exception is thrown if @p global_index is not stored - * by the current object. - * - * This function does the same thing as operator(). - */ - Number &operator [] (const size_type global_index); - - /** - * Instead of getting individual elements of a vector, this function - * allows to get 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. - */ - template - void extract_subvector_to (const std::vector &indices, - std::vector &values) const; - - /** - * Just as the above, but with pointers. Useful in minimizing copying of - * data around. - */ - template - void extract_subvector_to (ForwardIterator indices_begin, - const ForwardIterator indices_end, - OutputIterator values_begin) const; - - /** - * Read access to the data field specified by @p local_index. When you - * access elements in the order in which they are stored, it is necessary - * that you know in which they are stored. In other words, you need to - * know the map between the global indices of the elements this class - * stores, and the local indices into the contiguous array of these global - * elements. For this, see the general documentation of this class. - * - * Performance: Direct array access (fast). - */ - Number local_element (const size_type local_index) const; - - /** - * Read and write access to the data field specified by @p local_index. - * When you access elements in the order in which they are stored, it is - * necessary that you know in which they are stored. In other words, you - * need to know the map between the global indices of the elements this - * class stores, and the local indices into the contiguous array of these - * global elements. For this, see the general documentation of this class. - * - * Performance: Direct array access (fast). - */ - Number &local_element (const size_type local_index); - //@} - - - /** - * @name 3: Modification of vectors - */ - //@{ - - /** - * This function adds a whole set of values stored in @p values to the - * vector components specified by @p indices. - */ - template - void add (const std::vector &indices, - const std::vector &values); - - /** - * This function is similar to the previous one but takes a - * ReadWriteVector of values. - */ - template - void add (const std::vector &indices, - const ReadWriteVector &values); - - /** - * Take an address where n_elements are stored contiguously and - * add them into the vector. Handles all cases which are not covered by - * the other two add() functions above. - */ - template - void add (const size_type n_elements, - const size_type *indices, - const Number2 *values); - - /** - * Prints the vector to the output stream @p out. - */ - void print (std::ostream &out, - const unsigned int precision = 3, - const bool scientific = true) const; - - /** - * Returns the memory consumption of this class in bytes. - */ - std::size_t memory_consumption () const; - //@} - - protected: - - /** - * Return the local position of @p global_index. - */ - unsigned int - global_to_local (const types::global_dof_index global_index) const - { - // the following will throw an exception if the global_index is not - // in the remaining_elements - return static_cast(stored_elements.index_within_set(global_index)); - } - - /** - * A helper function that is used to resize the val array. - */ - void resize_val (const size_type new_allocated_size); - - /** - * Indices of the elements stored. - */ - IndexSet stored_elements; - - /** - * Pointer to the array of local elements of this vector. - */ - // TODO: use AlignedVector here for storage - Number *val; - }; - - /*@}*/ - - - /*----------------------- Inline functions ----------------------------------*/ - -#ifndef DOXYGEN - - template - inline - ReadWriteVector::ReadWriteVector () - : - val (NULL) - {} - - - - template - inline - ReadWriteVector::ReadWriteVector (const ReadWriteVector &v) - : - Subscriptor(), - val (NULL) - { - this->operator=(v); - } - - - - template - inline - ReadWriteVector::ReadWriteVector (const size_type size) - : - val (NULL) - { - reinit (size, false); - } - - - - template - inline - ReadWriteVector::ReadWriteVector (const IndexSet &locally_stored_indices) - : - val (NULL) - { - reinit (locally_stored_indices); - } - - - - template - inline - ReadWriteVector::~ReadWriteVector () - { - resize_val(0); - } - - - - template - inline - ReadWriteVector & - ReadWriteVector::operator= (const ReadWriteVector &in_vector) - { - resize_val(in_vector.n_elements()); - stored_elements = in_vector.get_stored_elements(); - std::copy(in_vector.begin(),in_vector.end(),begin()); - - return *this; - } - - - - template - template - inline - ReadWriteVector & - ReadWriteVector::operator= (const ReadWriteVector &in_vector) - { - resize_val(in_vector.n_elements()); - stored_elements = in_vector.get_stored_elements(); - std::copy(in_vector.begin(),in_vector.end(),begin()); - - return *this; - } - - - - template - inline - typename ReadWriteVector::size_type - ReadWriteVector::size() const - { - return stored_elements.size(); - } - - - - template - inline - typename ReadWriteVector::size_type - ReadWriteVector::n_elements() const - { - return stored_elements.n_elements(); - } - - - - template - inline - const IndexSet & - ReadWriteVector::get_stored_elements () const - { - return stored_elements; - } - - - - template - inline - typename ReadWriteVector::iterator - ReadWriteVector::begin () - { - return &val[0]; - } - - - - template - inline - typename ReadWriteVector::const_iterator - ReadWriteVector::begin () const - { - return &val[0]; - } - - - - template - inline - typename ReadWriteVector::iterator - ReadWriteVector::end () - { - return &val[this->n_elements()]; - } - - - - template - inline - typename ReadWriteVector::const_iterator - ReadWriteVector::end () const - { - return &val[this->n_elements()]; - } - - - - template - inline - Number - ReadWriteVector::operator() (const size_type global_index) const - { - return val[global_to_local(global_index)]; - } - - - - template - inline - Number & - ReadWriteVector::operator() (const size_type global_index) - { - return val[global_to_local (global_index)]; - } - - - - template - inline - Number - ReadWriteVector::operator[] (const size_type global_index) const - { - return operator()(global_index); - } - - - - template - inline - Number & - ReadWriteVector::operator[] (const size_type global_index) - { - return operator()(global_index); - } - - - - template - template - inline - void ReadWriteVector::extract_subvector_to (const std::vector &indices, - std::vector &values) const - { - for (size_type i = 0; i < indices.size(); ++i) - values[i] = operator()(indices[i]); - } - - - - template - template - inline - void ReadWriteVector::extract_subvector_to (ForwardIterator indices_begin, - const ForwardIterator indices_end, - OutputIterator values_begin) const - { - while (indices_begin != indices_end) - { - *values_begin = operator()(*indices_begin); - indices_begin++; - values_begin++; - } - } - - - - template - inline - Number - ReadWriteVector::local_element (const size_type local_index) const - { - AssertIndexRange (local_index, this->n_elements()); - - return val[local_index]; - } - - - - template - inline - Number & - ReadWriteVector::local_element (const size_type local_index) - { - AssertIndexRange (local_index, this->n_elements()); - - return val[local_index]; - } - - - - template - inline - ReadWriteVector & - ReadWriteVector::operator= (const Number s) - { - Assert(s==static_cast(0), ExcMessage("Only 0 can be assigned to a vector.")); - (void)s; - - std::fill(begin(),end(),Number()); - - return *this; - } - - - - template - template - inline - void - ReadWriteVector::add (const std::vector &indices, - const std::vector &values) - { - AssertDimension (indices.size(), values.size()); - add (indices.size(), &indices[0], &values[0]); - } - - - - template - template - inline - void - ReadWriteVector::add (const std::vector &indices, - const ReadWriteVector &values) - { - for (size_type i=0; ioperator()(indices[i]) += values[indices[i]]; - } - } - - - - template - template - inline - void - ReadWriteVector::add (const size_type n_indices, - const size_type *indices, - const Number2 *values) - { - for (size_type i=0; ioperator()(indices[i]) += values[i]; - } - } - -#endif // ifndef DOXYGEN - -} // end of 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. - * - * @relates Vector - */ -template -inline -void swap (LinearAlgebra::ReadWriteVector &u, - LinearAlgebra::ReadWriteVector &v) -{ - u.swap (v); -} - - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h deleted file mode 100644 index 9e334de71d..0000000000 --- a/include/deal.II/lac/read_write_vector.templates.h +++ /dev/null @@ -1,287 +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. -// -// --------------------------------------------------------------------- - -#ifndef dealii__parallel_vector_templates_h -#define dealii__parallel_vector_templates_h - - -#include -#include -#include -#include - -#ifdef DEAL_II_WITH_TRILINOS -# include -#endif - -DEAL_II_NAMESPACE_OPEN - - -namespace LinearAlgebra -{ - template - void - ReadWriteVector::resize_val (const size_type new_alloc_size) - { - if (new_alloc_size == 0) - { - if (val != NULL) - free(val); - val = NULL; - } - else - { - if (val != NULL) - free(val); - - Utilities::System::posix_memalign ((void **)&val, 64, sizeof(Number)*new_alloc_size); - } - } - - - - template - void - ReadWriteVector::reinit (const size_type size, - const bool omit_zeroing_entries) - { - // check whether we need to reallocate - resize_val(size); - - stored_elements = complete_index_set(size); - stored_elements.compress(); - - // set entries to zero if so requested - if (omit_zeroing_entries == false) - this->operator = (Number()); - } - - - - template - template - void - ReadWriteVector::reinit (const ReadWriteVector &v, - const bool omit_zeroing_entries) - { - resize_val(v.n_elements()); - - stored_elements = v.get_stored_elements(); - - if (omit_zeroing_entries == false) - this->operator= (Number()); - } - - - - template - void - ReadWriteVector::reinit (const IndexSet &locally_stored_indices, - const bool omit_zeroing_entries) - { - stored_elements = locally_stored_indices; - - // set vector size and allocate memory - resize_val(stored_elements.n_elements()); - - // initialize to zero - if (omit_zeroing_entries == false) - this->operator= (Number()); - } - - - -#ifdef DEAL_II_WITH_PETSC - namespace internal - { - template - 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 - void copy_petsc_vector (const std::complex *petsc_start_ptr, - const std::complex *petsc_end_ptr, - std::complex *ptr) - { - std::copy(petsc_start_ptr, petsc_end_ptr, ptr); - } - - template - void copy_petsc_vector (const std::complex *petsc_start_ptr, - const std::complex *petsc_end_ptr, - Number *ptr) - { - AssertThrow(false, ExcMessage("Tried to copy complex -> real")); - } - } - - - - template - ReadWriteVector & - ReadWriteVector::operator = (const PETScWrappers::MPI::Vector &petsc_vec) - { - //TODO: this works only if no communication is needed. - Assert(petsc_vec.locally_owned_elements() == stored_elements, - StandardExceptions::ExcInvalidState()); - - // get a representation of the vector and copy it - PetscScalar *start_ptr; - int ierr = VecGetArray (static_cast(petsc_vec), &start_ptr); - AssertThrow (ierr == 0, ExcPETScError(ierr)); - - const size_type vec_size = petsc_vec.local_size(); - internal::copy_petsc_vector (start_ptr, start_ptr + vec_size, begin()); - - // restore the representation of the vector - ierr = VecRestoreArray (static_cast(petsc_vec), &start_ptr); - AssertThrow (ierr == 0, ExcPETScError(ierr)); - - // return a pointer to this object per normal c++ operator overloading - // semantics - return *this; - } -#endif - - - -#ifdef DEAL_II_WITH_TRILINOS - template - ReadWriteVector & - ReadWriteVector::operator = (const TrilinosWrappers::MPI::Vector &trilinos_vec) - { - // Copy the local elements - std::vector trilinos_indices, readwrite_indices; - trilinos_vec.locally_owned_elements().fill_index_vector(trilinos_indices); - stored_elements.fill_index_vector(readwrite_indices); - std::vector intersection(trilinos_indices.size()); - std::vector::iterator end; - end = std::set_intersection(trilinos_indices.begin(), trilinos_indices.end(), - readwrite_indices.begin(), readwrite_indices.end(), - intersection.begin()); - const size_t intersection_size = end-intersection.begin(); - intersection.resize(intersection_size); - std::vector trilinos_position(intersection_size); - unsigned int j=0; - for (size_t i=0; i off_proc_indices(readwrite_indices.size()); - // Subtract the elements of readwrite that are in intersection. - end = std::set_difference(readwrite_indices.begin(), readwrite_indices.end(), - intersection.begin(), intersection.end(), off_proc_indices.begin()); - off_proc_indices.resize(end-off_proc_indices.begin()); - // Cast size_type to TrilinosWrappers::type::int_type - TrilinosWrappers::types::int_type *off_proc_array = - new TrilinosWrappers::types::int_type [off_proc_indices.size()]; - for (size_t i=0; i ( - off_proc_indices[i]); - Epetra_Map target_map(off_proc_indices.size(),off_proc_indices.size(), - off_proc_array,0,source_map.Comm()); - delete [] off_proc_array; - Epetra_Import import(target_map, source_map); - Epetra_FEVector target_vector(target_map); - target_vector.Import(trilinos_vec.trilinos_vector(), import, Insert); - values = target_vector.Values(); - for (unsigned int i=0; i - void - ReadWriteVector::swap (ReadWriteVector &v) - { - std::swap(stored_elements, v.stored_elements); - std::swap(val, v.val); - } - - - - template - std::size_t - ReadWriteVector::memory_consumption () const - { - std::size_t memory = sizeof(*this); - memory += sizeof (Number) * static_cast(this->n_elements()); - - memory += stored_elements.memory_consumption(); - - return memory; - } - - - - template - void - ReadWriteVector::print (std::ostream &out, - const unsigned int precision, - const bool scientific) const - { - 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); - - out << "IndexSet: "; - stored_elements.print(out); - out << std::endl; - for (unsigned int i=0; in_elements(); ++i) - out << val[i] << std::endl; - out << std::flush; - - - AssertThrow (out, ExcIO()); - // reset output format - out.flags (old_flags); - out.precision(old_precision); - } -} // end of namespace LinearAlgebra - - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/lac/vector_space_vector.h b/include/deal.II/lac/vector_space_vector.h deleted file mode 100644 index ff306146d8..0000000000 --- a/include/deal.II/lac/vector_space_vector.h +++ /dev/null @@ -1,191 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2015 - 2016 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. -// -// --------------------------------------------------------------------- - -#ifndef dealii__vector_space_vector_h -#define dealii__vector_space_vector_h - -#include -#include -#include - - -DEAL_II_NAMESPACE_OPEN - -class IndexSet; - -namespace LinearAlgebra -{ - template - class ReadWriteVector; - - - /*! @addtogroup Vectors - *@{ - */ - - /** - * VectorSpaceVector is an abstract class that is used to define the - * interface that vector classes need to implement when they want to - * implement global operations. This class is complementary of - * ReadWriteVector which allows the access of individual elements but does - * not allow global operations. - * - * @author Bruno Turcksin, 2015. - */ - template - class VectorSpaceVector - { - public: - typedef Number value_type; - typedef types::global_dof_index size_type; - typedef typename numbers::NumberTraits::real_type real_type; - - - /** - * Multiply the entire vector by a fixed factor. - */ - virtual VectorSpaceVector &operator*= (const Number factor) = 0; - - /** - * Divide the entire vector by a fixed factor. - */ - virtual VectorSpaceVector &operator/= (const Number factor) = 0; - - /** - * Add the vector @p V to the present one. - */ - virtual VectorSpaceVector &operator+= (const VectorSpaceVector &V) = 0; - - /** - * Subtract the vector @p V from the present one. - */ - virtual VectorSpaceVector &operator-= (const VectorSpaceVector &V) = 0; - - /** - * Return the scalar product of two vectors. - */ - virtual Number operator* (const VectorSpaceVector &V) = 0; - - /** - * Add @p a to all components. Note that @p a is a scalar not a vector. - */ - virtual void add(const Number a) = 0; - - /** - * Simple addition of a multiple of a vector, i.e. *this += a*V. - */ - virtual void add(const Number a, const VectorSpaceVector &V) = 0; - - /** - * Multiple addition of scaled vectors, i.e. *this += a*V+b*W. - */ - virtual void add(const Number a, const VectorSpaceVector &V, - const Number b, const VectorSpaceVector &W) = 0; - - /** - * Scaling and simple addition of a multiple of a vector, i.e. *this = - * s*(*this)+a*V. - */ - virtual void sadd(const Number s, const Number a, - const VectorSpaceVector &V) = 0; - - /** - * Scale each element of this vector by the corresponding element in the - * argument. This function is mostly meant to simulate multiplication (and - * immediate re-assignement) by a diagonal scaling matrix. - */ - virtual void scale(const VectorSpaceVector &scaling_factors) = 0; - - /** - * Assignement *this = a*V. - */ - virtual void equ(const Number a, const VectorSpaceVector &V) = 0; - - /** - * Return the l1 norm of the vector (i.e., the sum of the - * absolute values of all entries among all processors). - */ - virtual real_type l1_norm() = 0; - - /** - * Return the l2 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() = 0; - - /** - * 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() = 0; - - /** - * 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). - */ - virtual Number add_and_dot(const Number a, - const VectorSpaceVector &V, - const VectorSpaceVector &W) = 0; - - /** - * 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 = 0; - - /** - * 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 procesors if this is a distributed 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 const dealii::IndexSet &locally_owned_elements() const = 0; - - /** - * 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 = 0; - - /** - * Return the memory consumption of this class in bytes. - */ - virtual std::size_t memory_consumption() const = 0; - }; - /*@}*/ -} - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/source/lac/CMakeLists.txt b/source/lac/CMakeLists.txt index 648edcfc28..3290ed89a5 100644 --- a/source/lac/CMakeLists.txt +++ b/source/lac/CMakeLists.txt @@ -27,14 +27,12 @@ SET(_src constraint_matrix.cc full_matrix.cc lapack_full_matrix.cc - la_vector.cc matrix_lib.cc matrix_out.cc parallel_vector.cc precondition_block.cc precondition_block_ez.cc relaxation_block.cc - read_write_vector.cc solver.cc solver_control.cc sparse_decomposition.cc @@ -61,11 +59,9 @@ SET(_inst constraint_matrix.inst.in full_matrix.inst.in lapack_full_matrix.inst.in - la_vector.inst.in parallel_vector.inst.in precondition_block.inst.in relaxation_block.inst.in - read_write_vector.inst.in solver.inst.in sparse_matrix_ez.inst.in sparse_matrix.inst.in diff --git a/source/lac/la_vector.cc b/source/lac/la_vector.cc deleted file mode 100644 index ea6500ae7f..0000000000 --- a/source/lac/la_vector.cc +++ /dev/null @@ -1,25 +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 - -DEAL_II_NAMESPACE_OPEN - -namespace LinearAlgebra -{ -#include "la_vector.inst" -} - -DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/la_vector.inst.in b/source/lac/la_vector.inst.in deleted file mode 100644 index 77324a3d57..0000000000 --- a/source/lac/la_vector.inst.in +++ /dev/null @@ -1,29 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 1999 - 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. -// -// --------------------------------------------------------------------- - - - -for (SCALAR : REAL_SCALARS) - { - template class Vector; - } - - - -for (SCALAR : COMPLEX_SCALARS) - { - template class Vector; - } - diff --git a/source/lac/read_write_vector.cc b/source/lac/read_write_vector.cc deleted file mode 100644 index d92a39f040..0000000000 --- a/source/lac/read_write_vector.cc +++ /dev/null @@ -1,23 +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 -#include - -DEAL_II_NAMESPACE_OPEN - -#include "read_write_vector.inst" - -DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/read_write_vector.inst.in b/source/lac/read_write_vector.inst.in deleted file mode 100644 index dd83ecdbca..0000000000 --- a/source/lac/read_write_vector.inst.in +++ /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. -// -// --------------------------------------------------------------------- - - - -for (SCALAR : REAL_SCALARS) -{ - namespace LinearAlgebra - \{ - template class ReadWriteVector; - \} -} - -for (S1, S2 : REAL_SCALARS) -{ - namespace LinearAlgebra - \{ - template void ReadWriteVector::reinit (const ReadWriteVector&, - const bool); - \} -} - - -for (SCALAR : COMPLEX_SCALARS) -{ - namespace LinearAlgebra - \{ - template class ReadWriteVector; - \} -} - -for (S1, S2 : COMPLEX_SCALARS) -{ - namespace LinearAlgebra - \{ - template void ReadWriteVector::reinit (const ReadWriteVector&, - const bool); - \} -} - diff --git a/tests/lac/la_vector_accumulation_01.cc b/tests/lac/la_vector_accumulation_01.cc deleted file mode 100644 index c8b2b9d179..0000000000 --- a/tests/lac/la_vector_accumulation_01.cc +++ /dev/null @@ -1,66 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2012 - 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. -// -// --------------------------------------------------------------------- - - -// check that the l2 norm is exactly the same for many runs on random vector -// size with random vector entries. this is to ensure that the result is -// reproducible also when parallel evaluation is done - -#include "../tests.h" -#include -#include -#include -#include -#include - - - - -template -void check_norms () -{ - for (unsigned int test=0; test<20; ++test) - { - const unsigned int size = Testing::rand() % 100000; - LinearAlgebra::Vector vec (size); - for (unsigned int i=0; i(Testing::rand())/static_cast(RAND_MAX); - const typename LinearAlgebra::ReadWriteVector::real_type norm = vec.l2_norm(); - for (unsigned int i=0; i<30; ++i) - AssertThrow (vec.l2_norm() == norm, ExcInternalError()); - - LinearAlgebra::Vector vec2 (vec); - for (unsigned int i=0; i<10; ++i) - AssertThrow (vec2.l2_norm() == norm, ExcInternalError()); - } -} - - -int main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(2); - deallog.attach(logfile); - deallog.threshold_double(1.e-10); - - check_norms(); - check_norms(); - check_norms(); - check_norms >(); - deallog << "OK" << std::endl; -} - - diff --git a/tests/lac/la_vector_accumulation_01.output b/tests/lac/la_vector_accumulation_01.output deleted file mode 100644 index 0fd8fc12f0..0000000000 --- a/tests/lac/la_vector_accumulation_01.output +++ /dev/null @@ -1,2 +0,0 @@ - -DEAL::OK diff --git a/tests/lac/la_vector_add_and_dot.cc b/tests/lac/la_vector_add_and_dot.cc deleted file mode 100644 index d4f7524acc..0000000000 --- a/tests/lac/la_vector_add_and_dot.cc +++ /dev/null @@ -1,91 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2012 - 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. -// -// --------------------------------------------------------------------- - - -// check that LinearAlgebra::Vector::add_and_dot works correctly - -#include "../tests.h" -#include -#include -#include -#include -#include - - - - -template -void check () -{ - for (unsigned int test=0; test<5; ++test) - { - const unsigned int size = 17 + test*1101; - LinearAlgebra::Vector v1 (size), v2(size), v3(size), check(size); - // Check that the assignment works - v1 = 0.; - for (unsigned int i=0; i::value) - { - deallog << "Vector add reference: "; - for (unsigned int i=0; i(size) - << ", is " << prod_check/static_cast(size) - << std::endl; - } -} - - -int main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(2); - deallog.attach(logfile); - deallog.threshold_double(1.e-10); - - check(); - check(); - check(); - check >(); - deallog << "OK" << std::endl; -} - - diff --git a/tests/lac/la_vector_add_and_dot.output b/tests/lac/la_vector_add_and_dot.output deleted file mode 100644 index 0c6f95da73..0000000000 --- a/tests/lac/la_vector_add_and_dot.output +++ /dev/null @@ -1,25 +0,0 @@ - -DEAL::Add and dot should be 0.30, is 0.30 -DEAL::Add and dot should be 13.40, is 13.40 -DEAL::Add and dot should be 26.50, is 26.50 -DEAL::Add and dot should be 39.61, is 39.61 -DEAL::Add and dot should be 52.71, is 52.71 -DEAL::Vector add reference: 0.03 0.03 0.04 0.05 0.06 0.06 0.07 0.08 0.09 0.09 0.10 0.11 0.12 0.12 0.13 0.14 0.15 -DEAL::Vector check reference: 0.03 0.03 0.04 0.05 0.06 0.06 0.07 0.08 0.09 0.09 0.10 0.11 0.12 0.12 0.13 0.14 0.15 -DEAL::Vector add constant: 1.03 1.03 1.04 1.05 1.06 1.06 1.07 1.08 1.09 1.09 1.10 1.11 1.12 1.12 1.13 1.14 1.15 -DEAL::Add and dot should be 0.30, is 0.30 -DEAL::Add and dot should be 13.40, is 13.40 -DEAL::Add and dot should be 26.50, is 26.50 -DEAL::Add and dot should be 39.61, is 39.61 -DEAL::Add and dot should be 52.71, is 52.71 -DEAL::Add and dot should be 0.30, is 0.30 -DEAL::Add and dot should be 13.40, is 13.40 -DEAL::Add and dot should be 26.50, is 26.50 -DEAL::Add and dot should be 39.61, is 39.61 -DEAL::Add and dot should be 52.71, is 52.71 -DEAL::Add and dot should be (0.30,0.00), is (0.30,0.00) -DEAL::Add and dot should be (13.40,0.00), is (13.40,0.00) -DEAL::Add and dot should be (26.50,0.00), is (26.50,0.00) -DEAL::Add and dot should be (39.61,0.00), is (39.61,0.00) -DEAL::Add and dot should be (52.71,0.00), is (52.71,0.00) -DEAL::OK diff --git a/tests/lac/la_vector_large_numbers.cc b/tests/lac/la_vector_large_numbers.cc deleted file mode 100644 index 4f1d95229f..0000000000 --- a/tests/lac/la_vector_large_numbers.cc +++ /dev/null @@ -1,80 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2012 - 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 "../tests.h" -#include -#include -#include -#include -#include -#include - - - -void check_large_numbers() -{ - LinearAlgebra::Vector v(10); - v(0) = 1e13; - v(1) = 1e19; - v(2) = 1e21; - v(4) = 1.; - v(7) = 1e20; - - const double correct = std::sqrt(1e26 + 1e38 + 1e42 + 1e40 + 1.); - AssertThrow (std::abs(v.l2_norm() - correct) < 1e-6*correct, ExcInternalError()); - - for (unsigned int i=0; i w(7); - w(0) = 1e232; - w(6) = 1e231; - w(3) = 1e200; - const double correct3 = std::sqrt(100. + 1.) * 1e231; - AssertThrow (std::abs(w.l2_norm() - correct3) < 1e-13*correct3, ExcInternalError()); - - for (unsigned int i=0; i -#include -#include -#include -#include - - - - -template -void check_norms () -{ - const number acc = 1e1*std::numeric_limits::epsilon(); - unsigned int skip = 73; - for (unsigned int size=1; size<200000; size+=skip) - { - // test correctness - if (size > 10000) - skip += 17; - LinearAlgebra::Vector vec(size); - for (unsigned int i=0; i -void check_complex_norms () -{ - const number acc = 1e2*std::numeric_limits::epsilon(); - unsigned int skip = 73; - for (unsigned int size=1; size<100000; size+=skip) - { - // test correctness - if (size > 10000) - skip += 17; - LinearAlgebra::Vector > vec(size); - long double sum = 0.; - for (unsigned int i=0; i(i+1,i+2); - sum += std::sqrt((long double)(i+1)*(i+1) + (long double)(i+2)*(i+2)); - } - - const number l1_norm = vec.l1_norm(); - AssertThrow (std::abs(l1_norm-sum) < acc*sum, - ExcInternalError()); - - // test accuracy of summation - const std::complex value (3.14159265358979323846, 0.1); - for (unsigned int i=0; i(value); - const number l1_norma = vec.l1_norm(); - AssertThrow (std::abs(l1_norma-std::abs(value)*size) < - acc*size*std::abs(value), - ExcInternalError()); - const number l2_norma = vec.l2_norm(); - AssertThrow (std::abs(l2_norma-std::abs(value)*std::sqrt((number)size)) < - acc*std::sqrt((number)size)*std::abs(value), - ExcInternalError()); - } -} - - -int main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(2); - deallog.attach(logfile); - deallog.threshold_double(1.e-10); - - check_norms(); - check_norms(); - check_norms(); - check_complex_norms(); - check_complex_norms(); - deallog << "OK" << std::endl; -} - - diff --git a/tests/lac/la_vector_norms.output b/tests/lac/la_vector_norms.output deleted file mode 100644 index 0fd8fc12f0..0000000000 --- a/tests/lac/la_vector_norms.output +++ /dev/null @@ -1,2 +0,0 @@ - -DEAL::OK diff --git a/tests/lac/la_vector_output.cc b/tests/lac/la_vector_output.cc deleted file mode 100644 index 62c91f8445..0000000000 --- a/tests/lac/la_vector_output.cc +++ /dev/null @@ -1,87 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2012 - 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 "../tests.h" -#include -#include -#include -#include -#include -#include - - -void test() -{ - const char *filename = "test.txt"; - const unsigned int size(10); - LinearAlgebra::Vector vec(size); - for (unsigned int i=0; i> vec; - } - // Check the vector - for (unsigned int i=0; i -#include -#include -#include -#include - - - - -const unsigned int N=10; -unsigned int check_point = 0; - - - - -template -void print (const LinearAlgebra::Vector &v) -{ - for (unsigned int i=0; i -void check_vectors (LinearAlgebra::Vector &d1, LinearAlgebra::Vector &d2) -{ - deallog << "Fill & Swap" << std::endl; - LinearAlgebra::Vector d3(d1.size()); - print (d3); - - for (unsigned int i=0; i d4 (&array[0], &array[4]); - print (d4); - - deallog << "Extract number" << std::endl; - // Each line should contain two equal numbers - double sum = 0.; - for (unsigned int i=0; isum) sum = t; - } - deallog << d3.linfty_norm() << '\t' << sum << std::endl; - - deallog << "add & sub" << std::endl; - - d1 += d2; - print (d1); - - d2 -= d1; - print (d2); - - d1.add (2, d3); - print (d1); - - d1.add (2., d2, .5, d3); - print (d1); - - deallog << "sadd & scale" << std::endl; - - d2.sadd (2., 1., d1); - print (d2); - - d2.sadd (2., .5, d1); - print (d2); - - d1 *= 4.; - print (d1); - - deallog << "equ" << std::endl; - - d2.equ (.25, d1); - print (d2); -} - - -int main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(2); - deallog.attach(logfile); - deallog.threshold_double(1.e-10); - - LinearAlgebra::Vector d1(N), d2(N); - LinearAlgebra::Vector f1(N), f2(N); - LinearAlgebra::Vector l1(N), l2(N); - - // cross-tests with double/float - // vectors at the same time are not - // supported at present, - // unfortunately, as many functions - // don't accept other data types as - // arguments - check_vectors (d1, d2); - check_vectors (f1, f2); - check_vectors (l1, l2); -} - - diff --git a/tests/lac/la_vector_vector.output b/tests/lac/la_vector_vector.output deleted file mode 100644 index 46bd0872f1..0000000000 --- a/tests/lac/la_vector_vector.output +++ /dev/null @@ -1,76 +0,0 @@ - -DEAL::Fill & Swap -DEAL::0 0 0 0 0 0 0 0 0 0 -DEAL::0 2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 -DEAL::0 2.00 8.00 18.00 32.00 50.00 72.00 98.00 128.00 162.00 -DEAL::2.00 1.50 1.00 0.50 0 -0.50 -1.00 -1.50 -2.00 -2.50 -DEAL::0 2.00 8.00 18.00 32.00 50.00 72.00 98.00 128.00 162.00 -DEAL::0 2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 -DEAL::2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 -DEAL::0 1.10 2.20 3.30 -DEAL::Extract number -DEAL::-105.00 -105.00 -DEAL::33.76 33.76 -DEAL::12.50 12.50 -DEAL::2.50 2.50 -DEAL::add & sub -DEAL::2.50 4.50 6.50 8.50 10.50 12.50 14.50 16.50 18.50 20.50 -DEAL::-2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -DEAL::6.50 7.50 8.50 9.50 10.50 11.50 12.50 13.50 14.50 15.50 -DEAL::2.50 3.25 4.00 4.75 5.50 6.25 7.00 7.75 8.50 9.25 -DEAL::sadd & scale -DEAL::-2.50 -1.75 -1.00 -0.25 0.50 1.25 2.00 2.75 3.50 4.25 -DEAL::-3.75 -1.88 0 1.88 3.75 5.62 7.50 9.38 11.25 13.12 -DEAL::10.00 13.00 16.00 19.00 22.00 25.00 28.00 31.00 34.00 37.00 -DEAL::equ -DEAL::2.50 3.25 4.00 4.75 5.50 6.25 7.00 7.75 8.50 9.25 -DEAL::Fill & Swap -DEAL::0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 -DEAL::0.00 2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 -DEAL::0.00 2.00 8.00 18.00 32.00 50.00 72.00 98.00 128.00 162.00 -DEAL::2.00 1.50 1.00 0.50 0.00 -0.50 -1.00 -1.50 -2.00 -2.50 -DEAL::0.00 2.00 8.00 18.00 32.00 50.00 72.00 98.00 128.00 162.00 -DEAL::0.00 2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 -DEAL::2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 -DEAL::0.00 1.10 2.20 3.30 -DEAL::Extract number -DEAL::-105.00 -105.00 -DEAL::33.76 33.76 -DEAL::12.50 12.50 -DEAL::2.50 2.50 -DEAL::add & sub -DEAL::2.50 4.50 6.50 8.50 10.50 12.50 14.50 16.50 18.50 20.50 -DEAL::-2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -DEAL::6.50 7.50 8.50 9.50 10.50 11.50 12.50 13.50 14.50 15.50 -DEAL::2.50 3.25 4.00 4.75 5.50 6.25 7.00 7.75 8.50 9.25 -DEAL::sadd & scale -DEAL::-2.50 -1.75 -1.00 -0.25 0.50 1.25 2.00 2.75 3.50 4.25 -DEAL::-3.75 -1.88 0.00 1.88 3.75 5.62 7.50 9.38 11.25 13.12 -DEAL::10.00 13.00 16.00 19.00 22.00 25.00 28.00 31.00 34.00 37.00 -DEAL::equ -DEAL::2.50 3.25 4.00 4.75 5.50 6.25 7.00 7.75 8.50 9.25 -DEAL::Fill & Swap -DEAL::0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 -DEAL::0.00 2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 -DEAL::0.00 2.00 8.00 18.00 32.00 50.00 72.00 98.00 128.00 162.00 -DEAL::2.00 1.50 1.00 0.50 0.00 -0.50 -1.00 -1.50 -2.00 -2.50 -DEAL::0.00 2.00 8.00 18.00 32.00 50.00 72.00 98.00 128.00 162.00 -DEAL::0.00 2.00 4.00 6.00 8.00 10.00 12.00 14.00 16.00 18.00 -DEAL::2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 -DEAL::0.00 1.10 2.20 3.30 -DEAL::Extract number -DEAL::-105.00 -105.00 -DEAL::33.76 33.76 -DEAL::12.50 12.50 -DEAL::2.50 2.50 -DEAL::add & sub -DEAL::2.50 4.50 6.50 8.50 10.50 12.50 14.50 16.50 18.50 20.50 -DEAL::-2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -2.50 -DEAL::6.50 7.50 8.50 9.50 10.50 11.50 12.50 13.50 14.50 15.50 -DEAL::2.50 3.25 4.00 4.75 5.50 6.25 7.00 7.75 8.50 9.25 -DEAL::sadd & scale -DEAL::-2.50 -1.75 -1.00 -0.25 0.50 1.25 2.00 2.75 3.50 4.25 -DEAL::-3.75 -1.88 0.00 1.88 3.75 5.62 7.50 9.38 11.25 13.12 -DEAL::10.00 13.00 16.00 19.00 22.00 25.00 28.00 31.00 34.00 37.00 -DEAL::equ -DEAL::2.50 3.25 4.00 4.75 5.50 6.25 7.00 7.75 8.50 9.25 diff --git a/tests/lac/readwritevector_add.cc b/tests/lac/readwritevector_add.cc deleted file mode 100644 index be52267a3d..0000000000 --- a/tests/lac/readwritevector_add.cc +++ /dev/null @@ -1,83 +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. -// -// --------------------------------------------------------------------- - -// Check reinit and add - -#include "../tests.h" -#include -#include -#include - -#include -#include - -void test() -{ - unsigned int size(10); - std::vector indices(size); - std::vector float_values(size); - std::vector double_values(size); - for (unsigned int i=0; i float_vector; - LinearAlgebra::ReadWriteVector double_vector_1; - LinearAlgebra::ReadWriteVector double_vector_2; - IndexSet is(50); - is.add_range(0,10); - is.add_index(46); - is.add_range(11,25); - is.compress(); - - float_vector.reinit(size); - double_vector_1.reinit(float_vector); - double_vector_2.reinit(is); - - IndexSet is_copy(double_vector_2.get_stored_elements()); - is.print(std::cout); - is_copy.print(std::cout); - - float_vector.add(indices,float_values); - double_vector_1.add(indices,float_vector); - double_vector_2.add(size,&indices[0],&double_values[0]); - double_vector_1.swap(double_vector_2); - LinearAlgebra::ReadWriteVector::iterator val_1 = double_vector_1.begin(); - LinearAlgebra::ReadWriteVector::iterator end_1 = double_vector_1.end(); - LinearAlgebra::ReadWriteVector::iterator val_2 = double_vector_2.begin(); - LinearAlgebra::ReadWriteVector::iterator end_2 = double_vector_2.end(); - deallog << "double_vector_1" << std::endl; - for (; val_1 -#include -#include - -#include - - -void test() -{ - unsigned int double_size = 2; - unsigned int float_size = 10; - IndexSet is(50); - is.add_range(0,2); - is.add_index(46); - is.add_range(10,15); - LinearAlgebra::ReadWriteVector double_vector(is); - LinearAlgebra::ReadWriteVector float_vector(float_size); - deallog << "double_size " << double_vector.n_elements() < double_vector2(double_vector); - double_vector2.print(deallog.get_file_stream()); - - for (unsigned int i=0; i -#include -#include - -#include - - -void test() -{ - unsigned int double_size = 2; - unsigned int float_size = 10; - IndexSet is(50); - is.add_range(0,2); - is.add_index(46); - is.add_range(10,15); - LinearAlgebra::ReadWriteVector double_vector(is); - LinearAlgebra::ReadWriteVector float_vector(float_size); - deallog << "double_size " << double_vector.n_elements() < double_vector2(double_vector); - double_vector2.print(deallog.get_file_stream()); - - for (unsigned int i=0; i -#include -#include -#include -#include - -#include -#include - -void test() -{ - IndexSet is(8); - unsigned int rank = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); - if (rank==0) - is.add_range(0,4); - if (rank==1) - is.add_range(4,8); - is.compress(); - TrilinosWrappers::MPI::Vector tril_vector(is); - Vector tmp(8); - for (unsigned int i=0; i<8; ++i) - tmp[i] = i; - tril_vector = tmp; - - tril_vector.compress(VectorOperation::insert); - - MPI_Barrier(MPI_COMM_WORLD); - - IndexSet readwrite_is(8); - if (rank==0) - { - readwrite_is.add_range(0,2); - readwrite_is.add_range(6,8); - } - if (rank==1) - readwrite_is.add_range(2,6); - readwrite_is.compress(); - LinearAlgebra::ReadWriteVector readwrite(readwrite_is); - readwrite = tril_vector; - if (rank==0) - { - std::vector comp(4); - comp[0] = 0.; - comp[1] = 1.; - comp[2] = 6.; - comp[3] = 7.; - for (unsigned int i=0; i<4; ++i) - AssertThrow(readwrite.local_element(i)==comp[i], - ExcMessage("Element not copied correctly")); - } - MPI_Barrier(MPI_COMM_WORLD); - if (rank==1) - { - std::vector comp(4); - comp[0] = 2.; - comp[1] = 3.; - comp[2] = 4.; - comp[3] = 5.; - for (unsigned int i=0; i<4; ++i) - AssertThrow(readwrite.local_element(i)==comp[i], - ExcMessage("Element not copied correctly")); - } - deallog << "OK" <