]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Introduce more powerful constructors for the table classes.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 21 Dec 2013 14:53:37 +0000 (14:53 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 21 Dec 2013 14:53:37 +0000 (14:53 +0000)
git-svn-id: https://svn.dealii.org/trunk@32084 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/base/table.h
tests/base/table_04.cc [new file with mode: 0644]
tests/base/table_04.output [new file with mode: 0644]

index 273a5504d54ed8ad0a457f3d78df768ddc011dc7..2d21073c7e67c8fd19f46699e64de421cafa4b2d 100644 (file)
@@ -88,6 +88,12 @@ inconvenience this causes.
   order in which table elements are initialized, allowing to switch between
   C- and Fortran-style table layouts.
   <br>
+  Along with the TableBase::fill function, the Table classes of various
+  ranks have also gotten constructors that allow the in-place construction
+  not only of a table of correct size, but already initialized from
+  somewhere. This finally allows to mark Table objects as const by creating
+  them already with the correct content.
+  <br>
   (Wolfgang Bangerth, 2013/12/21)
   </li>
 
index 984446817b2e3a3fa4052a1d4c4b7dccf99d7f60..d898917f70c5d9d7e9233203d72850bc9ae8e3f1 100644 (file)
@@ -481,6 +481,18 @@ public:
    */
   TableBase (const TableIndices<N> &sizes);
 
+  /**
+   * Constructor. Initialize the
+   * array with the given
+   * dimensions in each index
+   * component, and then initialize the elements of the table using the
+   * second and third argument by calling fill(entries,C_style_indexing).
+   */
+  template <typename ForwardIterator>
+  TableBase (const TableIndices<N> &sizes,
+             ForwardIterator entries,
+             const bool      C_style_indexing = true);
+
   /**
    * Copy constructor. Performs a
    * deep copy.
@@ -806,6 +818,50 @@ public:
    */
   Table (const unsigned int size);
 
+  /**
+   * Constructor. Create a table with a given size and initialize it from
+   * a set of iterators.
+   *
+   * This function is entirely equivalent to creating a table <code>t</code>
+   * of the given size and then calling
+   * @code
+   *   t.fill (entries, C_style_indexing);
+   * @endcode
+   * on it, using the TableBase::fill() function where the arguments are
+   * explained in more detail. The point, however, is that that is only
+   * possible if the table can be changed after running the constructor,
+   * whereas calling the current constructor allows sizing and initializing
+   * an object right away so that it can be marked const.
+   *
+   * Using this constructor, you can do things like this:
+   * @code
+   *   const double values[] = { 1, 2, 3 };
+   *   const Table<1,double> t(3, entries, true);
+   * @endcode
+   * You can also initialize a table right from a file, using input iterators:
+   * @code
+   *   std::ifstream input ("myfile");
+   *   const Table<1,double> t(3,
+   *                           std::istream_iterator<double>(input),
+   *                           true);
+   * @endcode
+   *
+   *
+   * @param size The size of this one-dimensional table.
+   * @param entries An iterator to a set of elements from which to
+   *   initialize this table. It is assumed that iterator can be
+   *   incremented and dereferenced a sufficient number of times
+   *   to fill this table.
+   * @param C_style_indexing If true, run over elements of the
+   *   table with the last index changing fastest as we dereference
+   *   subsequent elements of the input range. If false, change
+   *   the first index fastest.
+   */
+  template <typename ForwardIterator>
+  Table (const unsigned int size,
+         ForwardIterator entries,
+         const bool      C_style_indexing = true);
+
   /**
    * Access operator. Since this is
    * a one-dimensional object, this
@@ -898,6 +954,52 @@ public:
   Table (const unsigned int size1,
          const unsigned int size2);
 
+  /**
+   * Constructor. Create a table with a given size and initialize it from
+   * a set of iterators.
+   *
+   * This function is entirely equivalent to creating a table <code>t</code>
+   * of the given size and then calling
+   * @code
+   *   t.fill (entries, C_style_indexing);
+   * @endcode
+   * on it, using the TableBase::fill() function where the arguments are
+   * explained in more detail. The point, however, is that that is only
+   * possible if the table can be changed after running the constructor,
+   * whereas calling the current constructor allows sizing and initializing
+   * an object right away so that it can be marked const.
+   *
+   * Using this constructor, you can do things like this:
+   * @code
+   *   const double values[] = { 1, 2, 3, 4, 5, 6 };
+   *   const Table<2,double> t(2, 3, entries, true);
+   * @endcode
+   * You can also initialize a table right from a file, using input iterators:
+   * @code
+   *   std::ifstream input ("myfile");
+   *   const Table<2,double> t(2, 3,
+   *                           std::istream_iterator<double>(input),
+   *                           true);
+   * @endcode
+   *
+   *
+   * @param size1 The size of this table in the first dimension.
+   * @param size2 The size of this table in the second dimension.
+   * @param entries An iterator to a set of elements from which to
+   *   initialize this table. It is assumed that iterator can be
+   *   incremented and dereferenced a sufficient number of times
+   *   to fill this table.
+   * @param C_style_indexing If true, run over elements of the
+   *   table with the last index changing fastest as we dereference
+   *   subsequent elements of the input range. If false, change
+   *   the first index fastest.
+   */
+  template <typename ForwardIterator>
+  Table (const unsigned int size1,
+         const unsigned int size2,
+         ForwardIterator entries,
+         const bool      C_style_indexing = true);
+
   /**
    * Reinitialize the object. This
    * function is mostly here for
@@ -1078,6 +1180,55 @@ public:
          const unsigned int size2,
          const unsigned int size3);
 
+  /**
+   * Constructor. Create a table with a given size and initialize it from
+   * a set of iterators.
+   *
+   * This function is entirely equivalent to creating a table <code>t</code>
+   * of the given size and then calling
+   * @code
+   *   t.fill (entries, C_style_indexing);
+   * @endcode
+   * on it, using the TableBase::fill() function where the arguments are
+   * explained in more detail. The point, however, is that that is only
+   * possible if the table can be changed after running the constructor,
+   * whereas calling the current constructor allows sizing and initializing
+   * an object right away so that it can be marked const.
+   *
+   * Using this constructor, you can do things like this (shown here for
+   * a two-dimensional table, but the same works for the current class):
+   * @code
+   *   const double values[] = { 1, 2, 3, 4, 5, 6 };
+   *   const Table<2,double> t(2, 3, entries, true);
+   * @endcode
+   * You can also initialize a table right from a file, using input iterators:
+   * @code
+   *   std::ifstream input ("myfile");
+   *   const Table<2,double> t(2, 3,
+   *                           std::istream_iterator<double>(input),
+   *                           true);
+   * @endcode
+   *
+   *
+   * @param size1 The size of this table in the first dimension.
+   * @param size2 The size of this table in the second dimension.
+   * @param size3 The size of this table in the third dimension.
+   * @param entries An iterator to a set of elements from which to
+   *   initialize this table. It is assumed that iterator can be
+   *   incremented and dereferenced a sufficient number of times
+   *   to fill this table.
+   * @param C_style_indexing If true, run over elements of the
+   *   table with the last index changing fastest as we dereference
+   *   subsequent elements of the input range. If false, change
+   *   the first index fastest.
+   */
+  template <typename ForwardIterator>
+  Table (const unsigned int size1,
+         const unsigned int size2,
+         const unsigned int size3,
+         ForwardIterator entries,
+         const bool      C_style_indexing = true);
+
   /**
    * Access operator. Generate an
    * object that accesses the
@@ -1751,6 +1902,20 @@ TableBase<N,T>::TableBase (const TableIndices<N> &sizes)
 
 
 
+template <int N, typename T>
+template <typename ForwardIterator>
+TableBase<N,T>::
+TableBase (const TableIndices<N> &sizes,
+           ForwardIterator entries,
+           const bool      C_style_indexing)
+{
+  reinit (sizes);
+  fill (entries, C_style_indexing);
+}
+
+
+
+
 template <int N, typename T>
 TableBase<N,T>::TableBase (const TableBase<N,T> &src)
   :
@@ -2289,6 +2454,20 @@ Table<1,T>::Table (const unsigned int size)
 
 
 
+template <typename T>
+template <typename ForwardIterator>
+inline
+Table<1,T>::Table (const unsigned int size,
+                   ForwardIterator entries,
+                   const bool C_style_indexing)
+  :
+  TableBase<1,T> (TableIndices<1> (size),
+                  entries,
+                  C_style_indexing)
+{}
+
+
+
 template <typename T>
 inline
 typename std::vector<T>::const_reference
@@ -2375,6 +2554,21 @@ Table<2,T>::Table (const unsigned int size1,
 
 
 
+template <typename T>
+template <typename ForwardIterator>
+inline
+Table<2,T>::Table (const unsigned int size1,
+                   const unsigned int size2,
+                   ForwardIterator entries,
+                   const bool C_style_indexing)
+  :
+  TableBase<2,T> (TableIndices<2> (size1, size2),
+                  entries,
+                  C_style_indexing)
+{}
+
+
+
 template <typename T>
 inline
 void
@@ -2629,6 +2823,22 @@ Table<3,T>::Table (const unsigned int size1,
 
 
 
+template <typename T>
+template <typename ForwardIterator>
+inline
+Table<3,T>::Table (const unsigned int size1,
+                   const unsigned int size2,
+                   const unsigned int size3,
+                   ForwardIterator entries,
+                   const bool C_style_indexing)
+  :
+  TableBase<3,T> (TableIndices<3> (size1, size2, size3),
+                  entries,
+                  C_style_indexing)
+{}
+
+
+
 template <typename T>
 inline
 dealii::internal::TableBaseAccessors::Accessor<3,T,true,2>
diff --git a/tests/base/table_04.cc b/tests/base/table_04.cc
new file mode 100644 (file)
index 0000000..940c354
--- /dev/null
@@ -0,0 +1,143 @@
+// ---------------------------------------------------------------------
+// $Id: table.cc 31349 2013-10-20 19:07:06Z maier $
+//
+// Copyright (C) 1998 - 2013 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.
+//
+// ---------------------------------------------------------------------
+
+
+// like_02, but don't just check TableBase::fill but instead use the
+// constructors that already initialize
+
+
+#include "../tests.h"
+#include <cmath>
+#include <fstream>
+#include <iomanip>
+#include <cstdlib>
+
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/table.h>
+
+
+int
+main ()
+{
+  std::ofstream logfile("output");
+  deallog << std::fixed;
+  deallog << std::setprecision(0);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  // rank=1
+  {
+    deallog << "rank=1" << std::endl;
+    
+    const double entries[] = { 1, 2, 3 };
+
+    {
+      Table<1,double> t (3, entries, true);
+
+      for (unsigned int i=0; i<t.size()[0]; ++i)
+       deallog << t[i] << ' ';
+      deallog << std::endl;
+    }
+    
+    // passing false as second argument shouldn't
+    // make a difference for rank-1 tables
+    {
+      Table<1,double> t (3, entries, false);
+      for (unsigned int i=0; i<t.size()[0]; ++i)
+       deallog << t[i] << ' ';
+      deallog << std::endl;
+    }
+  }
+
+
+  // rank=2
+  {
+    deallog << "rank=2" << std::endl;
+    
+    const double entries[] = { 1, 2, 3, 4, 5, 6 };
+
+    // create a 2x3 table from this
+    {
+      Table<2,double> t (2, 3, entries, true);
+
+      for (unsigned int i=0; i<t.size()[0]; ++i)
+       {
+         for (unsigned int j=0; j<t.size()[1]; ++j)
+           deallog << t[i][j] << ' ';
+         deallog << std::endl;
+       }
+    }
+    
+    // same data, same table, but filled in transpose ordering
+    {
+      Table<2,double> t (2, 3, entries, false);
+
+      for (unsigned int i=0; i<t.size()[0]; ++i)
+       {
+         for (unsigned int j=0; j<t.size()[1]; ++j)
+           deallog << t[i][j] << ' ';
+         deallog << std::endl;
+       }
+    }
+  }
+
+
+  // rank=3
+  {
+    deallog << "rank=3" << std::endl;
+    
+    const double entries[] = { 1, 2, 3, 4, 5, 6,
+                              7, 8, 9, 10, 11, 12 };
+
+    // create a 2x3x2 table from this
+    {
+      Table<3,double> t (2,3,2, entries, true);
+
+      for (unsigned int i=0; i<t.size()[0]; ++i)
+       {
+         for (unsigned int j=0; j<t.size()[1]; ++j)
+           {
+             deallog << '(';
+             for (unsigned int k=0; k<t.size()[2]; ++k)
+               deallog << t[i][j][k] << ' ';
+             deallog << ')';
+           }
+         deallog << std::endl;
+       }
+    }
+    
+    // same data, same table, but filled in transpose ordering
+    {
+      Table<3,double> t (2,3,2, entries, false);
+
+      for (unsigned int i=0; i<t.size()[0]; ++i)
+       {
+         for (unsigned int j=0; j<t.size()[1]; ++j)
+           {
+             deallog << '(';
+             for (unsigned int k=0; k<t.size()[2]; ++k)
+               deallog << t[i][j][k] << ' ';
+             deallog << ')';
+           }
+         deallog << std::endl;
+       }
+    }
+  }
+}
+
+
+
diff --git a/tests/base/table_04.output b/tests/base/table_04.output
new file mode 100644 (file)
index 0000000..220d255
--- /dev/null
@@ -0,0 +1,14 @@
+
+DEAL::rank=1
+DEAL::1. 2. 3. 
+DEAL::1. 2. 3. 
+DEAL::rank=2
+DEAL::1. 2. 3. 
+DEAL::4. 5. 6. 
+DEAL::1. 3. 5. 
+DEAL::2. 4. 6. 
+DEAL::rank=3
+DEAL::(1. 2. )(3. 4. )(5. 6. )
+DEAL::(7. 8. )(9. 10. )(11. 12. )
+DEAL::(1. 7. )(3. 9. )(5. 11. )
+DEAL::(2. 8. )(4. 10. )(6. 12. )

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.