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.
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
*/
virtual
Point<spacedim> project_to_manifold (const Point<spacedim> &candidate) const;
+private:
+ const Point<spacedim> periodicity;
};
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.
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;
};
template <int spacedim>
-FlatManifold<spacedim>::FlatManifold ()
+FlatManifold<spacedim>::FlatManifold (const Point<spacedim> periodicity) :
+ periodicity(periodicity)
{}
template <int spacedim>
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);
}
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>
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]);
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__)
{
return 2.0*numbers::PI-angle;
}
-
-
template <>
Point<3>
TorusBoundary<2,3>::push_forward (const Point<2> &surfP) const
--- /dev/null
+//---------------------------- 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;
+}
+
--- /dev/null
+
+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
projects / dealii-svn.git / commitdiff