]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Generalize `GridTools::rotate` for arbitrary rotation axes.
authorMarc Fehling <mafehling.git@gmail.com>
Tue, 12 Oct 2021 02:06:19 +0000 (20:06 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Tue, 19 Oct 2021 15:58:14 +0000 (09:58 -0600)
13 files changed:
doc/news/changes/incompatibilities/20211011Fehling [new file with mode: 0644]
doc/news/changes/minor/20211011Fehling [new file with mode: 0644]
include/deal.II/grid/grid_tools.h
source/grid/grid_generator.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
source/grid/grid_tools_nontemplates.cc
tests/fe/bdm_2.cc
tests/fe/nedelec.cc
tests/fe/rt_2.cc
tests/fe/rt_bubbles_2.cc
tests/grid/grid_tools.cc
tests/grid/rotate_01.cc

diff --git a/doc/news/changes/incompatibilities/20211011Fehling b/doc/news/changes/incompatibilities/20211011Fehling
new file mode 100644 (file)
index 0000000..7105a29
--- /dev/null
@@ -0,0 +1,5 @@
+Changed: The 3D implementation of GridTools::rotate() with an integer
+for a Cartesian coordinate direction has been superseded by the version
+that accepts unit vectors as rotation axes.
+<br>
+(Marc Fehling, 2021/10/11)
diff --git a/doc/news/changes/minor/20211011Fehling b/doc/news/changes/minor/20211011Fehling
new file mode 100644 (file)
index 0000000..020002b
--- /dev/null
@@ -0,0 +1,4 @@
+Changed: The 3D implementation of GridTools::rotate() now accepts unit
+vectors as rotation axes.
+<br>
+(Marc Fehling, 2021/10/11)
index b639f59453f3739bd35cd9a41bd795dd91a8f8b3..d43eb0435b332cc4e72564efc4838204b3fa06d1 100644 (file)
 // ---------------------------------------------------------------------
 
 #ifndef dealii_grid_tools_h
-#  define dealii_grid_tools_h
+#define dealii_grid_tools_h
 
 
-#  include <deal.II/base/config.h>
+#include <deal.II/base/config.h>
 
-#  include <deal.II/base/bounding_box.h>
-#  include <deal.II/base/geometry_info.h>
-#  include <deal.II/base/std_cxx17/optional.h>
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/geometry_info.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/std_cxx17/optional.h>
 
-#  include <deal.II/boost_adaptors/bounding_box.h>
+#include <deal.II/boost_adaptors/bounding_box.h>
 
-#  include <deal.II/distributed/shared_tria.h>
+#include <deal.II/distributed/shared_tria.h>
 
-#  include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_handler.h>
 
-#  include <deal.II/fe/fe_values.h>
-#  include <deal.II/fe/mapping.h>
-#  include <deal.II/fe/mapping_q1.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping.h>
+#include <deal.II/fe/mapping_q1.h>
 
-#  include <deal.II/grid/manifold.h>
-#  include <deal.II/grid/tria.h>
-#  include <deal.II/grid/tria_accessor.h>
-#  include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/manifold.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
 
-#  include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/dof_handler.h>
 
-#  include <deal.II/lac/la_parallel_vector.h>
-#  include <deal.II/lac/la_vector.h>
-#  include <deal.II/lac/petsc_vector.h>
-#  include <deal.II/lac/sparsity_tools.h>
-#  include <deal.II/lac/trilinos_vector.h>
+#include <deal.II/lac/la_parallel_vector.h>
+#include <deal.II/lac/la_vector.h>
+#include <deal.II/lac/petsc_vector.h>
+#include <deal.II/lac/sparsity_tools.h>
+#include <deal.II/lac/trilinos_vector.h>
 
-#  include <deal.II/numerics/rtree.h>
+#include <deal.II/numerics/rtree.h>
 
 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
-#  include <boost/archive/binary_iarchive.hpp>
-#  include <boost/archive/binary_oarchive.hpp>
-#  include <boost/geometry/index/rtree.hpp>
-#  include <boost/random/mersenne_twister.hpp>
-#  include <boost/serialization/array.hpp>
-#  include <boost/serialization/vector.hpp>
-
-#  ifdef DEAL_II_WITH_ZLIB
-#    include <boost/iostreams/device/back_inserter.hpp>
-#    include <boost/iostreams/filter/gzip.hpp>
-#    include <boost/iostreams/filtering_stream.hpp>
-#    include <boost/iostreams/stream.hpp>
-#  endif
+#include <boost/archive/binary_iarchive.hpp>
+#include <boost/archive/binary_oarchive.hpp>
+#include <boost/geometry/index/rtree.hpp>
+#include <boost/random/mersenne_twister.hpp>
+#include <boost/serialization/array.hpp>
+#include <boost/serialization/vector.hpp>
+
+#ifdef DEAL_II_WITH_ZLIB
+#  include <boost/iostreams/device/back_inserter.hpp>
+#  include <boost/iostreams/filter/gzip.hpp>
+#  include <boost/iostreams/filtering_stream.hpp>
+#  include <boost/iostreams/stream.hpp>
+#endif
 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
 
-#  include <bitset>
-#  include <list>
-#  include <set>
+#include <bitset>
+#include <list>
+#include <set>
 
 DEAL_II_NAMESPACE_OPEN
 
 // Forward declarations
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
 namespace parallel
 {
   namespace distributed
@@ -88,7 +89,7 @@ namespace hp
 }
 
 class SparsityPattern;
-#  endif
+#endif
 
 namespace internal
 {
@@ -96,14 +97,14 @@ namespace internal
   class ActiveCellIterator
   {
   public:
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
     using type = typename MeshType::active_cell_iterator;
-#  else
+#else
     using type = TriaActiveIterator<dealii::CellAccessor<dim, spacedim>>;
-#  endif
+#endif
   };
 
-#  ifdef _MSC_VER
+#ifdef _MSC_VER
   template <int dim, int spacedim>
   class ActiveCellIterator<dim, spacedim, dealii::DoFHandler<dim, spacedim>>
   {
@@ -111,7 +112,7 @@ namespace internal
     using type =
       TriaActiveIterator<dealii::DoFCellAccessor<dim, spacedim, false>>;
   };
-#  endif
+#endif
 } // namespace internal
 
 /**
@@ -555,6 +556,23 @@ namespace GridTools
   void
   rotate(const double angle, Triangulation<dim> &triangulation);
 
+  /**
+   * Rotate all vertices of the given @p triangulation in counter-clockwise
+   * direction around the @p axis described by a unit vector. Otherwise like the
+   * function above.
+   *
+   * @param[in] angle Angle in radians to rotate the Triangulation by.
+   * @param[in] axis A unit vector that defines the axis of rotation.
+   * @param[in,out] triangulation The Triangulation object to rotate.
+   *
+   * @note Implemented for spacedim=3 and dim=1, 2, and 3.
+   */
+  template <int dim>
+  void
+  rotate(const double            angle,
+         const Point<3, double> &axis,
+         Triangulation<dim, 3> & triangulation);
+
   /**
    * Rotate all vertices of the given @p triangulation in counter-clockwise
    * direction around the axis with the given index. Otherwise like the
@@ -566,9 +584,11 @@ namespace GridTools
    * @param[in,out] triangulation The Triangulation object to rotate.
    *
    * @note Implemented for dim=1, 2, and 3.
+   *
+   * @deprecated Use the alternative with the unit vector instead.
    */
   template <int dim>
-  void
+  DEAL_II_DEPRECATED_EARLY void
   rotate(const double           angle,
          const unsigned int     axis,
          Triangulation<dim, 3> &triangulation);
@@ -904,14 +924,14 @@ namespace GridTools
    * the performance if the function is called only once on few points.
    */
   template <int dim, int spacedim>
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
   std::tuple<
     std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>,
     std::vector<std::vector<Point<dim>>>,
     std::vector<std::vector<unsigned int>>>
-#  else
+#else
   return_type
-#  endif
+#endif
   compute_point_locations(
     const Cache<dim, spacedim> &        cache,
     const std::vector<Point<spacedim>> &points,
@@ -953,15 +973,15 @@ namespace GridTools
    * GridTools::compute_point_locations().
    */
   template <int dim, int spacedim>
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
   std::tuple<
     std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>,
     std::vector<std::vector<Point<dim>>>,
     std::vector<std::vector<unsigned int>>,
     std::vector<unsigned int>>
-#  else
+#else
   return_type
-#  endif
+#endif
   compute_point_locations_try_all(
     const Cache<dim, spacedim> &        cache,
     const std::vector<Point<spacedim>> &points,
@@ -1039,16 +1059,16 @@ namespace GridTools
    * of this page.
    */
   template <int dim, int spacedim>
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
   std::tuple<
     std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>,
     std::vector<std::vector<Point<dim>>>,
     std::vector<std::vector<unsigned int>>,
     std::vector<std::vector<Point<spacedim>>>,
     std::vector<std::vector<unsigned int>>>
-#  else
+#else
   return_type
-#  endif
+#endif
   distributed_compute_point_locations(
     const GridTools::Cache<dim, spacedim> &                cache,
     const std::vector<Point<spacedim>> &                   local_points,
@@ -1285,13 +1305,13 @@ namespace GridTools
    * for this case.
    */
   template <int dim, template <int, int> class MeshType, int spacedim>
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
   std::vector<typename MeshType<dim, spacedim>::active_cell_iterator>
-#  else
+#else
   std::vector<
     typename dealii::internal::
       ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type>
-#  endif
+#endif
   find_cells_adjacent_to_vertex(const MeshType<dim, spacedim> &container,
                                 const unsigned int             vertex_index);
 
@@ -1358,13 +1378,13 @@ namespace GridTools
    * processor will return a locally owned cell and the other one a ghost cell.
    */
   template <int dim, template <int, int> class MeshType, int spacedim>
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
   std::pair<typename MeshType<dim, spacedim>::active_cell_iterator, Point<dim>>
-#  else
+#else
   std::pair<typename dealii::internal::
               ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
             Point<dim>>
-#  endif
+#endif
   find_active_cell_around_point(const Mapping<dim, spacedim> & mapping,
                                 const MeshType<dim, spacedim> &mesh,
                                 const Point<spacedim> &        p,
@@ -1379,12 +1399,12 @@ namespace GridTools
    * @return An iterator into the mesh that points to the surrounding cell.
    */
   template <int dim, template <int, int> class MeshType, int spacedim>
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
   typename MeshType<dim, spacedim>::active_cell_iterator
-#  else
+#else
   typename dealii::internal::
     ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type
-#  endif
+#endif
   find_active_cell_around_point(const MeshType<dim, spacedim> &mesh,
                                 const Point<spacedim> &        p,
                                 const std::vector<bool> &marked_vertices = {},
@@ -1481,13 +1501,13 @@ namespace GridTools
    * call the function above with argument `cache` in this case.
    */
   template <int dim, template <int, int> class MeshType, int spacedim>
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
   std::pair<typename MeshType<dim, spacedim>::active_cell_iterator, Point<dim>>
-#  else
+#else
   std::pair<typename dealii::internal::
               ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
             Point<dim>>
-#  endif
+#endif
   find_active_cell_around_point(
     const Mapping<dim, spacedim> & mapping,
     const MeshType<dim, spacedim> &mesh,
@@ -1529,15 +1549,15 @@ namespace GridTools
    * @endcode
    */
   template <int dim, template <int, int> class MeshType, int spacedim>
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
   std::vector<std::pair<typename MeshType<dim, spacedim>::active_cell_iterator,
                         Point<dim>>>
-#  else
+#else
   std::vector<std::pair<
     typename dealii::internal::
       ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
     Point<dim>>>
-#  endif
+#endif
   find_all_active_cells_around_point(
     const Mapping<dim, spacedim> & mapping,
     const MeshType<dim, spacedim> &mesh,
@@ -1553,15 +1573,15 @@ namespace GridTools
    * function find_all_active_cells_around_point() above.
    */
   template <int dim, template <int, int> class MeshType, int spacedim>
-#  ifndef _MSC_VER
+#ifndef _MSC_VER
   std::vector<std::pair<typename MeshType<dim, spacedim>::active_cell_iterator,
                         Point<dim>>>
-#  else
+#else
   std::vector<std::pair<
     typename dealii::internal::
       ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
     Point<dim>>>
-#  endif
+#endif
   find_all_active_cells_around_point(
     const Mapping<dim, spacedim> & mapping,
     const MeshType<dim, spacedim> &mesh,
@@ -1908,13 +1928,13 @@ namespace GridTools
    * of this page.
    */
   template <int spacedim>
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
   std::tuple<std::vector<std::vector<unsigned int>>,
              std::map<unsigned int, unsigned int>,
              std::map<unsigned int, std::vector<unsigned int>>>
-#  else
+#else
   return_type
-#  endif
+#endif
   guess_point_owner(
     const std::vector<std::vector<BoundingBox<spacedim>>> &global_bboxes,
     const std::vector<Point<spacedim>> &                   points);
@@ -1955,13 +1975,13 @@ namespace GridTools
    * of this page.
    */
   template <int spacedim>
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
   std::tuple<std::map<unsigned int, std::vector<unsigned int>>,
              std::map<unsigned int, unsigned int>,
              std::map<unsigned int, std::vector<unsigned int>>>
-#  else
+#else
   return_type
-#  endif
+#endif
   guess_point_owner(
     const RTree<std::pair<BoundingBox<spacedim>, unsigned int>> &covering_rtree,
     const std::vector<Point<spacedim>> &                         points);
@@ -3254,7 +3274,7 @@ namespace GridTools
     void
     load(Archive &ar, const unsigned int version);
 
-#  ifdef DOXYGEN
+#ifdef DOXYGEN
     /**
      * Read or write the data of this object to or from a stream for the
      * purpose of serialization using the [BOOST serialization
@@ -3263,11 +3283,11 @@ namespace GridTools
     template <class Archive>
     void
     serialize(Archive &archive, const unsigned int version);
-#  else
+#else
     // This macro defines the serialize() method that is compatible with
     // the templated save() and load() method that have been implemented.
     BOOST_SERIALIZATION_SPLIT_MEMBER()
-#  endif
+#endif
   };
 
 
@@ -3455,7 +3475,7 @@ DeclExceptionMsg(ExcMeshNotOrientable,
 
 /* ----------------- Template function --------------- */
 
-#  ifndef DOXYGEN
+#ifndef DOXYGEN
 
 namespace GridTools
 {
@@ -4260,7 +4280,7 @@ namespace GridTools
                                           MeshType::space_dimension> &)>
         &compute_ghost_owners)
     {
-#    ifndef DEAL_II_WITH_MPI
+#  ifndef DEAL_II_WITH_MPI
       (void)mesh;
       (void)pack;
       (void)unpack;
@@ -4268,7 +4288,7 @@ namespace GridTools
       (void)process_cells;
       (void)compute_ghost_owners;
       Assert(false, ExcNeedsMPI());
-#    else
+#  else
       constexpr int dim      = MeshType::dimension;
       constexpr int spacedim = MeshType::space_dimension;
       auto          tria =
@@ -4490,7 +4510,7 @@ namespace GridTools
         }
 
 
-#    endif // DEAL_II_WITH_MPI
+#  endif // DEAL_II_WITH_MPI
     }
 
   } // namespace internal
@@ -4506,13 +4526,13 @@ namespace GridTools
     const std::function<bool(const typename MeshType::active_cell_iterator &)>
       &cell_filter)
   {
-#    ifndef DEAL_II_WITH_MPI
+#  ifndef DEAL_II_WITH_MPI
     (void)mesh;
     (void)pack;
     (void)unpack;
     (void)cell_filter;
     Assert(false, ExcNeedsMPI());
-#    else
+#  else
     internal::exchange_cell_data<DataType,
                                  MeshType,
                                  typename MeshType::active_cell_iterator>(
@@ -4526,7 +4546,7 @@ namespace GridTools
             process(cell, cell->subdomain_id());
       },
       [](const auto &tria) { return tria.ghost_owners(); });
-#    endif
+#  endif
   }
 
 
@@ -4542,13 +4562,13 @@ namespace GridTools
     const std::function<bool(const typename MeshType::level_cell_iterator &)>
       &cell_filter)
   {
-#    ifndef DEAL_II_WITH_MPI
+#  ifndef DEAL_II_WITH_MPI
     (void)mesh;
     (void)pack;
     (void)unpack;
     (void)cell_filter;
     Assert(false, ExcNeedsMPI());
-#    else
+#  else
     internal::exchange_cell_data<DataType,
                                  MeshType,
                                  typename MeshType::level_cell_iterator>(
@@ -4564,15 +4584,12 @@ namespace GridTools
             process(cell, cell->level_subdomain_id());
       },
       [](const auto &tria) { return tria.level_ghost_owners(); });
-#    endif
+#  endif
   }
 } // namespace GridTools
 
-#  endif
+#endif // DOXYGEN
 
 DEAL_II_NAMESPACE_CLOSE
 
-/*----------------------------   grid_tools.h     ---------------------------*/
-/* end of #ifndef dealii_grid_tools_h */
 #endif
-/*----------------------------   grid_tools.h     ---------------------------*/
index 814300b3b1d56d8e549f055687ee086a59e6e896..ff9868ae4411c9c78fe48ab468434470fb778f39 100644 (file)
@@ -4940,7 +4940,7 @@ namespace GridGenerator
                                          n_slices,
                                          2 * half_length,
                                          triangulation);
-    GridTools::rotate(numbers::PI / 2, 1, triangulation);
+    GridTools::rotate(numbers::PI_2, Point<3>({0., 1., 0.}), triangulation);
     GridTools::shift(Tensor<1, 3>({-half_length, 0.0, 0.0}), triangulation);
     // At this point we have a cylinder. Multiply the y and z coordinates by a
     // factor that scales (with x) linearly between radius_0 and radius_1 to fix
index 4b37c3eb1ffced1077be61d3b2c88be058c82e4a..9abd1b773362184e30b2b386fb66220e4003af1a 100644 (file)
@@ -56,6 +56,8 @@
 #include <deal.II/numerics/matrix_tools.h>
 #include <deal.II/numerics/vector_tools_integrate_difference.h>
 
+#include <deal.II/physics/transformations.h>
+
 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
 #include <boost/random/mersenne_twister.hpp>
 #include <boost/random/uniform_real_distribution.hpp>
@@ -1970,33 +1972,23 @@ namespace GridTools
     class Rotate3d
     {
     public:
-      Rotate3d(const double angle, const unsigned int axis)
-        : angle(angle)
-        , axis(axis)
+      Rotate3d(const double angle, const Point<3, double> &axis)
+        : rotation_matrix(
+            Physics::Transformations::Rotations::rotation_matrix_3d(axis,
+                                                                    angle))
       {}
 
       Point<3>
       operator()(const Point<3> &p) const
       {
-        if (axis == 0)
-          return {p(0),
-                  std::cos(angle) * p(1) - std::sin(angle) * p(2),
-                  std::sin(angle) * p(1) + std::cos(angle) * p(2)};
-        else if (axis == 1)
-          return {std::cos(angle) * p(0) + std::sin(angle) * p(2),
-                  p(1),
-                  -std::sin(angle) * p(0) + std::cos(angle) * p(2)};
-        else
-          return {std::cos(angle) * p(0) - std::sin(angle) * p(1),
-                  std::sin(angle) * p(0) + std::cos(angle) * p(1),
-                  p(2)};
+        return static_cast<Point<3>>(rotation_matrix * p);
       }
 
     private:
-      const double       angle;
-      const unsigned int axis;
+      const Tensor<2, 3, double> rotation_matrix;
     };
 
+
     template <int spacedim>
     class Scale
     {
@@ -2025,6 +2017,16 @@ namespace GridTools
   }
 
 
+  template <int dim>
+  void
+  rotate(const double            angle,
+         const Point<3, double> &axis,
+         Triangulation<dim, 3> & triangulation)
+  {
+    transform(internal::Rotate3d(angle, axis), triangulation);
+  }
+
+
   template <int dim>
   void
   rotate(const double           angle,
@@ -2033,9 +2035,13 @@ namespace GridTools
   {
     Assert(axis < 3, ExcMessage("Invalid axis given!"));
 
-    transform(internal::Rotate3d(angle, axis), triangulation);
+    Point<3, double> vector;
+    vector[axis] = 1.;
+
+    transform(internal::Rotate3d(angle, vector), triangulation);
   }
 
+
   template <int dim, int spacedim>
   void
   scale(const double                  scaling_factor,
index eb2215c3f8f51d2299536431bdec4e08e3b0fda7..193bfe7dfc705b6a2c594048c83d5b7398128cf8 100644 (file)
@@ -271,6 +271,12 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
         Triangulation<deal_II_dimension, deal_II_space_dimension> &);
 
 #  if deal_II_space_dimension == 3
+      template void
+      rotate<deal_II_dimension>(
+        const double,
+        const Point<deal_II_space_dimension, double> &,
+        Triangulation<deal_II_dimension, deal_II_space_dimension> &);
+
       template void
       rotate<deal_II_dimension>(
         const double,
index 158d76d113618170454a16f3a49b8531afa5164b..38f9d14cc2263b5048267f159553bb6d84345a43 100644 (file)
@@ -17,6 +17,8 @@
 
 #include <deal.II/grid/grid_tools.h>
 
+#include <deal.II/physics/transformations.h>
+
 #include <vector>
 
 // GridTools functions that are template specializations (i.e., only compiled
@@ -399,24 +401,23 @@ namespace GridTools
 
   namespace
   {
-    // the following class is only
-    // needed in 2d, so avoid trouble
-    // with compilers warning otherwise
+    // the following class is only needed in 2d, so avoid trouble with compilers
+    // warning otherwise
     class Rotate2d
     {
     public:
       explicit Rotate2d(const double angle)
-        : angle(angle)
+        : rotation_matrix(
+            Physics::Transformations::Rotations::rotation_matrix_2d(angle))
       {}
       Point<2>
       operator()(const Point<2> &p) const
       {
-        return {std::cos(angle) * p(0) - std::sin(angle) * p(1),
-                std::sin(angle) * p(0) + std::cos(angle) * p(1)};
+        return static_cast<Point<2>>(rotation_matrix * p);
       }
 
     private:
-      const double angle;
+      const Tensor<2, 2, double> rotation_matrix;
     };
   } // namespace
 
index 46703c566d43a29e98968262b64b9b7c77d7af69..352de8c28de9387db1dcf81f855b56570ceb37dd 100644 (file)
@@ -19,6 +19,7 @@
 // and on each of these cells each time we evaluate the shape
 // functions.
 
+#include <deal.II/base/numbers.h>
 #include <deal.II/base/quadrature_lib.h>
 
 #include <deal.II/dofs/dof_accessor.h>
@@ -69,7 +70,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // second round: rotate
       // triangulation
       case 1:
-        GridTools::rotate(3.14159265358 / 2, tria);
+        GridTools::rotate(numbers::PI_2, tria);
         break;
 
       // third round: inflate
@@ -83,7 +84,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // stretch
       case 3:
         GridTools::scale(.5, tria);
-        GridTools::rotate(-3.14159265358 / 2, tria);
+        GridTools::rotate(-numbers::PI_2, tria);
         GridTools::transform(&stretch_coordinates, tria);
 
         break;
index 7e074536cd51447ed1285a45afb80cf421b96749..4396c243176c8dc76b063a5bb458c65ff50c8b03 100644 (file)
@@ -25,6 +25,7 @@
 
 // perl -n -e 'print if s/DEAL:NedelecK-TransformN::value//' output
 
+#include <deal.II/base/numbers.h>
 #include <deal.II/base/quadrature_lib.h>
 
 #include <deal.II/dofs/dof_accessor.h>
@@ -75,7 +76,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // second round: rotate
       // triangulation
       case 1:
-        GridTools::rotate(3.14159265358 / 2, tria);
+        GridTools::rotate(numbers::PI_2, tria);
         break;
 
       // third round: inflate
@@ -89,7 +90,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // stretch
       case 3:
         GridTools::scale(.5, tria);
-        GridTools::rotate(-3.14159265358 / 2, tria);
+        GridTools::rotate(-numbers::PI_2, tria);
         GridTools::transform(&stretch_coordinates, tria);
 
         break;
index 92df4009d4d5536523fe6ffe120d87e8af7a8280..1757083cd72f01ee4ba16e404d72b550abbf6560 100644 (file)
@@ -19,6 +19,7 @@
 // and on each of these cells each time we evaluate the shape
 // functions.
 
+#include <deal.II/base/numbers.h>
 #include <deal.II/base/quadrature_lib.h>
 
 #include <deal.II/dofs/dof_accessor.h>
@@ -69,7 +70,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // second round: rotate
       // triangulation
       case 1:
-        GridTools::rotate(3.14159265358 / 2, tria);
+        GridTools::rotate(numbers::PI_2, tria);
         break;
 
       // third round: inflate
@@ -83,7 +84,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // stretch
       case 3:
         GridTools::scale(.5, tria);
-        GridTools::rotate(-3.14159265358 / 2, tria);
+        GridTools::rotate(-numbers::PI_2, tria);
         GridTools::transform(&stretch_coordinates, tria);
 
         break;
index e6462d54a08a5b67b79cb40e99ec624e35ae24e7..5d44fb2a027297f21129a0a1d8e8064658a72b5f 100644 (file)
@@ -19,6 +19,7 @@
 // and on each of these transformed cells we evaluate the shape
 // functions.
 
+#include <deal.II/base/numbers.h>
 #include <deal.II/base/quadrature_lib.h>
 
 #include <deal.II/dofs/dof_accessor.h>
@@ -67,7 +68,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
 
       // second test: rotate cell
       case 1:
-        GridTools::rotate(3.14159265358 / 2, tria);
+        GridTools::rotate(numbers::PI_2, tria);
         break;
 
       // third test: scale cell by a factor of 2
@@ -79,7 +80,7 @@ transform_grid(Triangulation<2> &tria, const unsigned int transform)
       // stretch
       case 3:
         GridTools::scale(.5, tria);
-        GridTools::rotate(-3.14159265358 / 2, tria);
+        GridTools::rotate(-numbers::PI_2, tria);
         GridTools::transform(&stretch_coordinates, tria);
 
         break;
index 0196fefcda029fb75a5b43ee12eba98f96078339..667fb32566702b772573becbd043fe16b1a48c33 100644 (file)
@@ -15,6 +15,8 @@
 
 
 
+#include <deal.II/base/numbers.h>
+
 #include <deal.II/grid/grid_generator.h>
 #include <deal.II/grid/grid_out.h>
 #include <deal.II/grid/grid_tools.h>
@@ -77,7 +79,7 @@ test2()
   GridOut().write_gnuplot(tria, logfile);
 
   logfile << "Rotated grid:" << std::endl;
-  GridTools::rotate(3.14159265258 / 4, tria);
+  GridTools::rotate(numbers::PI_4, tria);
   GridOut().write_gnuplot(tria, logfile);
 }
 
index f3f75a3f7f1c8c0d0c481bdbca939b9a7a35e092..a49faef88c76d2718335dfe381fbfd9101d6252b 100644 (file)
@@ -15,6 +15,8 @@
 
 // test GridTools::rotate(), output is checked for correctness visually
 
+#include <deal.II/base/numbers.h>
+
 #include <deal.II/grid/grid_generator.h>
 #include <deal.II/grid/grid_out.h>
 #include <deal.II/grid/grid_tools.h>
@@ -63,11 +65,11 @@ test()
   GridGenerator::subdivided_parallelepiped<dim, spacedim>(tria, origin, edges);
 
   // GridOut().write_gnuplot (tria, deallog.get_file_stream());
-  GridTools::rotate(numbers::PI / 3.0, 0, tria);
+  GridTools::rotate(numbers::PI / 3.0, Point<3>({1., 0., 0.}), tria);
   // GridOut().write_gnuplot (tria, deallog.get_file_stream());
-  GridTools::rotate(numbers::PI / 10.0, 1, tria);
+  GridTools::rotate(numbers::PI / 10.0, Point<3>({0., 1., 0.}), tria);
   // GridOut().write_gnuplot (tria, deallog.get_file_stream());
-  GridTools::rotate(-numbers::PI / 5.0, 2, tria);
+  GridTools::rotate(-numbers::PI / 5.0, Point<3>({0., 0., 1.}), tria);
   GridOut().write_gnuplot(tria, deallog.get_file_stream());
 }
 

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.