From c80ebe09b3b3430d4efe90b170fdd525b436d077 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Wed, 11 Nov 2015 08:39:40 -0600 Subject: [PATCH] Add read/write operations to LinearAlgebra::Vector. --- include/deal.II/base/index_set.h | 1 + include/deal.II/lac/la_vector.h | 118 ++++++++++++++++++ include/deal.II/lac/la_vector.templates.h | 63 ++++++++++ .../deal.II/lac/read_write_vector.templates.h | 1 + include/deal.II/lac/vector_space_vector.h | 5 + tests/lac/la_vector_add_and_dot.cc | 15 ++- tests/lac/la_vector_add_and_dot.output | 1 + tests/lac/la_vector_norms.cc | 4 +- tests/lac/la_vector_output.cc | 88 +++++++++++++ tests/lac/la_vector_output.output | 2 + 10 files changed, 293 insertions(+), 5 deletions(-) create mode 100644 tests/lac/la_vector_output.cc create mode 100644 tests/lac/la_vector_output.output diff --git a/include/deal.II/base/index_set.h b/include/deal.II/base/index_set.h index b0dc5a4016..0ce6ef03b1 100644 --- a/include/deal.II/base/index_set.h +++ b/include/deal.II/base/index_set.h @@ -19,6 +19,7 @@ #include #include #include +#include #include #include diff --git a/include/deal.II/lac/la_vector.h b/include/deal.II/lac/la_vector.h index 0e0ea0acb7..519fe6114f 100644 --- a/include/deal.II/lac/la_vector.h +++ b/include/deal.II/lac/la_vector.h @@ -23,6 +23,8 @@ #include #include #include +#include +#include #include #include @@ -88,6 +90,23 @@ namespace LinearAlgebra */ 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. */ @@ -113,6 +132,11 @@ namespace LinearAlgebra */ 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. */ @@ -204,6 +228,25 @@ namespace LinearAlgebra 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. @@ -233,6 +276,15 @@ namespace LinearAlgebra 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; }; /*@}*/ @@ -279,6 +331,56 @@ namespace LinearAlgebra + 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 @@ -309,6 +411,22 @@ namespace LinearAlgebra + 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 diff --git a/include/deal.II/lac/la_vector.templates.h b/include/deal.II/lac/la_vector.templates.h index f436504d48..92f8d0bbe7 100644 --- a/include/deal.II/lac/la_vector.templates.h +++ b/include/deal.II/lac/la_vector.templates.h @@ -17,6 +17,8 @@ #define dealii__la_vector_templates_h #include +#include +#include DEAL_II_NAMESPACE_OPEN @@ -346,6 +348,67 @@ namespace LinearAlgebra 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 diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index 35e2c301d6..9e334de71d 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -61,6 +61,7 @@ namespace LinearAlgebra resize_val(size); stored_elements = complete_index_set(size); + stored_elements.compress(); // set entries to zero if so requested if (omit_zeroing_entries == false) diff --git a/include/deal.II/lac/vector_space_vector.h b/include/deal.II/lac/vector_space_vector.h index 31ba1b1c72..1f09437bac 100644 --- a/include/deal.II/lac/vector_space_vector.h +++ b/include/deal.II/lac/vector_space_vector.h @@ -77,6 +77,11 @@ namespace LinearAlgebra */ 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. */ diff --git a/tests/lac/la_vector_add_and_dot.cc b/tests/lac/la_vector_add_and_dot.cc index 4441163888..07b68029d4 100644 --- a/tests/lac/la_vector_add_and_dot.cc +++ b/tests/lac/la_vector_add_and_dot.cc @@ -33,11 +33,13 @@ void check () { 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(size) diff --git a/tests/lac/la_vector_add_and_dot.output b/tests/lac/la_vector_add_and_dot.output index 329e9f8d54..0c6f95da73 100644 --- a/tests/lac/la_vector_add_and_dot.output +++ b/tests/lac/la_vector_add_and_dot.output @@ -6,6 +6,7 @@ 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 diff --git a/tests/lac/la_vector_norms.cc b/tests/lac/la_vector_norms.cc index 7f2faad389..e451006c2c 100644 --- a/tests/lac/la_vector_norms.cc +++ b/tests/lac/la_vector_norms.cc @@ -40,7 +40,7 @@ void check_norms () skip += 17; LinearAlgebra::Vector vec(size); for (unsigned int i=0; i void check_complex_norms () { - const number acc = 1e1*std::numeric_limits::epsilon(); + const number acc = 1e2*std::numeric_limits::epsilon(); unsigned int skip = 73; for (unsigned int size=1; size<100000; size+=skip) { diff --git a/tests/lac/la_vector_output.cc b/tests/lac/la_vector_output.cc new file mode 100644 index 0000000000..0aeb7d89d3 --- /dev/null +++ b/tests/lac/la_vector_output.cc @@ -0,0 +1,88 @@ +// --------------------------------------------------------------------- +// +// 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() +{ + std::string 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