]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Periodicity as Tensor<1,spacedim>.
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 6 Apr 2016 15:20:23 +0000 (17:20 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Mon, 11 Apr 2016 07:41:41 +0000 (09:41 +0200)
include/deal.II/grid/manifold.h
include/deal.II/grid/manifold_lib.h
source/grid/manifold.cc
source/grid/manifold_lib.cc
tests/manifold/chart_manifold_03.cc
tests/manifold/chart_manifold_03_embedded.cc
tests/manifold/chart_manifold_06.cc
tests/manifold/chart_manifold_06_embedded.cc
tests/manifold/chart_manifold_08.cc
tests/manifold/flat_manifold_03.cc
tests/manifold/flat_manifold_06.cc

index 617cd427d4e1571b848b7733ff67464582dc25f6..b76e8be958930174a9c085bbfc1789dfb4e3a0f4 100644 (file)
@@ -529,7 +529,7 @@ public:
    * guaranteed to lie in the periodicity box plus or minus
    * tolerance*periodicity.norm().
    */
-  FlatManifold (const Point<spacedim> &periodicity = Point<spacedim>(),
+  FlatManifold (const Tensor<1,spacedim> &periodicity = Tensor<1,spacedim>(),
                 const double tolerance=1e-10);
 
   /**
@@ -599,7 +599,7 @@ public:
   /**
    * Return the periodicity of this Manifold.
    */
-  const Point<spacedim> &get_periodicity() const;
+  const Tensor<1,spacedim> &get_periodicity() const;
 
 private:
   /**
@@ -615,12 +615,12 @@ private:
    * A periodicity 0 along one direction means no periodicity. This is the
    * default value for all directions.
    */
-  const Point<spacedim> periodicity;
+  const Tensor<1,spacedim> periodicity;
 
-  DeclException4(ExcPeriodicBox, int, Point<spacedim>, Point<spacedim>, double,
+  DeclException4(ExcPeriodicBox, int, Point<spacedim>, double, double,
                  << "The component number " << arg1 << " of the point [ " << arg2
                  << " ] is not in the interval [ " << -arg4
-                 << ", " << arg3[arg4] << "), bailing out.");
+                 << ", " << arg3 << "), bailing out.");
 
   /**
    * Relative tolerance. This tolerance is used to compute distances in double
@@ -736,7 +736,7 @@ public:
    * of (2*pi-eps) and (eps) is not pi, but 2*pi (or zero), since, on the
    * manifold, these two points are at distance 2*eps and not (2*pi-eps)
    */
-  ChartManifold(const Point<chartdim> &periodicity = Point<chartdim>());
+  ChartManifold(const Tensor<1,chartdim> &periodicity = Tensor<1,chartdim>());
 
   /**
    * Destructor. Does nothing here, but needs to be declared to make it
@@ -856,7 +856,7 @@ public:
   /**
    * Return the periodicity associated with the submanifold.
    */
-  const Point<chartdim> &get_periodicity() const;
+  const Tensor<1,chartdim> &get_periodicity() const;
 
 private:
   /**
index b366d1cc0050d24720bec699e3e97de42f1871ce..6be207102a3334d4aff584af842eae1d6c2c1420 100644 (file)
@@ -122,7 +122,7 @@ private:
    * Helper function which returns the periodicity associated with this
    * coordinate system, according to dim, chartdim, and spacedim.
    */
-  static Point<spacedim> get_periodicity();
+  static Tensor<1,spacedim> get_periodicity();
 };
 
 
@@ -227,7 +227,7 @@ public:
    */
   FunctionManifold(const Function<chartdim> &push_forward_function,
                    const Function<spacedim> &pull_back_function,
-                   const Point<chartdim> periodicity=Point<chartdim>(),
+                   const Tensor<1,chartdim> &periodicity=Tensor<1,chartdim>(),
                    const double tolerance=1e-10);
 
   /**
@@ -246,7 +246,7 @@ public:
    */
   FunctionManifold(const std::string push_forward_expression,
                    const std::string pull_back_expression,
-                   const Point<chartdim> periodicity=Point<chartdim>(),
+                   const Tensor<1,chartdim> &periodicity=Tensor<1,chartdim>(),
                    const typename FunctionParser<spacedim>::ConstMap = typename FunctionParser<spacedim>::ConstMap(),
                    const std::string chart_vars=FunctionParser<chartdim>::default_variable_names(),
                    const std::string space_vars=FunctionParser<spacedim>::default_variable_names(),
index ae5e6d1ae04b27b1a527785d3c213801412bc7cb..bd9174d31fb8caec0777ca762ac870c355126062 100644 (file)
@@ -399,7 +399,7 @@ Manifold<dim,spacedim>::get_tangent_vector(const Point<spacedim> &x1,
 
 
 template <int dim, int spacedim>
-FlatManifold<dim,spacedim>::FlatManifold (const Point<spacedim> &periodicity,
+FlatManifold<dim,spacedim>::FlatManifold (const Tensor<1,spacedim> &periodicity,
                                           const double tolerance)
   :
   periodicity(periodicity),
@@ -417,16 +417,16 @@ get_new_point (const Quadrature<spacedim> &quad) const
   const std::vector<Point<spacedim> > &surrounding_points = quad.get_points();
   const std::vector<double> &weights = quad.get_weights();
 
-  Point<spacedim> minP = periodicity;
+  Tensor<1,spacedim> minP = periodicity;
 
   for (unsigned int d=0; d<spacedim; ++d)
     if (periodicity[d] > 0)
       for (unsigned int i=0; i<surrounding_points.size(); ++i)
         {
           minP[d] = std::min(minP[d], surrounding_points[i][d]);
-          Assert( (surrounding_points[i][d] < periodicity[d]+tolerance*periodicity[d]) ||
-                  (surrounding_points[i][d] >= -tolerance*periodicity[d]),
-                  ExcPeriodicBox(d, surrounding_points[i], periodicity, tolerance*periodicity[d]));
+          Assert( (surrounding_points[i][d] < periodicity[d]+tolerance*periodicity.norm()) ||
+                  (surrounding_points[i][d] >= -tolerance*periodicity.norm()),
+                  ExcPeriodicBox(d, surrounding_points[i], periodicity[i], tolerance*periodicity.norm()));
         }
 
   // compute the weighted average point, possibly taking into account periodicity
@@ -464,7 +464,7 @@ FlatManifold<dim, spacedim>::project_to_manifold (const std::vector<Point<spaced
 
 
 template <int dim, int spacedim>
-const Point<spacedim> &
+const Tensor<1,spacedim> &
 FlatManifold<dim, spacedim>::get_periodicity() const
 {
   return periodicity;
@@ -506,12 +506,13 @@ ChartManifold<dim,spacedim,chartdim>::~ChartManifold ()
 
 
 template <int dim, int spacedim, int chartdim>
-ChartManifold<dim,spacedim,chartdim>::ChartManifold (const Point<chartdim> &periodicity)
+ChartManifold<dim,spacedim,chartdim>::ChartManifold (const Tensor<1,chartdim> &periodicity)
   :
   sub_manifold(periodicity)
 {}
 
 
+
 template <int dim, int spacedim, int chartdim>
 Point<spacedim>
 ChartManifold<dim,spacedim,chartdim>::
@@ -544,6 +545,7 @@ push_forward_gradient(const Point<chartdim> &) const
 }
 
 
+
 template <int dim, int spacedim, int chartdim>
 Tensor<1,spacedim>
 ChartManifold<dim,spacedim,chartdim>::
@@ -561,10 +563,11 @@ get_tangent_vector (const Point<spacedim> &x1,
   return result;
 }
 
+
+
 template <int dim, int spacedim, int chartdim>
-const Point<chartdim> &
-ChartManifold<dim,spacedim,chartdim>::
-get_periodicity() const
+const Tensor<1, chartdim> &
+ChartManifold<dim, spacedim, chartdim>::get_periodicity() const
 {
   return sub_manifold.get_periodicity();
 }
index ab23791e94f6b82794ac90465b5cd494ef55680a..16ab8980c729f10fc334a1a59276336eea873b2c 100644 (file)
@@ -32,10 +32,10 @@ SphericalManifold<dim,spacedim>::SphericalManifold(const Point<spacedim> center)
 
 
 template <int dim, int spacedim>
-Point<spacedim>
+Tensor<1,spacedim>
 SphericalManifold<dim,spacedim>::get_periodicity()
 {
-  Point<spacedim> periodicity;
+  Tensor<1,spacedim> periodicity;
   periodicity[spacedim-1] = 2*numbers::PI; // theta and phi period.
   return periodicity;
 }
@@ -247,7 +247,7 @@ template <int dim, int spacedim, int chartdim>
 FunctionManifold<dim,spacedim,chartdim>::FunctionManifold
 (const Function<chartdim> &push_forward_function,
  const Function<spacedim> &pull_back_function,
- const Point<chartdim> periodicity,
+ const Tensor<1,chartdim> &periodicity,
  const double tolerance):
   ChartManifold<dim,spacedim,chartdim>(periodicity),
   push_forward_function(&push_forward_function),
@@ -263,7 +263,7 @@ template <int dim, int spacedim, int chartdim>
 FunctionManifold<dim,spacedim,chartdim>::FunctionManifold
 (const std::string push_forward_expression,
  const std::string pull_back_expression,
- const Point<chartdim> periodicity,
+ const Tensor<1,chartdim> &periodicity,
  const typename FunctionParser<spacedim>::ConstMap const_map,
  const std::string chart_vars,
  const std::string space_vars,
index dde1cf10f38902bdca4d77ef4a7142989edc7445..219bcc226ca4a0a02964dbd4646ebf25d4af1d96 100644 (file)
@@ -31,7 +31,7 @@ template <int dim, int spacedim>
 class MyFlatManifold : public ChartManifold<dim,spacedim,spacedim>
 {
 public:
-  MyFlatManifold (const Point<spacedim> &periodicity)
+  MyFlatManifold (const Tensor<1,spacedim> &periodicity)
     :
     ChartManifold<dim,spacedim,spacedim> (periodicity)
   {}
@@ -72,7 +72,7 @@ void test(unsigned int ref=1)
   deallog << "Testing dim=" << dim
           << ", spacedim="<< spacedim << std::endl;
 
-  Point<spacedim> periodicity;
+  Tensor<1,spacedim> periodicity;
   periodicity[0] = 5.0;
 
   MyFlatManifold<dim,spacedim> manifold(periodicity);
index fb46426f1fdcde18bc8a21b366bd580eee9d36d8..df653af9fe82e644c636ad2d8b1a1d6fae980057 100644 (file)
@@ -33,7 +33,7 @@ template <int dim, int spacedim>
 class MyFlatManifold : public ChartManifold<dim,spacedim,spacedim+1>
 {
 public:
-  MyFlatManifold (const Point<spacedim+1> &periodicity)
+  MyFlatManifold (const Tensor<1,spacedim+1> &periodicity)
     :
     ChartManifold<dim,spacedim,spacedim+1> (periodicity)
   {}
@@ -80,7 +80,7 @@ void test(unsigned int ref=1)
   deallog << "Testing dim=" << dim
           << ", spacedim="<< spacedim << std::endl;
 
-  Point<spacedim+1> periodicity;
+  Tensor<1,spacedim+1> periodicity;
   periodicity[0] = 5.0;
 
   MyFlatManifold<dim,spacedim> manifold(periodicity);
index fcd443268f4c6562bca64c44c06bb61f313842db..487dfabde8718dafaae5b66a869f09bd61dfac6b 100644 (file)
@@ -24,7 +24,7 @@ template <int dim, int spacedim>
 class MyFlatManifold : public ChartManifold<dim,spacedim,spacedim>
 {
 public:
-  MyFlatManifold (const Point<spacedim> &periodicity)
+  MyFlatManifold (const Tensor<1,spacedim> &periodicity)
     :
     ChartManifold<dim,spacedim,spacedim>(periodicity)
   {}
@@ -67,7 +67,7 @@ void test()
           << ", spacedim="<< spacedim << std::endl;
 
   // make the domain periodic in the first direction with periodicity 1.1
-  Point<spacedim> periodicity;
+  Tensor<1,spacedim> periodicity;
   periodicity[0] = 1.1;
   MyFlatManifold<dim,spacedim> manifold(periodicity);
 
index 9d127b3b29a43240a26afd1360b35ea3950b7672..96a3789d06bd0c2f5f2cfd5a785927272674270b 100644 (file)
@@ -26,7 +26,7 @@ template <int dim, int spacedim>
 class MyFlatManifold : public ChartManifold<dim,spacedim,spacedim+1>
 {
 public:
-  MyFlatManifold (const Point<spacedim+1> &periodicity)
+  MyFlatManifold (const Tensor<1,spacedim+1> &periodicity)
     :
     ChartManifold<dim,spacedim,spacedim+1>(periodicity)
   {}
@@ -74,7 +74,7 @@ void test()
           << ", spacedim="<< spacedim << std::endl;
 
   // make the domain periodic in the first direction with periodicity 1.1
-  Point<spacedim+1> periodicity;
+  Tensor<1,spacedim+1> periodicity;
   periodicity[0] = 1.1;
   MyFlatManifold<dim,spacedim> manifold(periodicity);
 
index cbf25744085e63c85bc7ecb42815dcf2c44e9c8a..3d0ca0383867cf69751ecd172d40b52de25f3d98 100644 (file)
@@ -19,9 +19,7 @@
 #include <deal.II/grid/manifold.h>
 
 
-Point<3> periodicity (/*r=*/0,
-                            /*phi=*/2*numbers::PI,
-                            /*z=*/0);
+Tensor<1,3> periodicity(static_cast<Tensor<1,3> >(Point<3>(0,2*numbers::PI,0)));
 
 class MyCylinderManifold : public ChartManifold<2,3,3>
 {
index 2621541c2a9105749775fee0bb6fa18443a15aff..946142ac5657253d6a69c74dcf4c0b062319d729 100644 (file)
@@ -31,7 +31,7 @@ void test(unsigned int ref=1)
   deallog << "Testing dim=" << dim
           << ", spacedim="<< spacedim << std::endl;
 
-  Point<spacedim> periodicity;
+  Tensor<1,spacedim> periodicity;
   periodicity[0] = 5.0;
 
   FlatManifold<dim,spacedim> manifold(periodicity);
index 6e7f8531690b0a507805de8ee24b2c00f070bdd7..770552420cbe2e9afba009e0f1e48f7b94533e36 100644 (file)
@@ -25,7 +25,7 @@ void test()
           << ", spacedim="<< spacedim << std::endl;
 
   // make the domain periodic in the first direction with periodicity 1.1
-  Point<spacedim> periodicity;
+  Tensor<1,spacedim> periodicity;
   periodicity[0] = 1.1;
   FlatManifold<dim,spacedim> manifold(periodicity);
 

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.