]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Moved CompositionManifold to its own directory.
authorLuca Heltai <luca.heltai@sissa.it>
Sat, 9 Apr 2016 12:23:54 +0000 (14:23 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Mon, 11 Apr 2016 07:41:41 +0000 (09:41 +0200)
include/deal.II/grid/composition_manifold.h [new file with mode: 0644]
include/deal.II/grid/manifold.h
source/grid/CMakeLists.txt
source/grid/manifold.cc
source/grid/manifold_tools.cc [deleted file]
source/grid/manifold_tools.inst.in [deleted file]
tests/manifold/composition_manifold_01.cc
tests/manifold/composition_manifold_02.cc
tests/manifold/composition_manifold_03.cc
tests/manifold/composition_manifold_04.cc

diff --git a/include/deal.II/grid/composition_manifold.h b/include/deal.II/grid/composition_manifold.h
new file mode 100644 (file)
index 0000000..c9809cd
--- /dev/null
@@ -0,0 +1,183 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii__composition_manifold_h
+#define dealii__composition_manifold_h
+
+
+/*----------------------------   composition_manifold.h     ------------*/
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/subscriptor.h>
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/thread_management.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/derivative_form.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/manifold.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * CompositionManifold.  Take two ChartManifold objects, and make
+ * their composition. The CompositionManifold object is a
+ * ChartManifold going from the chart of the first ChartManifold to
+ * the embedding space of the second ChartManifold. If the first
+ * ChartManifold is periodic, so is the resulting ChartManifold, with
+ * the same periodicity. Periodicity on the second ChartManifold is not
+ * allowed, and the constructor will throw an axception if the second
+ * Manifold is periodic.
+ *
+ * This class only works for dim <= chartdim <= intermediate_spacedim
+ * <= spacedim. If you try to instantiate anything different, an
+ * Exception will be thrown in one of the ChartManifold classes that
+ * violates this condition.
+ *
+ * Given the ChartManifold F and the ChartManifold G, this class
+ * represents the composition of G after F.
+ *
+ * The template parameters have the following meaning:
+ *
+ * @tparam dim The dimension of the resulting ChartManifold
+ * @tparam spacedim The space dimension of the resulting ChartManifold
+ * @tparam chartdim The chart dimension of the resulting ChartManifold
+ * @tparam intermediate_dim The space dimension of the first ChartManifold
+ * @tparam dim1 The dimension of the first ChartManifold, which coincides also
+ * with the chart dimension of the second ChartManifold
+ * @tparam dim2 The dimension of the second ChartManifold
+ *
+ * @ingroup manifold
+ * @author Luca Heltai, Timo Heister, 2016
+ */
+template <int dim, int spacedim=dim, int chartdim=dim,
+          int intermediate_dim=dim, int dim1=dim, int dim2=dim>
+class CompositionManifold : public ChartManifold<dim, spacedim, chartdim>
+{
+public:
+
+  /**
+   * Construct the composition of the two given manifolds.
+   */
+  CompositionManifold(const ChartManifold<dim1, intermediate_dim, chartdim> &F,
+                      const ChartManifold<dim2, spacedim, intermediate_dim> &G);
+
+  /**
+   * Pull back the given point in spacedim to the Euclidean chartdim
+   * dimensional space. This function calls the pull_back() function
+   * of G, and then the pull_back() function of F.
+   */
+  virtual
+  Point<chartdim>
+  pull_back(const Point<spacedim> &space_point) const;
+
+
+  /**
+   * Push forward the chartdim dimensional point to a spacedim
+   * Euclidean point. The function calls first the push_forward() of
+   * F, and then the push_foward() of G.
+   */
+  virtual
+  Point<spacedim>
+  push_forward(const Point<chartdim> &chart_point) const;
+
+  /**
+   * Return the derivative of the composition of G after F.
+   */
+  virtual
+  DerivativeForm<1,chartdim,spacedim>
+  push_forward_gradient(const Point<chartdim> &chart_point) const;
+
+private:
+  /**
+   * The first ChartManifold.
+   */
+  SmartPointer<const ChartManifold<dim1, intermediate_dim, chartdim>,
+               CompositionManifold<dim,spacedim,chartdim,dim1,dim2,intermediate_dim> > F;
+
+
+  /**
+   * The second ChartManifold.
+   */
+  SmartPointer<const ChartManifold<dim2, spacedim, intermediate_dim>,
+               CompositionManifold<dim,spacedim,chartdim,dim1,dim2,intermediate_dim> > G;
+};
+
+
+/*------------------Template Implementations------------------------*/
+
+template <int dim, int spacedim, int chartdim, int intermediate_dim,
+          int dim1, int dim2>
+CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::CompositionManifold
+(const ChartManifold<dim1, intermediate_dim, chartdim> &F,
+ const ChartManifold<dim2, spacedim, intermediate_dim> &G) :
+  ChartManifold<dim,spacedim,chartdim>(F.get_periodicity()),
+  F(&F),
+  G(&G)
+{
+  // We don't know what to do with a periodicity in the second manifold, so
+  // throw an assertion if the second manifold is periodic
+  Assert(G.get_periodicity().norm() == 0.0,
+         ExcMessage("The second manifold cannot be periodic."));
+}
+
+
+
+template <int dim, int spacedim, int chartdim, int intermediate_dim,
+          int dim1, int dim2>
+Point<chartdim>
+CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::pull_back(const Point<spacedim> &space_point) const
+{
+  return F->pull_back(G->pull_back(space_point));
+}
+
+
+
+template <int dim, int spacedim, int chartdim, int intermediate_dim,
+          int dim1, int dim2>
+Point<spacedim>
+CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::push_forward(const Point<chartdim> &chart_point) const
+{
+  return G->push_forward(F->push_forward(chart_point));
+}
+
+
+
+template <int dim, int spacedim, int chartdim, int intermediate_dim,
+          int dim1, int dim2>
+DerivativeForm<1,chartdim,spacedim>
+CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::push_forward_gradient(const Point<chartdim> &chart_point) const
+{
+  DerivativeForm<1,chartdim,intermediate_dim> DF =
+    F->push_forward_gradient(chart_point);
+
+  DerivativeForm<1,intermediate_dim,spacedim> DG =
+    G->push_forward_gradient(F->push_forward(chart_point));
+
+  DerivativeForm<1,chartdim,spacedim> DF_DG;
+
+  for (unsigned int d=0; d<spacedim; ++d)
+    for (unsigned int c=0; c<chartdim; ++c)
+      for (unsigned int s=0; s<intermediate_dim; ++s)
+        DF_DG[d][c] += DG[d][s]*DF[s][c];
+
+  return DF_DG;
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index b76e8be958930174a9c085bbfc1789dfb4e3a0f4..c82e737dc7f2504203d45c3c17bd81bf3eb2de59 100644 (file)
@@ -617,10 +617,9 @@ private:
    */
   const Tensor<1,spacedim> periodicity;
 
-  DeclException4(ExcPeriodicBox, int, Point<spacedim>, double, double,
+  DeclException3(ExcPeriodicBox, int, Point<spacedim>, double,
                  << "The component number " << arg1 << " of the point [ " << arg2
-                 << " ] is not in the interval [ " << -arg4
-                 << ", " << arg3 << "), bailing out.");
+                 << " ] is not in the interval [ 0, " << arg3 << "), bailing out.");
 
   /**
    * Relative tolerance. This tolerance is used to compute distances in double
index 18abd372114d0f7b99f2e88fa45cc820442c5a77..8279350bf27ab47b9fd7994a23b8f81f983d743e 100644 (file)
@@ -26,7 +26,6 @@ SET(_src
   intergrid_map.cc
   manifold.cc
   manifold_lib.cc
-  manifold_tools.cc
   persistent_tria.cc
   tria_accessor.cc
   tria_boundary.cc
@@ -46,7 +45,6 @@ SET(_inst
   intergrid_map.inst.in
   manifold.inst.in
   manifold_lib.inst.in
-  manifold_tools.inst.in
   tria_accessor.inst.in
   tria_boundary.inst.in
   tria_boundary_lib.inst.in
index bd9174d31fb8caec0777ca762ac870c355126062..07ff568521122682292e8bf3919e7dc4d9eed6b3 100644 (file)
@@ -424,9 +424,9 @@ get_new_point (const Quadrature<spacedim> &quad) const
       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.norm()) ||
-                  (surrounding_points[i][d] >= -tolerance*periodicity.norm()),
-                  ExcPeriodicBox(d, surrounding_points[i], periodicity[i], tolerance*periodicity.norm()));
+          Assert( (surrounding_points[i][d] < periodicity[d]+tolerance*periodicity[d]) ||
+                  (surrounding_points[i][d] >= -tolerance*periodicity[d]),
+                  ExcPeriodicBox(d, surrounding_points[i], periodicity[i]));
         }
 
   // compute the weighted average point, possibly taking into account periodicity
diff --git a/source/grid/manifold_tools.cc b/source/grid/manifold_tools.cc
deleted file mode 100644 (file)
index ac98682..0000000
+++ /dev/null
@@ -1,93 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 1998 - 2016 by the deal.II authors
-//
-// This file is part of the deal.II library.
-//
-// The deal.II library is free software; you can use it, redistribute
-// it, and/or modify it under the terms of the GNU Lesser General
-// Public License as published by the Free Software Foundation; either
-// version 2.1 of the License, or (at your option) any later version.
-// The full text of the license can be found in the file LICENSE at
-// the top level of the deal.II distribution.
-//
-// ---------------------------------------------------------------------
-
-#include <deal.II/base/tensor.h>
-#include <deal.II/grid/tria.h>
-#include <deal.II/grid/tria_iterator.h>
-#include <deal.II/grid/tria_accessor.h>
-#include <deal.II/grid/manifold_tools.h>
-#include <deal.II/fe/fe_q.h>
-#include <cmath>
-
-DEAL_II_NAMESPACE_OPEN
-
-
-/* -------------------------- CompositionManifold --------------------- */
-
-template <int dim, int spacedim, int chartdim, int intermediate_dim,
-          int dim1, int dim2>
-CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::CompositionManifold
-(const ChartManifold<dim1, intermediate_dim, chartdim> &F,
- const ChartManifold<dim2, spacedim, intermediate_dim> &G) :
-  ChartManifold<dim,spacedim,chartdim>(F.get_periodicity()),
-  F(&F),
-  G(&G)
-{
-  // We don't know what to do with a periodicity in the second manifold, so
-  // throw an assertion if the second manifold is periodic
-  Assert(G.get_periodicity().norm() == 0.0,
-         ExcMessage("The second manifold cannot be periodic."));
-}
-
-
-
-template <int dim, int spacedim, int chartdim, int intermediate_dim,
-          int dim1, int dim2>
-Point<chartdim>
-CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::pull_back(const Point<spacedim> &space_point) const
-{
-  return F->pull_back(G->pull_back(space_point));
-}
-
-
-
-template <int dim, int spacedim, int chartdim, int intermediate_dim,
-          int dim1, int dim2>
-Point<spacedim>
-CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::push_forward(const Point<chartdim> &chart_point) const
-{
-  return G->push_forward(F->push_forward(chart_point));
-}
-
-
-
-template <int dim, int spacedim, int chartdim, int intermediate_dim,
-          int dim1, int dim2>
-DerivativeForm<1,chartdim,spacedim>
-CompositionManifold<dim,spacedim,chartdim,intermediate_dim,dim1,dim2>::push_forward_gradient(const Point<chartdim> &chart_point) const
-{
-  DerivativeForm<1,chartdim,intermediate_dim> DF =
-    F->push_forward_gradient(chart_point);
-
-  DerivativeForm<1,intermediate_dim,spacedim> DG =
-    G->push_forward_gradient(F->push_forward(chart_point));
-
-  DerivativeForm<1,chartdim,spacedim> DF_DG;
-
-  for (unsigned int d=0; d<spacedim; ++d)
-    for (unsigned int c=0; c<chartdim; ++c)
-      for (unsigned int s=0; s<intermediate_dim; ++s)
-        DF_DG[d][c] += DG[d][s]*DF[s][c];
-
-  return DF_DG;
-}
-
-
-
-// explicit instantiations
-#include "manifold_tools.inst"
-
-DEAL_II_NAMESPACE_CLOSE
-
diff --git a/source/grid/manifold_tools.inst.in b/source/grid/manifold_tools.inst.in
deleted file mode 100644 (file)
index 54684d5..0000000
+++ /dev/null
@@ -1,24 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 1999 - 2014 by the deal.II authors
-//
-// This file is part of the deal.II library.
-//
-// The deal.II library is free software; you can use it, redistribute
-// it, and/or modify it under the terms of the GNU Lesser General
-// Public License as published by the Free Software Foundation; either
-// version 2.1 of the License, or (at your option) any later version.
-// The full text of the license can be found in the file LICENSE at
-// the top level of the deal.II distribution.
-//
-// ---------------------------------------------------------------------
-
-
-for (d : DIMENSIONS; d1 : DIMENSIONS; d2: DIMENSIONS; 
-     c : SPACE_DIMENSIONS ; s : SPACE_DIMENSIONS; s1 : SPACE_DIMENSIONS)
-  {
-#if d <= s && c <= s1 && s1 <= s && d1 <= s1 && d2 <= s
-    template class CompositionManifold<d,s,c,s1,d1,d2>;
-#endif
-  }
-
index 19af4633b4c8ba15fbed500feddf10043933ae4a..433b7291ae1b0e92d9c31b4471c43ae7996977b6 100644 (file)
@@ -18,7 +18,7 @@
 
 // all include files you need here
 #include <deal.II/grid/manifold_lib.h>
-#include <deal.II/grid/manifold_tools.h>
+#include <deal.II/grid/composition_manifold.h>
 
 
 int main ()
index 38cfb6c18ae38286db7b7fe501cced0799a9325e..fc9626184532120c51e5add75d2414f1063ebfe8 100644 (file)
@@ -19,7 +19,7 @@
 
 // all include files you need here
 #include <deal.II/grid/manifold_lib.h>
-#include <deal.II/grid/manifold_tools.h>
+#include <deal.II/grid/composition_manifold.h>
 
 
 int main ()
index 38bb5945e8892f36abe1440055973b4983636ef4..08eeef5609aaa803c5cc87508da9f91c2a6c4702 100644 (file)
@@ -19,7 +19,7 @@
 
 // all include files you need here
 #include <deal.II/grid/manifold_lib.h>
-#include <deal.II/grid/manifold_tools.h>
+#include <deal.II/grid/composition_manifold.h>
 
 
 int main ()
index 107a5a2d55d64f86f67df53835f8e5fded106955..c0d3c33b63c59a835fa58b5baccde4ed91df90ae 100644 (file)
@@ -19,7 +19,7 @@
 
 // all include files you need here
 #include <deal.II/grid/manifold_lib.h>
-#include <deal.II/grid/manifold_tools.h>
+#include <deal.II/grid/composition_manifold.h>
 
 
 int main ()

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.