]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement a random-access iterator for TransposeTable. 5959/head
authorDavid Wells <wellsd2@rpi.edu>
Sun, 18 Mar 2018 22:41:09 +0000 (18:41 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Thu, 22 Mar 2018 19:48:51 +0000 (15:48 -0400)
doc/news/changes/minor/20180225DavidWells [new file with mode: 0644]
include/deal.II/base/table.h
tests/base/transpose_table_iterators_0.cc [new file with mode: 0644]
tests/base/transpose_table_iterators_0.output [new file with mode: 0644]
tests/base/transpose_table_iterators_1.cc [new file with mode: 0644]
tests/base/transpose_table_iterators_1.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20180225DavidWells b/doc/news/changes/minor/20180225DavidWells
new file mode 100644 (file)
index 0000000..503a5fe
--- /dev/null
@@ -0,0 +1,4 @@
+New: TransposeTable (the base class of LAPACKFullMatrix and ScaLAPACKMatrix) now
+has a random access iterator implementation similar to the one provided by
+SparseMatrix.
+<br> (David Wells, 2018/02/25)
index f4c8b05baeabb42ba58250f252a8b970efed8094..cf263b1059c027241dab20191b1a902f325529f8 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/base/table_indices.h>
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/aligned_vector.h>
+#include <deal.II/base/linear_index_iterator.h>
 
 #include <cstddef>
 #include <algorithm>
@@ -31,6 +32,7 @@ DEAL_II_NAMESPACE_OPEN
 // forward declaration
 template <int N, typename T> class TableBase;
 template <int N, typename T> class Table;
+template <typename T> class TransposeTable;
 template <typename T> class Table<1,T>;
 template <typename T> class Table<2,T>;
 template <typename T> class Table<3,T>;
@@ -1499,6 +1501,202 @@ public:
 
 
 
+/**
+ * A namespace for iterators for TransposeTable. TransposeTable is unique in
+ * that it stores entries in column-major order.
+ */
+namespace TransposeTableIterators
+{
+  // Forward declaration of the iterator class.
+  template <typename T, bool Constness>
+  class Iterator;
+
+  /**
+   * Accessor class template. This will be partially specialized for both
+   * values of <code>Constness</code>.
+   */
+  template <typename T, bool Constness>
+  class Accessor;
+
+  /**
+   * Accessor base class. Stores a const pointer when <code>Constness ==
+   * true</code> and a non-const pointer when <code>Constness == false</code>.
+   */
+  template <typename T, bool Constness>
+  class AccessorBase
+  {
+  public:
+    /**
+     * Type of the stored pointer to the TransposeTable.
+     */
+    typedef typename std::conditional<Constness, const TransposeTable<T>*, TransposeTable<T>*>::type
+    container_pointer_type;
+
+    /**
+     * Default constructor.
+     */
+    AccessorBase();
+
+    /**
+     * Constructor setting up the end iterator.
+     */
+    AccessorBase(const container_pointer_type table);
+
+    /**
+     * Copy constructor from a non-const Accessor.
+     */
+    AccessorBase(const AccessorBase<T, false> &);
+
+    /**
+     * Constructor taking an array index. Note that TransposeTable stores data
+     * in column-major order.
+     */
+    AccessorBase(const container_pointer_type table,
+                 const std::ptrdiff_t linear_index);
+
+    /**
+     * Get a constant reference to the value of the element represented by
+     * this object.
+     */
+    const T &value() const;
+
+    /**
+     * Numerical type of the row and column indices of the TransposeTable.
+     */
+    typedef typename TransposeTable<T>::size_type size_type;
+
+    /**
+     * Access the row of the current entry.
+     */
+    size_type row() const;
+
+    /**
+     * Access the column of the current entry.
+     */
+    size_type column() const;
+
+  protected:
+    /**
+     * Pointer to the TransposeTable.
+     */
+    container_pointer_type container;
+
+    /**
+     * Current index.
+     */
+    std::ptrdiff_t linear_index;
+
+    /**
+     * Make the const version a friend for copying.
+     */
+    friend class AccessorBase<T, true>;
+
+    /**
+     * Make the underlying iterator class a friend.
+     */
+    friend class LinearIndexIterator<Iterator<T, Constness>, Accessor<T, Constness>>;
+  };
+
+  /**
+   * Constant Accessor. Same as the base class.
+   */
+  template <typename T>
+  class Accessor<T, true> : public AccessorBase<T, true>
+  {
+  public:
+    /**
+     * Inherit the base class constructors.
+     */
+    using AccessorBase<T, true>::AccessorBase;
+  };
+
+  /**
+   * Accessor class for non-constant tables. Defines a few extra operations
+   * for modifying values.
+   */
+  template <typename T>
+  class Accessor<T, false> : public AccessorBase<T, false>
+  {
+  public:
+    /**
+     * Inherit the base class constructors.
+     */
+    using AccessorBase<T, false>::AccessorBase;
+
+    /**
+     * Assignment operator. This assigns a new value to the table entry at the
+     * current row and column coordinates.
+     */
+    const Accessor<T, false> &operator = (const T &) const;
+
+    /**
+     * Move assignment operator. This assigns a new value to the table entry at the
+     * current row and column coordinates.
+     */
+    const Accessor<T, false> &operator = (T&&) const;
+
+    /**
+     * Since we overload value() we have to explicitly use the base class
+     * version.
+     */
+    using AccessorBase<T, false>::value;
+
+    /**
+     * Get a reference to the value of the element represented by
+     * this object.
+     */
+    T &value() const;
+  };
+
+  /**
+   * Iterator class for TransposeTable.
+   */
+  template <typename T, bool Constness>
+  class Iterator : public LinearIndexIterator<Iterator<T, Constness>, Accessor<T, Constness>>
+  {
+  public:
+    /**
+     * Size type used by TransposeTable.
+     */
+    typedef typename TransposeTable<T>::size_type size_type;
+
+    /**
+     * Type of the stored pointer to the TransposeTable.
+     */
+    typedef typename std::conditional<Constness, const TransposeTable<T>*, TransposeTable<T>*>::type
+    container_pointer_type;
+
+    /**
+     * Constructor from an accessor.
+     */
+    Iterator (const Accessor<T, Constness> &accessor);
+
+    /**
+     * Constructor. Create the end iterator for a table.
+     */
+    Iterator (const container_pointer_type object);
+
+    /**
+     * Constructor for a particular table entry.
+     */
+    Iterator (const container_pointer_type object,
+              const size_type row,
+              const size_type column);
+
+    /**
+     * Copy constructor from a non-const iterator.
+     */
+    Iterator(const Iterator<T, false> &i);
+
+    /**
+     * Constructor for an entry with a particular linear index.
+     */
+    Iterator (const container_pointer_type container,
+              const std::ptrdiff_t         linear_index);
+  };
+}
+
+
 /**
  * A class representing a transpose two-dimensional table, i.e. a matrix of
  * objects (not necessarily only numbers) in column first numbering (FORTRAN
@@ -1536,6 +1734,17 @@ public:
    */
   typedef typename AlignedVector<T>::const_reference const_reference;
 
+  /**
+   * Typedef for a constant iterator that traverses the table in column-major
+   * order.
+   */
+  typedef TransposeTableIterators::Iterator<T, true> const_iterator;
+
+  /**
+   * Typedef for an iterator that traverses the table in column-major order.
+   */
+  typedef TransposeTableIterators::Iterator<T, false> iterator;
+
   /**
    * Default constructor. Set all dimensions to zero.
    */
@@ -1586,6 +1795,26 @@ public:
    */
   size_type n_cols () const;
 
+  /**
+   * Return an iterator pointing to the first entry.
+   */
+  iterator begin ();
+
+  /**
+   * Return a constant iterator pointing to the first entry.
+   */
+  const_iterator begin () const;
+
+  /**
+   * Return an iterator pointing to one past the last entry.
+   */
+  iterator end ();
+
+  /**
+   * Return a constant iterator pointing to one past the last entry.
+   */
+  const_iterator end () const;
+
 protected:
   /**
    * Return a read-write reference to the element <tt>(i,j)</tt>.
@@ -1616,6 +1845,19 @@ protected:
    */
   const_reference el (const size_type i,
                       const size_type j) const;
+
+  /**
+   * Make the AccessorBase class a friend so that it may directly index into
+   * the array.
+   */
+  friend class TransposeTableIterators::AccessorBase<T, true>;
+  friend class TransposeTableIterators::AccessorBase<T, false>;
+
+  /**
+   * Make the mutable accessor class a friend so that we can write to array
+   * entries with iterators.
+   */
+  friend class TransposeTableIterators::Accessor<T, false>;
 };
 
 
@@ -2382,6 +2624,178 @@ Table<2,T>::n_cols () const
 
 
 
+//---------------------------------------------------------------------------
+namespace TransposeTableIterators
+{
+  template <typename T, bool Constness>
+  inline
+  AccessorBase<T, Constness>::AccessorBase()
+    :
+    container(nullptr),
+    linear_index(std::numeric_limits<decltype(linear_index)>::max())
+  {}
+
+
+
+  template <typename T, bool Constness>
+  inline
+  AccessorBase<T, Constness>::AccessorBase(const container_pointer_type table)
+    :
+    container(table),
+    linear_index(container->values.size())
+  {}
+
+
+
+  template <typename T, bool Constness>
+  inline
+  AccessorBase<T, Constness>::AccessorBase(const AccessorBase<T, false> &a)
+    :
+    container(a.container),
+    linear_index(a.linear_index)
+  {}
+
+
+
+  template <typename T, bool Constness>
+  inline
+  AccessorBase<T, Constness>::AccessorBase(const container_pointer_type table,
+                                           const std::ptrdiff_t         index)
+    :
+    container(table),
+    linear_index(index)
+  {
+    Assert(0 <= linear_index && linear_index < container->values.size() + 1,
+           ExcMessage("The current iterator points outside of the table and is "
+                      "not the end iterator."));
+  }
+
+
+
+  template <typename T, bool Constness>
+  inline
+  const T &
+  AccessorBase<T, Constness>::value() const
+  {
+    Assert(0 <= linear_index && linear_index < container->values.size(),
+           ExcMessage("The current iterator points outside of the table."));
+    return this->container->values[linear_index];
+  }
+
+
+
+  template <typename T, bool Constness>
+  inline
+  typename AccessorBase<T, Constness>::size_type
+  AccessorBase<T, Constness>::AccessorBase::row() const
+  {
+    Assert(!container->empty(), ExcMessage("An empty table has no rows or columns."));
+    const auto row_n = linear_index % container->n_rows();
+    Assert(0 <= row_n && row_n < container->n_rows(),
+           ExcMessage("The current iterator points outside the table."));
+    return row_n;
+  }
+
+
+
+  template <typename T, bool Constness>
+  inline
+  typename AccessorBase<T, Constness>::size_type
+  AccessorBase<T, Constness>::AccessorBase::column() const
+  {
+    Assert(!container->empty(), ExcMessage("An empty table has no rows or columns."));
+    const auto column_n = linear_index/container->n_rows();
+    Assert(0 <= column_n && column_n < container->n_cols(),
+           ExcMessage("The current iterator points outside the table."));
+    return column_n;
+  }
+
+
+  template <typename T>
+  inline
+  const Accessor<T, false> &
+  Accessor<T, false>::operator = (const T &t) const
+  {
+    Assert(0 <= this->linear_index && this->linear_index < this->container->values.size(),
+           ExcMessage("The current iterator points outside of the table."));
+    this->container->values[this->linear_index] = t;
+    return *this;
+  }
+
+
+
+  template <typename T>
+  inline
+  const Accessor<T, false> &
+  Accessor<T, false>::operator = (T &&t) const
+  {
+    Assert(0 <= this->linear_index && this->linear_index < this->container->values.size(),
+           ExcMessage("The current iterator points outside of the table."));
+    this->container->values[this->linear_index] = t;
+    return *this;
+  }
+
+
+
+  template <typename T>
+  inline
+  T &
+  Accessor<T, false>::value() const
+  {
+    Assert(0 <= this->linear_index && this->linear_index < this->container->values.size(),
+           ExcMessage("The current iterator points outside of the table."));
+    return this->container->values[this->linear_index];
+  }
+
+
+
+  template <typename T, bool Constness>
+  Iterator<T, Constness>::Iterator(const Accessor<T, Constness> &a)
+    :
+    LinearIndexIterator<Iterator<T, Constness>, Accessor<T, Constness>>(a)
+  {}
+
+
+
+  template <typename T, bool Constness>
+  Iterator<T, Constness>::Iterator(const container_pointer_type table)
+    :
+    LinearIndexIterator<Iterator<T, Constness>, Accessor<T, Constness>>(Accessor<T, Constness>(table))
+  {}
+
+
+
+  template <typename T, bool Constness>
+  Iterator<T, Constness>::Iterator(const Iterator<T, false> &i)
+    :
+    LinearIndexIterator<Iterator<T, Constness>, Accessor<T, Constness>>(*i)
+  {}
+
+
+
+
+  template <typename T, bool Constness>
+  Iterator<T, Constness>::Iterator(const container_pointer_type table,
+                                   const size_type              row_n,
+                                   const size_type              col_n)
+    :
+    Iterator(table, table->n_rows()*col_n + row_n)
+  {}
+
+
+
+  template <typename T, bool Constness>
+  Iterator<T, Constness>::Iterator
+  (const container_pointer_type table,
+   const std::ptrdiff_t         linear_index)
+    :
+    LinearIndexIterator<Iterator<T, Constness>, Accessor<T, Constness>>
+        (Accessor<T, Constness>(table, linear_index))
+  {}
+}
+
+
+
 //---------------------------------------------------------------------------
 template <typename T>
 inline
@@ -2473,6 +2887,43 @@ TransposeTable<T>::n_cols () const
 
 
 
+template <typename T>
+inline
+typename TransposeTable<T>::iterator
+TransposeTable<T>::begin ()
+{
+  return typename TransposeTable<T>::iterator(this, 0, 0);
+}
+
+
+
+template <typename T>
+inline
+typename TransposeTable<T>::const_iterator
+TransposeTable<T>::begin () const
+{
+  return typename TransposeTable<T>::const_iterator(this, 0, 0);
+}
+
+
+
+template <typename T>
+inline
+typename TransposeTable<T>::iterator
+TransposeTable<T>::end ()
+{
+  return typename TransposeTable<T>::iterator(this);
+}
+
+
+
+template <typename T>
+inline
+typename TransposeTable<T>::const_iterator
+TransposeTable<T>::end () const
+{
+  return typename TransposeTable<T>::const_iterator(this);
+}
 //---------------------------------------------------------------------------
 
 
diff --git a/tests/base/transpose_table_iterators_0.cc b/tests/base/transpose_table_iterators_0.cc
new file mode 100644 (file)
index 0000000..e3ccf2f
--- /dev/null
@@ -0,0 +1,102 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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 some basic behavior of the TransposeTable iterators
+#include "../tests.h"
+
+#include <deal.II/base/table.h>
+
+#include <boost/core/demangle.hpp>
+
+template <typename Container>
+void test()
+{
+  Container table(3, 4);
+
+  auto it = table.begin();
+  auto t2 = it;
+  deallog << std::boolalpha;
+  deallog << boost::core::demangle(typeid(Container).name()) << std::endl;
+  deallog << "++it == it:           " << (++table.begin() == table.begin()) << std::endl;
+  deallog << "++it != it:           " << (++table.begin() != table.begin()) << std::endl;
+  deallog << "it++ == it:           " << (table.begin()++ == table.begin()) << std::endl;
+  deallog << "it++ != it + 1:       " << (table.begin()++ != table.begin() + 1) << std::endl;
+  deallog << "++(++it) == it + 2:   " << (++(++table.begin()) == table.begin() + 2) << std::endl;
+  deallog << "it++ == it + 1:       " << (table.begin()++ == table.begin() + 1) << std::endl;
+  deallog << "--(++it) == it:       " << (--(++table.begin()) == table.begin()) << std::endl;
+  deallog << "(it--)-- == it:       " << ((table.end()--)-- == table.end()) << std::endl;
+  deallog << "end - 5 <= end:       " << ((table.end() - 5) <= table.end()) << std::endl;
+  deallog << "end <= end:           " << (table.end() <= table.end()) << std::endl;
+  deallog << "end >= end:           " << (table.end() >= table.end()) << std::endl;
+  deallog << "begin + 5 >= begin:   " << ((table.begin() + 5) >= table.begin()) << std::endl;
+  deallog << "begin + 5 <= begin:   " << ((table.begin() + 5) <= table.begin()) << std::endl;
+  deallog << "end - begin:          " << table.end() - table.begin() << std::endl;
+  deallog << "begin - end:          " << table.begin() - table.end() << std::endl;
+  deallog << "(begin + 2) - begin:  " << (table.begin() + 2) - table.begin() << std::endl;
+  deallog << "begin + 6 < begin:    " << ((table.begin() + 6) < table.begin()) << std::endl;
+  deallog << "begin < begin + 1:    " << (table.begin() < (table.begin() + 6)) << std::endl;
+  deallog << "end - 5 < begin:      " << ((table.end() - 6) < table.begin()) << std::endl;
+  t2 = it + 5;
+  deallog << "it+5 == (t2 = (it+5)):" << (it + 5 == t2) << std::endl;
+  const auto it2 = table.end() - 5;
+  deallog << "end - 5 position:     " << it2->row() << ", " << it2->column() << std::endl;
+  const auto it3 = table.begin() + 5;
+  deallog << "begin + 5 position:   " << it3->row() << ", " << it3->column() << std::endl;
+}
+
+int main()
+{
+  initlog();
+
+  // Test an empty table
+  {
+    TransposeTable<double> transpose_table;
+    Assert(transpose_table.begin() == transpose_table.end(),
+           ExcMessage("The beginning and end iterators should be equal for an "
+                      "empty table."));
+    TransposeTable<double>::const_iterator begin = transpose_table.begin();
+    TransposeTable<double>::const_iterator end = transpose_table.end();
+    Assert(begin == end, ExcMessage("The beginning and end const iterators should"
+                                    " be equal for an empty table."));
+  }
+
+  // test some things with accessors
+  {
+    TransposeTable<double> table(2, 2);
+
+    TransposeTableIterators::Accessor<double, true> a3(&table, 2);
+    TransposeTableIterators::Accessor<double, false> a4(&table, 2);
+    deallog << "Accessors refer to the same entry: "
+            << (&(TransposeTableIterators::Accessor<double, true>(a4).value()) == &(a3.value()))
+            << std::endl;
+  }
+
+  // test a non-empty rectangular table
+  test<TransposeTable<double>>();
+  test<const TransposeTable<TransposeTable<double>>>();
+
+  {
+    // check that we can convert
+    TransposeTable<double> table(3, 4);
+    TransposeTable<double>::iterator it0 = table.begin();
+    TransposeTable<double>::const_iterator it1 = it0;
+    deallog << "converted iterators are equal: "
+            << (it1 == it0)
+            << std::endl;
+  }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/base/transpose_table_iterators_0.output b/tests/base/transpose_table_iterators_0.output
new file mode 100644 (file)
index 0000000..161b855
--- /dev/null
@@ -0,0 +1,50 @@
+
+DEAL::Accessors refer to the same entry: 1
+DEAL::dealii::TransposeTable<double>
+DEAL::++it == it:           false
+DEAL::++it != it:           true
+DEAL::it++ == it:           true
+DEAL::it++ != it + 1:       true
+DEAL::++(++it) == it + 2:   true
+DEAL::it++ == it + 1:       false
+DEAL::--(++it) == it:       true
+DEAL::(it--)-- == it:       true
+DEAL::end - 5 <= end:       true
+DEAL::end <= end:           true
+DEAL::end >= end:           true
+DEAL::begin + 5 >= begin:   true
+DEAL::begin + 5 <= begin:   false
+DEAL::end - begin:          12
+DEAL::begin - end:          -12
+DEAL::(begin + 2) - begin:  2
+DEAL::begin + 6 < begin:    false
+DEAL::begin < begin + 1:    true
+DEAL::end - 5 < begin:      false
+DEAL::it+5 == (t2 = (it+5)):true
+DEAL::end - 5 position:     1, 2
+DEAL::begin + 5 position:   2, 1
+DEAL::dealii::TransposeTable<dealii::TransposeTable<double> >
+DEAL::++it == it:           false
+DEAL::++it != it:           true
+DEAL::it++ == it:           true
+DEAL::it++ != it + 1:       true
+DEAL::++(++it) == it + 2:   true
+DEAL::it++ == it + 1:       false
+DEAL::--(++it) == it:       true
+DEAL::(it--)-- == it:       true
+DEAL::end - 5 <= end:       true
+DEAL::end <= end:           true
+DEAL::end >= end:           true
+DEAL::begin + 5 >= begin:   true
+DEAL::begin + 5 <= begin:   false
+DEAL::end - begin:          12
+DEAL::begin - end:          -12
+DEAL::(begin + 2) - begin:  2
+DEAL::begin + 6 < begin:    false
+DEAL::begin < begin + 1:    true
+DEAL::end - 5 < begin:      false
+DEAL::it+5 == (t2 = (it+5)):true
+DEAL::end - 5 position:     1, 2
+DEAL::begin + 5 position:   2, 1
+DEAL::converted iterators are equal: true
+DEAL::OK
diff --git a/tests/base/transpose_table_iterators_1.cc b/tests/base/transpose_table_iterators_1.cc
new file mode 100644 (file)
index 0000000..25e663f
--- /dev/null
@@ -0,0 +1,95 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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 usage of TransposeTable iterators with the table itself
+
+#include "../tests.h"
+
+#include <deal.II/base/table.h>
+
+#include <numeric>
+
+int main()
+{
+  initlog();
+
+  deallog << std::boolalpha;
+
+  // test a non-empty rectangular table
+  TransposeTable<double> table(3, 4);
+  std::iota(table.begin(), table.end(), 1.0);
+  for (const auto entry : table)
+    {
+      deallog << entry.row()
+              << ", "
+              << entry.column()
+              << ", "
+              << entry.value()
+              << std::endl;
+    }
+
+  deallog << "backwards order:" << std::endl;
+  auto it = table.end() - 1;
+  for (; it >= table.begin(); --it)
+    {
+      deallog << it->row()
+              << ", "
+              << it->column()
+              << ", "
+              << it->value()
+              << std::endl;
+    }
+  deallog << "iterator is one before the beginning: "
+          << (it == table.begin() - 1) << std::endl;
+
+  deallog << "every other entry:" << std::endl;
+  it = table.begin();
+  for (; it < table.end(); it += 2)
+    {
+      deallog << it->row()
+              << ", "
+              << it->column()
+              << ", "
+              << it->value()
+              << std::endl;
+    }
+
+  // print every other entry
+  it = table.end() - 1;
+  deallog << "every other entry:" << std::endl;
+  for (; it >= table.begin(); it -= 2)
+    {
+      deallog << it->row()
+              << ", "
+              << it->column()
+              << ", "
+              << it->value()
+              << std::endl;
+    }
+
+  // test some type equalities
+  static_assert(std::is_same<decltype(table.begin()->value()), double &>::value,
+                "The iterator value for a non-const table should not be const.");
+  static_assert(std::is_same<decltype(table.end()->value()), double &>::value,
+                "The iterator value for a non-const table should not be const.");
+
+  const TransposeTable<double> &ref = table;
+  static_assert(std::is_same<decltype(ref.begin()->value()), const double &>::value,
+                "The iterator value for a constant table should be const.");
+  static_assert(std::is_same<decltype(ref.end()->value()), const double &>::value,
+                "The iterator value for a constant table should be const.");
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/base/transpose_table_iterators_1.output b/tests/base/transpose_table_iterators_1.output
new file mode 100644 (file)
index 0000000..4614620
--- /dev/null
@@ -0,0 +1,42 @@
+
+DEAL::0, 0, 1.00000
+DEAL::1, 0, 2.00000
+DEAL::2, 0, 3.00000
+DEAL::0, 1, 4.00000
+DEAL::1, 1, 5.00000
+DEAL::2, 1, 6.00000
+DEAL::0, 2, 7.00000
+DEAL::1, 2, 8.00000
+DEAL::2, 2, 9.00000
+DEAL::0, 3, 10.0000
+DEAL::1, 3, 11.0000
+DEAL::2, 3, 12.0000
+DEAL::backwards order:
+DEAL::2, 3, 12.0000
+DEAL::1, 3, 11.0000
+DEAL::0, 3, 10.0000
+DEAL::2, 2, 9.00000
+DEAL::1, 2, 8.00000
+DEAL::0, 2, 7.00000
+DEAL::2, 1, 6.00000
+DEAL::1, 1, 5.00000
+DEAL::0, 1, 4.00000
+DEAL::2, 0, 3.00000
+DEAL::1, 0, 2.00000
+DEAL::0, 0, 1.00000
+DEAL::iterator is one before the beginning: true
+DEAL::every other entry:
+DEAL::0, 0, 1.00000
+DEAL::2, 0, 3.00000
+DEAL::1, 1, 5.00000
+DEAL::0, 2, 7.00000
+DEAL::2, 2, 9.00000
+DEAL::1, 3, 11.0000
+DEAL::every other entry:
+DEAL::2, 3, 12.0000
+DEAL::0, 3, 10.0000
+DEAL::1, 2, 8.00000
+DEAL::2, 1, 6.00000
+DEAL::0, 1, 4.00000
+DEAL::1, 0, 2.00000
+DEAL::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.