]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update documentation in Point and Tensor. 552/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 15 Feb 2015 21:13:53 +0000 (15:13 -0600)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 16 Feb 2015 21:54:17 +0000 (15:54 -0600)
In particular, document the intent of the various template arguments.

include/deal.II/base/point.h
include/deal.II/base/tensor.h
include/deal.II/base/tensor_base.h
include/deal.II/numerics/data_postprocessor.h

index ebb244e8195ff0cf17e61414e0e96837e5d908c7..d2e3925c7f376648fefec134767f6e4c695a6b9a 100644 (file)
 DEAL_II_NAMESPACE_OPEN
 
 /**
- * The <tt>Point</tt> class represents a point in a space with
- * arbitrary dimension <tt>dim</tt>.
+ * A class that represents a point in a space with arbitrary dimension
+ * <tt>dim</tt>.
  *
- * It is the preferred object to be passed to functions which operate on
- * points in spaces of a priori fixed dimension: rather than using functions
- * like <tt>double f(double x)</tt> and <tt>double f(double x, double y)</tt>,
- * you should use <tt>double f(Point<dim> &p)</tt> instead as it allows writing
- * dimension independent code.
+ * Objects of this class are used to represent points, i.e., vectors
+ * anchored at the origin of a Cartesian vector space. They are, among
+ * other uses, passed to functions that operate on points in spaces of
+ * a priori fixed dimension: rather than using functions like
+ * <tt>double f(double x)</tt> and <tt>double f(double x, double
+ * y)</tt>, you should use <tt>double f(Point<dim> &p)</tt> instead as
+ * it allows writing dimension independent code.
  *
  *
  * <h3>What's a <code>Point@<dim@></code> and what is a <code>Tensor@<1,dim@></code>?</h3>
@@ -67,10 +69,25 @@ DEAL_II_NAMESPACE_OPEN
  * class. Alternatively, as in the case of vector-valued functions,
  * you can use objects of type Vector or <code>std::vector<code>.
  *
+ *
+ * @tparam dim An integer that denotes the dimension of the space in which
+ *   a point lies. This of course equals the number of coordinates that
+ *   identify a point.
+ * @tparam Number The data type in which the coordinates values are
+ *   to be stored. This will, in almost all cases, simply be the default
+ *   @p double, but there are cases where one may want to store coordinates
+ *   in a different (and always scalar) type. An example would be an interval
+ *   type that can store the value of a coordinate as well as its uncertainty.
+ *   Another example would be a type that allows for Automatic Differentiation
+ *   (see, for example, the Sacado type used in step-33) and thereby can
+ *   generate analytic (spatial) derivatives of a function when passed a
+ *   Point object whose coordinates are stored in such a type.
+ *
+ *
  * @ingroup geomprimitives
  * @author Wolfgang Bangerth, 1997
  */
-template <int dim, typename Number>
+template <int dim, typename Number = double>
 class Point : public Tensor<1,dim,Number>
 {
 public:
@@ -122,7 +139,9 @@ public:
          const Number z);
 
   /**
-   * Return a unit vector in coordinate direction <tt>i</tt>.
+   * Return a unit vector in coordinate direction <tt>i</tt>, i.e., a
+   * vector that is zero in all coordinates except for a single 1 in
+   * the <tt>i</tt>th coordinate.
    */
   static Point<dim,Number> unit_vector(const unsigned int i);
 
@@ -137,13 +156,12 @@ public:
   Number &operator () (const unsigned int index);
 
   /*
-   * Plus and minus operators are re-implemented from Tensor<1,dim>
-   * to avoid additional casting.
+   * @name Addition and subtraction of points.
+   * @{
    */
 
   /**
-   * Add two point vectors. If possible, use <tt>operator +=</tt> instead
-   * since this does not need to copy a point at least once.
+   * Add two point vectors.
    */
   Point<dim,Number>   operator + (const Tensor<1,dim,Number> &) const;
 
@@ -171,37 +189,46 @@ public:
   Point<dim,Number>   operator - () const;
 
   /**
-   * Multiply by a factor. If possible, use <tt>operator *=</tt> instead since
-   * this does not need to copy a point at least once.
-   *
-   * There is a commutative complement to this function also
+   * @}
+   */
+
+  /*
+   * @name Multiplication and scaling of points. Dot products. Norms.
+   * @{
    */
-  Point<dim,Number>   operator * (const Number) const;
 
   /**
-   * Returns the scalar product of two vectors.
+   * Multiply the current point by a factor.
    */
-  Number       operator * (const Tensor<1,dim,Number> &) const;
+  Point<dim,Number>   operator * (const Number) const;
 
   /**
-   * Divide by a factor. If possible, use <tt>operator /=</tt> instead since
-   * this does not need to copy a point at least once.
+   * Divide the current point by a factor.
    */
   Point<dim,Number>   operator / (const Number) const;
 
   /**
-   * Returns the scalar product of this point vector with itself, i.e. the
+   * Return the scalar product of the vectors representing two points.
+   */
+  Number              operator * (const Tensor<1,dim,Number> &p) const;
+
+  /**
+   * Return the scalar product of this point vector with itself, i.e. the
    * square, or the square of the norm.
    */
   Number              square () const;
 
   /**
-   * Returns the Euclidean distance of <tt>this</tt> point to the point
+   * Return the Euclidean distance of <tt>this</tt> point to the point
    * <tt>p</tt>, i.e. the <tt>l_2</tt> norm of the difference between the
    * vectors representing the two points.
    */
   Number distance (const Point<dim,Number> &p) const;
 
+  /**
+   * @}
+   */
+
   /**
    * Read or write the data of this object to or from a stream for the purpose
    * of serialization
index df573e876e050c2a45f299c884e39457dde64ac5..87e578240a27630e22eeafc8195d22911d3a1d66 100644 (file)
@@ -26,21 +26,42 @@ template <int rank_, int dim, typename Number> class Tensor;
 template <int dim, typename Number> class Tensor<1,dim,Number>;
 
 /**
- * Provide a general tensor class with an arbitrary rank, i.e. with an
+ * A general tensor class with an arbitrary rank, i.e. with an
  * arbitrary number of indices. The Tensor class provides an indexing operator
  * and a bit of infrastructure, but most functionality is recursively handed
  * down to tensors of rank 1 or put into external templated functions, e.g.
  * the <tt>contract</tt> family.
  *
- * Using this tensor class for objects of rank 2 has advantages over matrices
- * in many cases since the dimension is known to the compiler as well as the
- * location of the data. It is therefore possible to produce far more
- * efficient code than for matrices with runtime-dependent dimension.
- *
- * This class provides an optional template argument for the type of the
- * underlying data. It defaults to @p double values. It can be used to base
- * tensors on @p float or @p complex numbers or any other data type that
- * implements basic arithmetic operations.
+ * Using this tensor class for objects of rank 2 has advantages over
+ * matrices in many cases since the dimension is known to the compiler
+ * as well as the location of the data. It is therefore possible to
+ * produce far more efficient code than for matrices with
+ * runtime-dependent dimension. It also makes the code easier to read
+ * because of the semantic difference between a tensor (an object that
+ * relates to a coordinate system and has transformation properties
+ * with regard to coordinate rotations and transforms) and matrices
+ * (which we consider as operators on arbitrary vector spaces related
+ * to linear algebra things).
+ *
+ * @tparam rank_ An integer that denotes the rank of this tensor. A
+ *   rank-0 tensor is a scalar, a rank-1 tensor is a vector with @p dim
+ *   components, a rank-2 tensor is a matrix with dim-by-dim components,
+ *   etc. There are specializations of this class for rank-0 and rank-1
+ *   tensors. There is also a related class SymmetricTensor for
+ *   tensors of even rank whose elements are symmetric.
+ * @tparam dim An integer that denotes the dimension of the space in which
+ *   this tensor operates. This of course equals the number of coordinates that
+ *   identify a point and rank-1 tensor.
+ * @tparam Number The data type in which the tensor elements are
+ *   to be stored. This will, in almost all cases, simply be the default
+ *   @p double, but there are cases where one may want to store elements
+ *   in a different (and always scalar) type. It can be used to base
+ *   tensors on @p float or @p complex numbers or any other data type that
+ *   implements basic arithmetic operations.
+ *   Another example would be a type that allows for Automatic Differentiation
+ *   (see, for example, the Sacado type used in step-33) and thereby can
+ *   generate analytic (spatial) derivatives of a function that takes a
+ *   tensor as argument.
  *
  * @ingroup geomprimitives
  * @author Wolfgang Bangerth, 1998-2005
index 30cc1f1f49f6a3240c3810d36401e4338861539d..5b52ecf3fe7cdc964079349ac162df643e337554 100644 (file)
@@ -41,7 +41,7 @@ template <typename number> class Vector;
 // this file must be included when using something like Tensor<1,dim>, and
 // Point and Tensor must not be forward declared without the number type
 // specified)
-template <int dim, typename Number=double> class Point;
+template <int dim, typename Number> class Point;
 
 // general template; specialized for rank==1; the general template is in
 // tensor.h
@@ -63,6 +63,21 @@ template <int dim, typename Number> class Tensor<1,dim,Number>;
  * (i.e. @p Number) for all purposes but is part of the Tensor template
  * family.
  *
+ * @tparam dim An integer that denotes the dimension of the space in which
+ *   this tensor operates. This of course equals the number of coordinates that
+ *   identify a point and rank-1 tensor. Since the current object is a rank-0
+ *   tensor (a scalar), this template argument has no meaning for this class.
+ * @tparam Number The data type in which the tensor elements are
+ *   to be stored. This will, in almost all cases, simply be the default
+ *   @p double, but there are cases where one may want to store elements
+ *   in a different (and always scalar) type. It can be used to base
+ *   tensors on @p float or @p complex numbers or any other data type that
+ *   implements basic arithmetic operations.
+ *   Another example would be a type that allows for Automatic Differentiation
+ *   (see, for example, the Sacado type used in step-33) and thereby can
+ *   generate analytic (spatial) derivatives of a function that takes a
+ *   tensor as argument.
+ *
  * @ingroup geomprimitives
  * @author Wolfgang Bangerth, 2009
  */
@@ -179,7 +194,7 @@ public:
   Tensor<0,dim,Number> &operator -= (const Tensor<0,dim,Number> &rhs);
 
   /**
-   * Scale the vector by <tt>factor</tt>, i.e. multiply all coordinates by
+   * Scale the vector by <tt>factor</tt>, i.e. multiply all elements by
    * <tt>factor</tt>.
    */
   Tensor<0,dim,Number> &operator *= (const Number factor);
@@ -287,6 +302,20 @@ private:
  * rank 1, or vector, with as many elements as a point object, but with
  * different physical units), we use the <tt>Tensor<1,dim,Number></tt> class.
  *
+ * @tparam dim An integer that denotes the dimension of the space in which
+ *   this tensor operates. This of course equals the number of coordinates that
+ *   identify a point and rank-1 tensor.
+ * @tparam Number The data type in which the tensor elements are
+ *   to be stored. This will, in almost all cases, simply be the default
+ *   @p double, but there are cases where one may want to store elements
+ *   in a different (and always scalar) type. It can be used to base
+ *   tensors on @p float or @p complex numbers or any other data type that
+ *   implements basic arithmetic operations.
+ *   Another example would be a type that allows for Automatic Differentiation
+ *   (see, for example, the Sacado type used in step-33) and thereby can
+ *   generate analytic (spatial) derivatives of a function that takes
+ *   a tensor as argument.
+ *
  * @ingroup geomprimitives
  * @author Wolfgang Bangerth, 1998-2005
  */
index 5f73960ff2deaf35f0bfcd792f1465ff140210ae..7a52721750fb6c8a3b2d673ac13041c23c9d2051 100644 (file)
@@ -20,6 +20,7 @@
 
 #include <deal.II/base/subscriptor.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/point.h>
 #include <deal.II/lac/vector.h>
 #include <deal.II/fe/fe_update_flags.h>
 #include <deal.II/numerics/data_component_interpretation.h>

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.