]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Added periodicity to FlatManifold
authorheltai <heltai@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 28 Jan 2014 10:17:57 +0000 (10:17 +0000)
committerheltai <heltai@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 28 Jan 2014 10:17:57 +0000 (10:17 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_manifold_id@32331 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/grid/manifold.h
deal.II/source/grid/manifold.cc
deal.II/source/grid/tria_boundary_lib.cc
tests/manifold/flat_manifold_01.cc [new file with mode: 0644]
tests/manifold/flat_manifold_01.output [new file with mode: 0644]

index ed85ca36786f4c879e96326b94bd7823c4901259..955962b69c840feaee62640181f6180d0405993e 100644 (file)
@@ -151,9 +151,22 @@ class FlatManifold: public Manifold<spacedim>
 public:
   /**
    * Default constructor. Some compilers require this for some
-   * reasons.
+   * reasons. The The optional argument can be used to specify the
+   * periodicity of the spacedim-dimensional manifold (one period per
+   * direction). A peridicity value of zero means that along that
+   * direction there is no peridicity. By default no periodicity is
+   * assumed.
+   *
+   * Periodicity affects the way a middle point is computed. It is
+   * assumed that if two points are more than half period distant,
+   * then the distance should be computed by crossing the periodicity
+   * boundary, i.e., the average is computed by adding a full period
+   * to the sum of the two. For example, if along direction 0 we have
+   * 2*pi periodicity, then the average of (2*pi-eps) and (eps) is not
+   * pi, but 2*pi (or zero), since, on a periodic manifold, these two
+   * points are at distance 2*eps and not (2*pi-eps)
    */
-  FlatManifold ();
+  FlatManifold (const Point<spacedim> periodicity=Point<spacedim>());
 
   /**
    * Let the new point be the average sum of surrounding vertices.
@@ -171,7 +184,7 @@ public:
     virtual Point<spacedim>
     get_new_point(const std::vector<Point<spacedim> > &surrounding_points,
                  const std::vector<double> &weights) const;
-
+  
   /**
    *  Project to FlatManifold: do nothing. Note however that this
    *  function can be overloaded by derived classes, which will then
@@ -181,6 +194,8 @@ public:
    */
   virtual
   Point<spacedim> project_to_manifold (const Point<spacedim> &candidate) const;
+private:
+  const Point<spacedim> periodicity;
 };
 
 
@@ -227,7 +242,24 @@ template <int dim, int spacedim>
 class ManifoldChart: public Manifold<spacedim>
 {
 public:
-
+  /**
+   * Constructor. The optional argument can be used to specify the
+   * periodicity of the dim-dimensional manifold (one period per
+   * direction). A peridicity value of zero means that along that
+   * direction there is no peridicity. By default no periodicity is
+   * assumed.
+   *
+   * Periodicity affects the way a middle point is computed. It is
+   * assumed that if two points are more than half period distant,
+   * then the distance should be computed by crossing the periodicity
+   * boundary, i.e., then the average is computed by adding a full
+   * period to the sum of the two. For example, if along direction 0
+   * we have 2*pi periodicity, then the average 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)
+   */
+  ManifoldChart(const Point<dim> periodicity=Point<dim>());
+  
   /**
    * Destructor. Does nothing here, but needs to be declared to make
    * it virtual.
@@ -264,10 +296,12 @@ public:
     virtual Point<spacedim>
     push_forward(const Point<dim> &chart_point) const = 0;
 
-
   private:
+  /**
+   * The sub_manifold object is used to compute the average of the
+   * points in the chart coordinates system. 
+   */
     FlatManifold<dim> sub_manifold;
-    mutable std::vector<Point<dim> > chart_points;
 };
 
 
index 1d598217a8bea86b0d87326cc905481b3aa8477c..60f6fb8f25e3d6b8cb407e77e604168221ca5abd 100644 (file)
@@ -70,7 +70,8 @@ Manifold<spacedim>::normal_vector(const Point<spacedim> &point) const
 
 
 template <int spacedim>
-FlatManifold<spacedim>::FlatManifold ()
+FlatManifold<spacedim>::FlatManifold (const Point<spacedim> periodicity) :
+  periodicity(periodicity)
 {}
 
 template <int spacedim>
@@ -90,8 +91,23 @@ get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
   
   
   Point<spacedim> p;
-  for(unsigned int i=0; i<surrounding_points.size(); ++i)
-    p += surrounding_points[i]*weights[i];
+  Point<spacedim> dp;
+  bool check_period = (periodicity.norm() != 0);
+  for(unsigned int i=0; i<surrounding_points.size(); ++i) {
+    if(check_period) {
+      for(unsigned int d=0; d<spacedim; ++d) 
+       if(periodicity[d] > 0)
+         dp[d] = std::abs(surrounding_points[i][d]-
+                          surrounding_points[0][d]) > periodicity[d]/2.0 ?
+           periodicity[d] : 0.0;
+    }
+    p += (surrounding_points[i]+dp)*weights[i];
+  }
+  if(check_period) 
+    for(unsigned int d=0; d<spacedim; ++d) 
+      if(periodicity[d] > 0) 
+       p[d] = std::fmod(p[d], periodicity[d]);
+
   return project_to_manifold(p);
 }
 
@@ -109,6 +125,11 @@ template <int dim, int spacedim>
 ManifoldChart<dim,spacedim>::~ManifoldChart ()
 {}
 
+template <int dim, int spacedim>
+ManifoldChart<dim,spacedim>::ManifoldChart (const Point<dim> periodicity):
+  sub_manifold(periodicity)
+{}
+
 
 template <int dim, int spacedim>
 Point<spacedim>
@@ -117,7 +138,7 @@ get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
               const std::vector<double> &weights) const
 {
   AssertDimension(surrounding_points.size(), weights.size());
-  chart_points.resize(surrounding_points.size());
+  std::vector<Point<dim> > chart_points(surrounding_points.size());
   
   for(unsigned int i=0; i<surrounding_points.size(); ++i) 
     chart_points[i] = pull_back(surrounding_points[i]);
index 61c86b0e8e6ca86cb55394a3bb092261482526e0..920e3e8840eea7948052aad6b1239809c6269d25 100644 (file)
@@ -420,6 +420,7 @@ template <>
 TorusBoundary<2,3>::TorusBoundary (const double R__,
                                    const double r__)
   :
+  ManifoldChart<2,3>(Point<2>(2*numbers::PI, 2*numbers::PI)),
   R(R__),
   r(r__)
 {
@@ -448,8 +449,6 @@ TorusBoundary<dim,spacedim>::get_correct_angle(const double angle,
   return 2.0*numbers::PI-angle;
 }
 
-
-
 template <>
 Point<3>
 TorusBoundary<2,3>::push_forward (const Point<2> &surfP) const
diff --git a/tests/manifold/flat_manifold_01.cc b/tests/manifold/flat_manifold_01.cc
new file mode 100644 (file)
index 0000000..9d463df
--- /dev/null
@@ -0,0 +1,55 @@
+//----------------------------  manifold_id_01.cc  ---------------------------
+//    Copyright (C) 2011, 2013 by the mathLab team.
+//
+//    This file is subject to LGPL 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.
+//
+//----------------------------  flat_manifold_01.cc  ---------------------------
+
+
+// Test periodicity conditions using FlatManifold
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+
+// all include files you need here
+#include <deal.II/grid/manifold.h>
+#include <deal.II/grid/manifold_lib.h>
+
+template <int spacedim> 
+void test() {
+  Point<spacedim> period;
+  for(unsigned int i=0; i<spacedim; ++i)
+    period[i] = 1.0;
+  FlatManifold<spacedim> manifold(period);
+  deallog << "Testing spacedim = " << spacedim << std::endl;
+  deallog << "Period: " << period << std::endl;
+  std::vector<Point<spacedim> > p(2);
+  std::vector<double> w(2, 0.5);
+  double eps = .1;
+  p[0] = period*eps*2;
+  p[1] = period*(1.0-eps);
+  Point<spacedim> middle = manifold.get_new_point(p, w);
+  deallog << "P(0): " << p[0] << std::endl;
+  deallog << "P(1): " << p[1] << std::endl;
+  deallog << "Middle: " << middle << std::endl;
+}
+
+int main () 
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-5);
+  
+  test<1>();
+  test<2>();
+  test<3>();
+
+  return 0;
+}
+                  
diff --git a/tests/manifold/flat_manifold_01.output b/tests/manifold/flat_manifold_01.output
new file mode 100644 (file)
index 0000000..f6f885f
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL::Testing spacedim = 1
+DEAL::Period: 1.00000
+DEAL::P(0): 0.200000
+DEAL::P(1): 0.900000
+DEAL::Middle: 0.0500000
+DEAL::Testing spacedim = 2
+DEAL::Period: 1.00000 1.00000
+DEAL::P(0): 0.200000 0.200000
+DEAL::P(1): 0.900000 0.900000
+DEAL::Middle: 0.0500000 0.0500000
+DEAL::Testing spacedim = 3
+DEAL::Period: 1.00000 1.00000 1.00000
+DEAL::P(0): 0.200000 0.200000 0.200000
+DEAL::P(1): 0.900000 0.900000 0.900000
+DEAL::Middle: 0.0500000 0.0500000 0.0500000

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.