]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Improve documentation. 1514/head
authorBruno Turcksin <bruno.turcksin@gmail.com>
Wed, 23 Sep 2015 00:27:06 +0000 (19:27 -0500)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Wed, 23 Sep 2015 00:27:06 +0000 (19:27 -0500)
include/deal.II/base/index_set.h
include/deal.II/lac/read_write_vector.h
tests/lac/readwritevector_assignment.cc [new file with mode: 0644]
tests/lac/readwritevector_assignment.output [new file with mode: 0644]
tests/trilinos/readwritevector.cc
tests/trilinos/readwritevector.mpirun=2.output

index 1c48dcd25f4b22943fa6de64d5e39fd9569f70bf..bfca45fa7162a1a19ef603686a4e2cc0df29119c 100644 (file)
@@ -225,7 +225,7 @@ public:
   /**
    * This function returns the local index of the beginning of the largest range.
    */
-  unsigned int largest_range_index() const;
+  unsigned int largest_range_starting_index() const;
 
   /**
    * Compress the internal representation by merging individual elements with
@@ -1414,12 +1414,12 @@ IndexSet::n_intervals () const
 
 inline
 unsigned int
-IndexSet::largest_range_index() const
+IndexSet::largest_range_starting_index() const
 {
   Assert(ranges.empty()==false, ExcMessage("IndexSet cannot be empty."));
 
   compress();
-  std::vector<Range>::const_iterator main_range=ranges.begin()+largest_range;
+  const std::vector<Range>::const_iterator main_range=ranges.begin()+largest_range;
 
   return main_range->nth_index_in_set;
 }
index 44bfa02a056cc3cc12f0d5b958b83f482a6c4f86..d0e389aacaaaa1ea88cb3d8c4dfa63f009269a06 100644 (file)
@@ -59,11 +59,38 @@ namespace LinearAlgebra
    */
 
   /**
-   * ReadWriteVector stores a subset of elements. It 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.
+   * 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.
+   *
+   * <h3>Storing elements</h3>
+   * 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 <i>global</i> index into a <i>local</i> 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 <i>global</i> 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 <i>in which order</i> 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.
    */
@@ -204,20 +231,24 @@ namespace LinearAlgebra
     ReadWriteVector<Number> &operator = (const Number s);
 
     /**
-     * Returns the size of the vector.
-     *
-     * Note that the result is not equal to the number of element stored locally.
-     * The latter information is returned by n_elements()
+     * 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;
 
     /**
-     * Return the number of elements stored locally.
+     * 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 stores the indices of the elements stored.
+     * Return the IndexSet that represents the indices of the elements stored.
      */
     const IndexSet &get_stored_elements () const;
 
@@ -255,19 +286,22 @@ namespace LinearAlgebra
 
     /**
      * Read access to the data in the position corresponding to @p
-     * global_index.
+     * 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.
+     * 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.
+     * 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().
      */
@@ -275,7 +309,8 @@ namespace LinearAlgebra
 
     /**
      * Read and write access to the data in the position corresponding to @p
-     * global_index.
+     * 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().
      */
@@ -301,16 +336,24 @@ namespace LinearAlgebra
                                OutputIterator           values_begin) const;
 
     /**
-     * Read access to the data field specified by @p local_index. The largest
-     * range of the IndexSet is the first one.
+     * 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. The
-     * largest range of the IndexSet is the first one.
+     * 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).
      */
@@ -649,7 +692,7 @@ namespace LinearAlgebra
   ReadWriteVector<Number>::operator= (const Number s)
   {
     Assert(s==static_cast<Number>(0), ExcMessage("Only 0 can be assigned to a vector."));
-    std::fill(begin(),end(),s);
+    std::fill(begin(),end(),Number());
 
     return *this;
   }
diff --git a/tests/lac/readwritevector_assignment.cc b/tests/lac/readwritevector_assignment.cc
new file mode 100644 (file)
index 0000000..0d647dc
--- /dev/null
@@ -0,0 +1,68 @@
+// ---------------------------------------------------------------------
+//
+// 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 constructor and operator=
+
+#include "../tests.h"
+#include <deal.II/base/index_set.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/read_write_vector.h>
+
+#include <fstream>
+
+
+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> double_vector(is);
+  LinearAlgebra::ReadWriteVector<float> float_vector(float_size);
+  deallog << "double_size " << double_vector.n_elements() <<std::endl;
+  deallog << "float_size " << float_vector.n_elements() <<std::endl;
+
+  double_vector = 0.;
+  for (unsigned int i=0; i<double_vector.n_elements(); ++i)
+    double_vector.local_element(i) += i;
+  for (unsigned int i=0; i<float_vector.n_elements(); ++i)
+    float_vector[i] = i;
+
+  double_vector.print(deallog.get_file_stream());
+  float_vector.print(deallog.get_file_stream());
+
+  float_vector = double_vector;
+  float_vector.print(deallog.get_file_stream());
+
+  LinearAlgebra::ReadWriteVector<double> double_vector2(double_vector);
+  double_vector2.print(deallog.get_file_stream());
+
+  for (unsigned int i=0; i<double_vector.n_elements(); ++i)
+    double_vector2.local_element(i) += i;
+  double_vector = double_vector2;
+  double_vector.print(deallog.get_file_stream());
+}
+
+int main()
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  test();
+
+  return 0;
+}
diff --git a/tests/lac/readwritevector_assignment.output b/tests/lac/readwritevector_assignment.output
new file mode 100644 (file)
index 0000000..2e8a7ec
--- /dev/null
@@ -0,0 +1,55 @@
+
+DEAL::double_size 8
+DEAL::float_size 10
+IndexSet: {[0,1], [10,14], 46}
+
+0.000e+00
+1.000e+00
+2.000e+00
+3.000e+00
+4.000e+00
+5.000e+00
+6.000e+00
+7.000e+00
+IndexSet: {[0,9]}
+
+0.000e+00
+1.000e+00
+2.000e+00
+3.000e+00
+4.000e+00
+5.000e+00
+6.000e+00
+7.000e+00
+8.000e+00
+9.000e+00
+IndexSet: {[0,1], [10,14], 46}
+
+0.000e+00
+1.000e+00
+2.000e+00
+3.000e+00
+4.000e+00
+5.000e+00
+6.000e+00
+7.000e+00
+IndexSet: {[0,1], [10,14], 46}
+
+0.000e+00
+1.000e+00
+2.000e+00
+3.000e+00
+4.000e+00
+5.000e+00
+6.000e+00
+7.000e+00
+IndexSet: {[0,1], [10,14], 46}
+
+0.000e+00
+2.000e+00
+4.000e+00
+6.000e+00
+8.000e+00
+1.000e+01
+1.200e+01
+1.400e+01
index 2a06fb3f39ef7063360b79e53e65eef6f3320781..dd2b9a0a3e2b10a6179d7bd88f9dd63bb5e2eb60 100644 (file)
@@ -78,6 +78,7 @@ void test()
         AssertThrow(readwrite.local_element(i)==comp[i],
                     ExcMessage("Element not copied correctly"));
     }
+  deallog << "OK" <<std::endl;
 }
 
 int main (int argc, char **argv)
index 8b137891791fe96927ad78e64b0aad7bded08bdc..0fd8fc12f0b442283edd8868867114c242b04d11 100644 (file)
@@ -1 +1,2 @@
 
+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.