]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Allow access to vector2d objects using x[i][j] instead of only x(i,j).
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Mar 2002 10:03:41 +0000 (10:03 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Mar 2002 10:03:41 +0000 (10:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@5551 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/vector2d.h
deal.II/doc/news/2002/c-3-3.html
tests/base/Makefile
tests/base/vector2d.cc [new file with mode: 0644]
tests/base/vector2d.checked [new file with mode: 0644]

index ba5826cbb9856373ca3eeae455a8a32f2ac4e50c..61b4a1557175231b5ea1768d68d40e07d0cbd090 100644 (file)
 #include <base/subscriptor.h>
 
 
+template <typename T> class vector2d;
+
+
 /**
  * Two-dimensional array of arbitrary data.
  *
  * This is an implementation of a two-dimensional array with access by
- * pairs of coordinates.
+ * pairs of coordinates. Access is either by @p{x(i,j)} to keep with
+ * matrix notation, or @p{x[i][j]} in accordance to C-style arrays.
  *
  * The name is chosen to be @p{vector2d} to be conformant with the
  * standard C++ classes. This, although array would have been
@@ -33,6 +37,108 @@ template<typename T>
 class vector2d : public Subscriptor
 {
   public:
+                                    /**
+                                     * Object representing one row of
+                                     * a @ref{vector2d} object. It
+                                     * allows to access the elements
+                                     * through the
+                                     * @p{operator[]}. Since objects
+                                     * of this type are also
+                                     * generated through the
+                                     * @p{vector2d::operator[]}, this
+                                     * allows accessing the elements
+                                     * of @p{vector2d} objects just
+                                     * as those of two-dimensional
+                                     * C-style arrays.
+                                     *
+                                     * This class is used for
+                                     * constant @p{vector2d}
+                                     * objects. It only allows
+                                     * read-only access to the
+                                     * elements.
+                                     */
+    class ConstRowAccessor
+    {
+      public:
+                                        /**
+                                         * Constructor.
+                                         */
+       ConstRowAccessor (const vector2d<T>  &parent,
+                         const unsigned int  row);
+
+                                        /**
+                                         * Access operator.
+                                         */
+       T  operator [] (const unsigned int column) const;
+       
+      protected:
+                                        /**
+                                         * Pointer to the parent
+                                         * object. Used only to check
+                                         * access indices for
+                                         * validity.
+                                         */
+       const vector2d<T>  &parent;
+
+                                        /**
+                                         * Pointer to the start of
+                                         * the row pointed to by this
+                                         * object.
+                                         */
+       const T * const     row_start;
+    };    
+
+
+                                    /**
+                                     * Object representing one row of
+                                     * a @ref{vector2d} object. It
+                                     * allows to access the elements
+                                     * through the
+                                     * @p{operator[]}. Since objects
+                                     * of this type are also
+                                     * generated through the
+                                     * @p{vector2d::operator[]}, this
+                                     * allows accessing the elements
+                                     * of @p{vector2d} objects just
+                                     * as those of two-dimensional
+                                     * C-style arrays.
+                                     *
+                                     * This class is used for
+                                     * non-constant @p{vector2d}
+                                     * objects. It allows read and
+                                     * write access to the elements.
+                                     */
+    class NonConstRowAccessor
+    {
+      public:
+                                        /**
+                                         * Constructor.
+                                         */
+       NonConstRowAccessor (vector2d<T>        &parent,
+                            const unsigned int  row);
+
+
+                                        /**
+                                         * Access operator.
+                                         */
+       T & operator [] (const unsigned int column) const;
+       
+      private:
+                                        /**
+                                         * Pointer to the parent
+                                         * object. Used only to check
+                                         * access indices for
+                                         * validity.
+                                         */
+       const vector2d<T>  &parent;
+
+                                        /**
+                                         * Pointer to the start of
+                                         * the row pointed to by this
+                                         * object.
+                                         */
+       T * const           row_start;
+    };    
                                     /**
                                      * Constructor for a quadratic
                                      * @p{rows} by @p{rows} array. The
@@ -145,6 +251,42 @@ class vector2d : public Subscriptor
     T & operator() (const unsigned int i,
                    const unsigned int j);
 
+                                    /**
+                                     * Return an object representing
+                                     * a certain row of this
+                                     * array. This object in turn has
+                                     * an @p{operator[]}, so that the
+                                     * elements of this array can be
+                                     * accessed through @p{x[i][j]}
+                                     * as well as through @p{x(i,j)}.
+                                     *
+                                     * This function is for constant
+                                     * objects. The returned row
+                                     * representing object only
+                                     * allows read access to the
+                                     * elements of the row pointed
+                                     * to.
+                                     */
+    ConstRowAccessor    operator [] (const unsigned int row) const;
+
+                                    /**
+                                     * Return an object representing
+                                     * a certain row of this
+                                     * array. This object in turn has
+                                     * an @p{operator[]}, so that the
+                                     * elements of this array can be
+                                     * accessed through @p{x[i][j]}
+                                     * as well as through @p{x(i,j)}.
+                                     *
+                                     * This function is for
+                                     * non-constant objects. The
+                                     * returned row representing
+                                     * object allows read and write
+                                     * access to the elements of the
+                                     * row pointed to.
+                                     */
+    NonConstRowAccessor operator [] (const unsigned int row);
+    
 
                                     /**
                                      * Set all entries to their
@@ -234,12 +376,65 @@ class vector2d : public Subscriptor
                                     /**
                                      * This is unfortunately needed.
                                      */
-    template <typename T2> friend class FullMatrix;    
+    template <typename T2> friend class FullMatrix;
+
+    friend class ConstRowAccessor;
+    friend class NonConstRowAccessor;    
 };
 
 
 /* --------------------- Template and inline functions ---------------- */
 
+template <typename T>
+inline
+vector2d<T>::ConstRowAccessor::
+ConstRowAccessor (const vector2d<T>  &parent,
+                 const unsigned int  row)
+                       :
+               parent (parent),
+               row_start(&parent.val[row*parent.n_cols()])
+{
+  Assert (row < parent.n_rows(), ExcInternalError());
+};
+
+
+template <typename T>
+inline
+T
+vector2d<T>::ConstRowAccessor::
+operator [] (const unsigned int column) const
+{
+  Assert (column < parent.n_cols(), ExcInternalError());
+  return *(row_start+column);
+};
+
+
+
+template <typename T>
+inline
+vector2d<T>::NonConstRowAccessor::
+NonConstRowAccessor (vector2d<T>        &parent,
+                    const unsigned int  row)
+                       :
+               parent (parent),
+               row_start(&parent.val[row*parent.n_cols()])
+{
+  Assert (row < parent.n_rows(), ExcInternalError());
+};
+
+
+template <typename T>
+inline
+T &
+vector2d<T>::NonConstRowAccessor::
+operator [] (const unsigned int column) const
+{
+  Assert (column < parent.n_cols(), ExcInternalError());
+  return *(row_start+column);
+};
+
+
+
 template <typename T>
 inline
 vector2d<T>::~vector2d ()
@@ -477,7 +672,8 @@ vector2d<T>::operator() (const unsigned int i, const unsigned int j) const
 
 
 template <typename T>
-inline T &
+inline
+T &
 vector2d<T>::operator() (const unsigned int i, const unsigned int j)
 {
   Assert (i<num_rows, ExcIndexRange (i, 0, num_rows));
@@ -487,6 +683,36 @@ vector2d<T>::operator() (const unsigned int i, const unsigned int j)
 
 
 
+template <typename T>
+inline
+typename vector2d<T>::ConstRowAccessor 
+vector2d<T>::operator [] (const unsigned int row) const 
+{
+                                  // Note: check for validity of row
+                                  // number is done in the
+                                  // constructor of the created
+                                  // object
+  return ConstRowAccessor(*this, row);
+};
+
+
+
+template <typename T>
+inline
+typename vector2d<T>::NonConstRowAccessor 
+vector2d<T>::operator [] (const unsigned int row)
+{
+                                  // Note: check for validity of row
+                                  // number is done in the
+                                  // constructor of the created
+                                  // object
+  return NonConstRowAccessor(*this, row);
+};
+
+
+
+
+
 template <typename T>
 inline
 const T *
index 87114332067b74bd58d523d5c4071c592eac85a0..7506939c154b4cd5d88676a4c85ae64adddfe90d 100644 (file)
@@ -67,6 +67,17 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>base</h3>
 
 <ol>
+  <li> <p>
+       New: The <code class="class">vector2d</code> class now not only
+       allows access to elements through the <code
+       class="member">operator()(unsingned int,unsigned int)</code>
+       (i.e. matrix or Fortran style access), but also through nested
+       brackets via an <code class="member">operator[]</code>
+       (i.e. like to a two-dimensional C-style array).
+       <br>
+       (WB 2002/03/08)
+       </p> 
+
   <li> <p>
        Changed: The function <code class="class">MultithreadInfo</code>::
        <code class="member">get_n_cpus</code> now reports the proper number
index 6e99a2c61a83336ba8ae05959aca2bebcdebf148..b60f5bcd694a1e4150d35a8263f7f635c436a16e 100644 (file)
@@ -2,7 +2,7 @@
 # Generated automatically from Makefile.in by configure.
 ############################################################
 # $Id$
-# Copyright (C) 2000, 2001 by the deal.II authors
+# Copyright (C) 2000, 2001, 2002 by the deal.II authors
 ############################################################
 
 ############################################################
@@ -27,10 +27,11 @@ quadrature_test.exe: quadrature_test.go    $(libraries)
 reference.exe      : reference.go abort.go $(libraries)
 tensor.exe         : tensor.go             $(libraries)
 timer.exe          : timer.go              $(libraries)
+vector2d.exe       : vector2d.go           $(libraries)
 auto_derivative_function.exe        : auto_derivative_function.go            $(libraries)
 
 
-tests = logtest reference quadrature_test tensor timer polynomial1d polynomial_test auto_derivative_function
+tests = logtest reference quadrature_test tensor timer polynomial1d polynomial_test auto_derivative_function vector2d
 
 ############################################################
 
diff --git a/tests/base/vector2d.cc b/tests/base/vector2d.cc
new file mode 100644 (file)
index 0000000..467399b
--- /dev/null
@@ -0,0 +1,56 @@
+//----------------------------  $RCSfile$  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  $RCSfile$  ---------------------------
+
+
+#include <cmath>
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <cstdlib>
+
+#include <base/logstream.h>
+#include <base/vector2d.h>
+
+const int entries[9] = { 11,12,13,21,22,23,31,32,33 };
+
+int
+main ()
+{
+  std::ofstream logfile("vector2d.output");
+  logfile.setf(std::ios::fixed);
+  logfile.precision(3);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  vector2d<double> Td(3,3);
+  vector2d<int> Ti(3,3);
+  
+  for (unsigned int i=0; i<9; ++i)
+    {
+      Td[i/3][i%3] = entries[i];
+      Ti[i/3][i%3] = entries[i];
+    };
+
+  for (unsigned int i=0; i<3; ++i)
+    for (unsigned int j=0; j<3; ++j)
+      {
+       Assert (Td[i][j] == Td(i,j), ExcInternalError());
+       Assert (Ti[i][j] == Ti(i,j), ExcInternalError());
+       Assert (Ti[i][j] == Td(i,j), ExcInternalError());
+
+       logfile << i << " " << j << " " << Td[i][j] << " ok" << std::endl;
+      };
+};
+
+
+      
diff --git a/tests/base/vector2d.checked b/tests/base/vector2d.checked
new file mode 100644 (file)
index 0000000..64879da
--- /dev/null
@@ -0,0 +1,10 @@
+
+0 0 11.000 ok
+0 1 12.000 ok
+0 2 13.000 ok
+1 0 21.000 ok
+1 1 22.000 ok
+1 2 23.000 ok
+2 0 31.000 ok
+2 1 32.000 ok
+2 2 33.000 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.