]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Move the point class from deal.II/include/grid to here.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 6 Nov 1998 12:09:10 +0000 (12:09 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 6 Nov 1998 12:09:10 +0000 (12:09 +0000)
git-svn-id: https://svn.dealii.org/trunk@640 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/function.h
deal.II/base/include/base/point.h [new file with mode: 0644]
deal.II/base/include/base/quadrature.h
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/include/grid/tria_boundary.h
deal.II/deal.II/include/grid/tria_iterator.h
deal.II/deal.II/include/grid/tria_levels.h

index 934d016c33d57a1257e385ecf02fe8054cfe11ca..c7a568791736b02403a918a92028649112e87f1a 100644 (file)
@@ -8,7 +8,7 @@
 
 #include <base/exceptions.h>
 #include <vector>
-#include <grid/point.h>
+#include <base/point.h>
 
 
 
diff --git a/deal.II/base/include/base/point.h b/deal.II/base/include/base/point.h
new file mode 100644 (file)
index 0000000..a546fe6
--- /dev/null
@@ -0,0 +1,330 @@
+/*----------------------------   point.h     ---------------------------*/
+/*      $Id$                 */
+/*      Copyright W. Bangerth, University of Heidelberg, 1998 */
+#ifndef __point_H
+#define __point_H
+/*----------------------------   point.h     ---------------------------*/
+
+
+#include <base/exceptions.h>
+#include <base/tensor_base.h>
+
+
+/**
+ * The #Point# class provides for a point or vector in a space with arbitrary
+ * dimension #dim#.
+ *
+ * It is the preferred object to be passed to functions which
+ * operate on points in spaces of a priori unknown dimension: rather than
+ * using functions like #double f(double x)# and #double f(double x, double y)#,
+ * you use double #f(Point<dim> &p)#.
+ *
+ * #Point# also serves as a starting point for the implementation of the
+ * geometrical primitives like #Polyhedron#, #Triangle#, etc.
+ *
+ * #Point#s can also be thought of as vectors, i.e. points in a vector space
+ * without an obvious meaning. For instance, it may be suitable to let the
+ * gradient of a function be a #point# vector:
+ * #Point<dim> gradient_of_f (const Point<dim> &x)#. #Point#s have all
+ * functionality for this, e.g. scalar products, addition etc. However, you
+ * should also consider using the simpler #Tensor<1,dim># class, which seems
+ * more suited to gradients.
+ *
+ * @author Wolfgang Bangerth, 1997
+ */
+template <int dim>
+class Point : public Tensor<1,dim> {
+  public:
+                                    /**
+                                     * Default constructor.
+                                     */
+                                    /**
+                                     * Constructor. Initialize all entries
+                                     * to zero.
+                                     */
+    explicit Point ();
+    
+                                    /**
+                                     *  Constructor for one dimensional points. This
+                                     *  function is only implemented for #dim==1#
+                                     *  since the usage is considered unsafe
+                                     *  for points with #dim!=1#.
+                                     */
+    explicit Point (const double x);
+
+                                    /**
+                                     *  Constructor for two dimensional points. This
+                                     *  function is only implemented for #dim==2#
+                                     *  since the usage is considered unsafe
+                                     *  for points with #dim!=2#.
+                                     */
+    Point (const double x, const double y);
+    
+                                    /**
+                                     *  Constructor for three dimensional points. This
+                                     *  function is only implemented for #dim==3#
+                                     *  since the usage is considered unsafe
+                                     *  for points with #dim!=3#.
+                                     */
+    Point (const double x, const double y, const double z);
+
+                                    /**
+                                     *  Read access to the #index#th coordinate.
+                                     */
+    double   operator () (const unsigned int index) const;
+
+                                    /**
+                                     *  Read and write access to the #index#th
+                                     *  coordinate.
+                                     */
+    double & operator () (const unsigned int index);
+
+                                    /**
+                                     *  Add two point vectors. If possible, use
+                                     *  #operator +=# instead since this does not
+                                     *  need to copy a point at least once.
+                                     */
+    Point<dim>   operator + (const Point<dim> &) const;
+
+                                    /**
+                                     *  Subtract two point vectors. If possible, use
+                                     *  #operator +=# instead since this does not
+                                     *  need to copy a point at least once.
+                                     */
+    Point<dim>   operator - (const Point<dim> &) const;
+
+                                    /**
+                                     *  Multiply by a factor. If possible, use
+                                     *  #operator *=# instead since this does not
+                                     *  need to copy a point at least once.
+                                     *
+                                     * There is a commutative complement to this
+                                     * function also
+                                     */
+    Point<dim>   operator * (const double) const;
+
+                                    /**
+                                     *  Divide by a factor. If possible, use
+                                     *  #operator /=# instead since this does not
+                                     *  need to copy a point at least once.
+                                     */
+    Point<dim>   operator / (const double) const;
+
+                                    /**
+                                     *  Add another vector, i.e. move this point by
+                                     *  the given offset.
+                                     */
+    Point<dim> & operator += (const Point<dim> &);
+                                    /**
+                                     *  Subtract another vector.
+                                     */
+    Point<dim> & operator -= (const Point<dim> &);
+
+                                    /**
+                                     *  Scale the vector by #factor#, i.e. multiply
+                                     *  all coordinates by #factor#.
+                                     */
+    Point<dim> & operator *= (const double &factor);
+
+                                    /**
+                                     *  Scale the vector by #1/factor#.
+                                     */
+    Point<dim> & operator /= (const double &factor);
+
+                                    /**
+                                     *  Returns the scalar product of two vectors.
+                                     */
+    double              operator * (const Point<dim> &) const;
+
+                                    /**
+                                     *  Returns the scalar product of this point
+                                     *  vector with itself, i.e. the square, or
+                                     *  the square of the norm.
+                                     */
+    double              square () const;
+
+
+                                    /**
+                                     *  Exception
+                                     */
+    DeclException1 (ExcDimTooSmall,
+                   int,
+                   << "Given dimensions must be >= 1, but was " << arg1);
+                                    /**
+                                     *  Exception
+                                     */
+    DeclException1 (ExcInvalidIndex,
+                   int,
+                   << "Invalid index " << arg1);
+};
+
+
+
+
+
+
+
+/*------------------------------- Inline functions: Point ---------------------------*/
+
+
+
+template <int dim>
+inline
+Point<dim>::Point () :
+               Tensor<1,dim>() {};
+
+
+
+template <>
+inline
+Point<1>::Point (const double x) {
+  values[0] = x;
+};
+
+
+
+template <>
+inline
+Point<2>::Point (const double x, const double y) {
+  values[0] = x;
+  values[1] = y;
+};
+
+
+
+template <>
+inline
+Point<3>::Point (const double x, const double y, const double z) {
+  values[0] = x;
+  values[1] = y;
+  values[2] = z;
+};
+
+
+
+template <int dim>
+inline
+double Point<dim>::operator () (const unsigned int index) const {
+  Assert (index<dim, ExcInvalidIndex (index));
+  return values[index];
+};
+
+
+
+template <int dim>
+inline
+double & Point<dim>::operator () (const unsigned int index) {
+  Assert (index<dim, ExcInvalidIndex (index));
+  return values[index];
+};
+
+
+
+template <int dim>
+inline
+Point<dim> Point<dim>::operator + (const Point<dim> &p) const {
+  return (Point<dim>(*this) += p);
+};
+
+
+
+template <int dim>
+inline
+Point<dim> Point<dim>::operator - (const Point<dim> &p) const {
+  return (Point<dim>(*this) -= p);
+};
+
+
+
+template <int dim>
+inline
+Point<dim> Point<dim>::operator * (const double factor) const {
+  return (Point<dim>(*this) *= factor);
+};
+
+
+
+template <int dim>
+inline
+Point<dim> operator * (const double factor, const Point<dim> &p) {
+  return p*factor;
+};
+
+
+
+template <int dim>
+inline
+Point<dim> Point<dim>::operator / (const double factor) const {
+  return (Point<dim>(*this) /= factor);
+};
+
+
+  
+template <int dim>
+inline
+Point<dim> & Point<dim>::operator += (const Point<dim> &p) {
+  for (unsigned int i=0; i<dim; ++i)
+    values[i] += p.values[i];
+  return *this;
+};
+
+
+
+template <int dim>
+inline
+Point<dim> & Point<dim>::operator -= (const Point<dim> &p) {
+  for (unsigned int i=0; i<dim; ++i)
+    values[i] -= p.values[i];
+  return *this;
+};
+
+
+
+template <int dim>
+inline
+Point<dim> & Point<dim>::operator *= (const double &s) {
+  for (unsigned int i=0; i<dim; ++i)
+    values[i] *= s;
+  return *this;
+};
+
+
+
+template <int dim>
+inline
+Point<dim> & Point<dim>::operator /= (const double &s) {
+  for (unsigned int i=0; i<dim; ++i)
+    values[i] /= s;
+  return *this;
+};
+
+
+
+template <int dim>
+inline
+double Point<dim>::operator * (const Point<dim> &p) const {
+  double q=0;
+  for (unsigned int i=0; i<dim; ++i)
+    q += values[i] * p.values[i];
+  return q;
+};
+
+
+
+template <int dim>
+inline
+double Point<dim>::square () const {
+  double q=0;
+  for (unsigned int i=0; i<dim; ++i)
+    q += values[i] * values[i];
+  return q;
+};
+
+
+
+
+
+/*----------------------------   point.h     ---------------------------*/
+/* end of #ifndef __point_H */
+#endif
+/*----------------------------   point.h     ---------------------------*/
index 7abc27cb0325a9ec57059ad2afa4749f74357a9a..f42c50f56e7eaff2befc3a790474da1e5fc4c1ff 100644 (file)
@@ -6,7 +6,7 @@
 /*----------------------------   quadrature.h     ---------------------------*/
 
 
-#include <grid/point.h>
+#include <base/point.h>
 #include <vector>
 
 
index 492b427b814f266cd0e9266c6bc72342e98b7f04..525275ad1d32d9dd8f8a8962a9eb507c7b1a50d9 100644 (file)
@@ -7,7 +7,7 @@
 
 #include <base/exceptions.h>
 #include <base/subscriptor.h>
-#include <grid/point.h>
+#include <base/point.h>
 #include <grid/dof.h>
 #include <grid/geometry_info.h>
 #include <lac/dfmatrix.h>
index 35ce803a3a0065e5cbeb7f214b26e3f0f2492054..bf0768e08ec70cb2361658162851f469a1078e2e 100644 (file)
@@ -9,7 +9,7 @@
 #include <base/exceptions.h>
 #include <base/subscriptor.h>
 #include <lac/dfmatrix.h>
-#include <grid/point.h>
+#include <base/point.h>
 #include <grid/tria.h>
 #include <fe/fe_update_flags.h>
 
index b17b6020bf860f8ac1f64bdb8bbfd0e96833d030..f67067033334ac13f3f8b149ff28a51d1ea474b3 100644 (file)
@@ -6,7 +6,7 @@
 /*----------------------------   tria.h     ---------------------------*/
 
 #include <vector>
-#include <grid/point.h>
+#include <base/point.h>
 #include <grid/geometry_info.h>
 
 
index 0dd06b3172750a63e9d4f84409b8e7bc7c7e5f89..f67f7b4515666f015fddf59ba38e0a4d075c47dd 100644 (file)
@@ -5,7 +5,7 @@
 #define __tria_boundary_H
 /*----------------------------   boundary-function.h     ---------------------------*/
 
-#include <grid/point.h>
+#include <base/point.h>
 #include <grid/geometry_info.h>
 
 
index ef6cc9c87b371bfd270aa2b19dab1c9b0f019eaf..65dea7aae8c5cef0aaefba426e29da9f2d7a4226 100644 (file)
@@ -10,7 +10,7 @@
 #include <iterator>
 #include <iostream>
 #include <base/exceptions.h>
-#include <grid/point.h>
+#include <base/point.h>
 #include <grid/tria_accessor.h>
 
 
index 397ba597cb793b2a8048ddbe7cecc971022cabd2..4a212e07ba1534be365037052903c713718530c5 100644 (file)
@@ -10,7 +10,7 @@
 #include <vector>
 #include <grid/tria_line.h>
 #include <grid/tria_quad.h>
-#include <grid/point.h>
+#include <base/point.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.