]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove deprecated function GridGenerator::laplace_transformation. Move the implementa...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 30 Dec 2014 14:18:05 +0000 (08:18 -0600)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 4 Jan 2015 21:13:44 +0000 (15:13 -0600)
doc/news/changes.h
include/deal.II/grid/grid_generator.h
source/grid/grid_generator.cc
source/grid/grid_generator.inst.in
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/deal.II/grid_transform.cc
tests/deal.II/grid_transform_3d.cc
tests/deal.II/grid_transform_coefficient.cc

index 7ab87361e945f82d3a5152a084d87f4363ea81e6..147bcaf23045f5981206d9bc5e5a451c291f1530 100644 (file)
@@ -55,6 +55,7 @@ inconvenience this causes.
   - Deprecated variants of MeshWorker::loop and MeshWorker::integration_loop.
   - ThreadManagement::spawn.
   - Threads::ThreadCondition and Threads::ThreadMutex.
+  - GridGenerator::laplace_transformation.
 </ol>
 
   <li> Removed: The config.h file no longer exports HAVE_* definitions.
index a4e7fd64a4b336b2c7f95db4a4c19227d89123cb..46745406ab915d26123459c17cd31b2337d2ec93 100644 (file)
@@ -42,11 +42,6 @@ template <typename number> class SparseMatrix;
  * allowing them to be distinguished for the purpose of attaching geometry
  * objects and evaluating different boundary conditions.
  *
- * This namespace also provides a function
- * GridGenerator::laplace_transformation that smoothly transforms a domain
- * into another one. This can be used to transform basic geometries to more
- * complicated ones, like a shell to a grid of an airfoil, for example.
- *
  * @ingroup grid
  */
 namespace GridGenerator
@@ -851,42 +846,6 @@ namespace GridGenerator
    * @}
    */
 
-  /**
-   * @name Deprecated functions
-   * @{
-   */
-
-  /**
-   * This function transforms the @p Triangulation @p tria smoothly to a
-   * domain that is described by the boundary points in the map @p new_points.
-   * This map maps the point indices to the boundary points in the transformed
-   * domain.
-   *
-   * Note, that the @p Triangulation is changed in-place, therefore you don't
-   * need to keep two triangulations, but the given triangulation is changed
-   * (overwritten).
-   *
-   * In 1d, this function is not currently implemented.
-   *
-   * An optional @p coefficient for the Laplace problem an be used to control
-   * the amount of mesh deformation in different parts of the domain. Larger
-   * values make cells less prone to deformation (effectively increasing their
-   * stiffness). The coefficient is evaluated in the coordinate system of the
-   * old, undeformed configuration of the triangulation as input, i.e., before
-   * the transformation is applied. Should this function be provided, sensible
-   * results can only be expected if all coefficients are positive.
-   *
-   * @deprecated This function has been moved to GridTools::laplace_transform
-   */
-  template <int dim>
-  void laplace_transformation (Triangulation<dim> &tria,
-                               const std::map<unsigned int,Point<dim> > &new_points,
-                               const Function<dim> *coefficient = 0) DEAL_II_DEPRECATED;
-
-  /*
-   * @}
-   */
-
   /**
    * @name Exceptions
    * @{
index 147d5d6357523eef1dcbd3c6892a9ddf284b1c73..e27fe2e1fc3179840185ded923e0a529e5ef7a6a 100644 (file)
 //
 // ---------------------------------------------------------------------
 
-#include <deal.II/base/quadrature_lib.h>
-#include <deal.II/base/thread_management.h>
-#include <deal.II/lac/vector.h>
-#include <deal.II/lac/vector_memory.h>
-#include <deal.II/lac/filtered_matrix.h>
-#include <deal.II/lac/precondition.h>
-#include <deal.II/lac/solver_cg.h>
-#include <deal.II/lac/sparse_matrix.h>
-#include <deal.II/lac/compressed_sparsity_pattern.h>
-#include <deal.II/lac/constraint_matrix.h>
 #include <deal.II/grid/grid_generator.h>
 #include <deal.II/grid/grid_reordering.h>
 #include <deal.II/grid/grid_tools.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_boundary_lib.h>
 #include <deal.II/grid/intergrid_map.h>
-#include <deal.II/dofs/dof_handler.h>
-#include <deal.II/dofs/dof_accessor.h>
-#include <deal.II/dofs/dof_tools.h>
-#include <deal.II/fe/mapping_q1.h>
-#include <deal.II/fe/fe_q.h>
-#include <deal.II/numerics/matrix_tools.h>
 
 #include <deal.II/distributed/tria.h>
 
@@ -3573,135 +3557,6 @@ namespace GridGenerator
   }
 
 
-  /**
-   * Solve the Laplace equation for @p laplace_transformation function for one
-   * of the @p dim space dimensions. Factorized into a function of its own
-   * in order to allow parallel execution.
-   */
-  void laplace_solve (const SparseMatrix<double> &S,
-                      const std::map<unsigned int,double> &m,
-                      Vector<double> &u)
-  {
-    const unsigned int n_dofs=S.n();
-    FilteredMatrix<Vector<double> > SF (S);
-    PreconditionJacobi<SparseMatrix<double> > prec;
-    prec.initialize(S, 1.2);
-    FilteredMatrix<Vector<double> > PF (prec);
-
-    SolverControl control (n_dofs, 1.e-10, false, false);
-    GrowingVectorMemory<Vector<double> > mem;
-    SolverCG<Vector<double> > solver (control, mem);
-
-    Vector<double> f(n_dofs);
-
-    SF.add_constraints(m);
-    SF.apply_constraints (f, true);
-    solver.solve(SF, u, f, PF);
-  }
-
-
-// Implementation for 1D only
-  template <>
-  void laplace_transformation (Triangulation<1> &,
-                               const std::map<unsigned int,Point<1> > &,
-                               const Function<1> *)
-  {
-    Assert(false, ExcNotImplemented());
-  }
-
-
-// Implementation for dimensions except 1
-  template <int dim>
-  void laplace_transformation (Triangulation<dim> &tria,
-                               const std::map<unsigned int,Point<dim> > &new_points,
-                               const Function<dim> *coefficient)
-  {
-    // first provide everything that is
-    // needed for solving a Laplace
-    // equation.
-    MappingQ1<dim> mapping_q1;
-    FE_Q<dim> q1(1);
-
-    DoFHandler<dim> dof_handler(tria);
-    dof_handler.distribute_dofs(q1);
-
-    CompressedSparsityPattern c_sparsity_pattern (dof_handler.n_dofs (),
-                                                  dof_handler.n_dofs ());
-    DoFTools::make_sparsity_pattern (dof_handler, c_sparsity_pattern);
-    c_sparsity_pattern.compress ();
-
-    SparsityPattern sparsity_pattern;
-    sparsity_pattern.copy_from (c_sparsity_pattern);
-    sparsity_pattern.compress ();
-
-    SparseMatrix<double> S(sparsity_pattern);
-
-    QGauss<dim> quadrature(4);
-
-    MatrixCreator::create_laplace_matrix(mapping_q1, dof_handler, quadrature, S,coefficient);
-
-    // set up the boundary values for
-    // the laplace problem
-    std::vector<std::map<unsigned int,double> > m(dim);
-    typename std::map<unsigned int,Point<dim> >::const_iterator map_end=new_points.end();
-
-    // fill these maps using the data
-    // given by new_points
-    typename DoFHandler<dim>::cell_iterator cell=dof_handler.begin_active(),
-                                            endc=dof_handler.end();
-    for (; cell!=endc; ++cell)
-      {
-        for (unsigned int face_no=0; face_no<GeometryInfo<dim>::faces_per_cell; ++face_no)
-          {
-            const typename DoFHandler<dim>::face_iterator face=cell->face(face_no);
-
-            // loop over all vertices of the cell and see if it is listed in the map
-            // given as first argument of the function
-            for (unsigned int vertex_no=0;
-                 vertex_no<GeometryInfo<dim>::vertices_per_face; ++vertex_no)
-              {
-                const unsigned int vertex_index=face->vertex_index(vertex_no);
-
-                const typename std::map<unsigned int,Point<dim> >::const_iterator map_iter
-                  = new_points.find(vertex_index);
-
-                if (map_iter!=map_end)
-                  for (unsigned int i=0; i<dim; ++i)
-                    m[i].insert(std::pair<unsigned int,double> (
-                                  face->vertex_dof_index(vertex_no, 0), map_iter->second(i)));
-              }
-          }
-      }
-
-    // solve the dim problems with
-    // different right hand sides.
-    Vector<double> us[dim];
-    for (unsigned int i=0; i<dim; ++i)
-      us[i].reinit (dof_handler.n_dofs());
-
-    // solve linear systems in parallel
-    Threads::TaskGroup<> tasks;
-    for (unsigned int i=0; i<dim; ++i)
-      tasks += Threads::new_task (&laplace_solve,
-                                  S, m[i], us[i]);
-    tasks.join_all ();
-
-    // change the coordinates of the
-    // points of the triangulation
-    // according to the computed values
-    for (cell=dof_handler.begin_active(); cell!=endc; ++cell)
-      for (unsigned int vertex_no=0;
-           vertex_no<GeometryInfo<dim>::vertices_per_cell; ++vertex_no)
-        {
-          Point<dim> &v=cell->vertex(vertex_no);
-          const unsigned int dof_index=cell->vertex_dof_index(vertex_no, 0);
-          for (unsigned int i=0; i<dim; ++i)
-            v(i)=us[i](dof_index);
-        }
-  }
-
-
-
   template <>
   void hyper_cube_with_cylindrical_hole (Triangulation<1> &,
                                          const double,
index c74fb017ea700b31fdf069f4f181304cb927650f..017dd2e5f288ef298d78a208b7dcbb5138bd8fb0 100644 (file)
@@ -127,13 +127,6 @@ namespace GridGenerator \{
        const bool);
   #endif
   
-#if deal_II_dimension > 1
-  template void
-    laplace_transformation<deal_II_dimension> (Triangulation<deal_II_dimension> &,
-                                              const std::map<unsigned int,Point<deal_II_dimension> > &,
-                                              const Function<deal_II_dimension> *);
-#endif
-  
 \}  
  }
 
index 95ee495eee291cfb825300ecf1c0a897f60433fe..0e7b73969395b8152275ca6588351335454fc578 100644 (file)
 
 #include <deal.II/base/std_cxx11/array.h>
 #include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/thread_management.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/vector_memory.h>
+#include <deal.II/lac/filtered_matrix.h>
+#include <deal.II/lac/precondition.h>
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/compressed_sparsity_pattern.h>
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/lac/sparsity_tools.h>
+#include <deal.II/lac/compressed_sparsity_pattern.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/distributed/tria.h>
 #include <deal.II/grid/tria_accessor.h>
@@ -22,9 +34,6 @@
 #include <deal.II/grid/tria_boundary.h>
 #include <deal.II/grid/grid_generator.h>
 #include <deal.II/grid/grid_tools.h>
-#include <deal.II/lac/sparsity_pattern.h>
-#include <deal.II/lac/sparsity_tools.h>
-#include <deal.II/lac/compressed_sparsity_pattern.h>
 #include <deal.II/dofs/dof_handler.h>
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/dofs/dof_tools.h>
@@ -34,6 +43,7 @@
 #include <deal.II/fe/fe_values.h>
 #include <deal.II/hp/mapping_collection.h>
 #include <deal.II/multigrid/mg_dof_handler.h>
+#include <deal.II/numerics/matrix_tools.h>
 
 #include <boost/random/uniform_real_distribution.hpp>
 #include <boost/random/mersenne_twister.hpp>
@@ -627,15 +637,136 @@ namespace GridTools
   }
 
 
+  namespace
+  {
+    /**
+     * Solve the Laplace equation for the @p laplace_transform function for one
+     * of the @p dim space dimensions. Factorized into a function of its own
+     * in order to allow parallel execution.
+     */
+    void laplace_solve (const SparseMatrix<double> &S,
+                        const std::map<unsigned int,double> &m,
+                        Vector<double> &u)
+    {
+      const unsigned int n_dofs=S.n();
+      FilteredMatrix<Vector<double> > SF (S);
+      PreconditionJacobi<SparseMatrix<double> > prec;
+      prec.initialize(S, 1.2);
+      FilteredMatrix<Vector<double> > PF (prec);
+
+      SolverControl control (n_dofs, 1.e-10, false, false);
+      GrowingVectorMemory<Vector<double> > mem;
+      SolverCG<Vector<double> > solver (control, mem);
+
+      Vector<double> f(n_dofs);
+
+      SF.add_constraints(m);
+      SF.apply_constraints (f, true);
+      solver.solve(SF, u, f, PF);
+    }
+  }
+
+
+
+  // Implementation for 1D only
+  template <>
+  void laplace_transform (const std::map<unsigned int,Point<1> > &,
+                          Triangulation<1> &,
+                          const Function<1> *)
+  {
+    Assert(false, ExcNotImplemented());
+  }
+
+
+  // Implementation for dimensions except 1
   template <int dim>
   void
   laplace_transform (const std::map<unsigned int,Point<dim> > &new_points,
                      Triangulation<dim> &triangulation,
                      const Function<dim> *coefficient)
   {
-    //TODO: Move implementation of this function into the current
-    // namespace
-    GridGenerator::laplace_transformation(triangulation, new_points, coefficient);
+    // first provide everything that is
+    // needed for solving a Laplace
+    // equation.
+    MappingQ1<dim> mapping_q1;
+    FE_Q<dim> q1(1);
+
+    DoFHandler<dim> dof_handler(triangulation);
+    dof_handler.distribute_dofs(q1);
+
+    CompressedSparsityPattern c_sparsity_pattern (dof_handler.n_dofs (),
+                                                  dof_handler.n_dofs ());
+    DoFTools::make_sparsity_pattern (dof_handler, c_sparsity_pattern);
+    c_sparsity_pattern.compress ();
+
+    SparsityPattern sparsity_pattern;
+    sparsity_pattern.copy_from (c_sparsity_pattern);
+    sparsity_pattern.compress ();
+
+    SparseMatrix<double> S(sparsity_pattern);
+
+    QGauss<dim> quadrature(4);
+
+    MatrixCreator::create_laplace_matrix(mapping_q1, dof_handler, quadrature, S,coefficient);
+
+    // set up the boundary values for
+    // the laplace problem
+    std::vector<std::map<unsigned int,double> > m(dim);
+    typename std::map<unsigned int,Point<dim> >::const_iterator map_end=new_points.end();
+
+    // fill these maps using the data
+    // given by new_points
+    typename DoFHandler<dim>::cell_iterator cell=dof_handler.begin_active(),
+                                            endc=dof_handler.end();
+    for (; cell!=endc; ++cell)
+      {
+        for (unsigned int face_no=0; face_no<GeometryInfo<dim>::faces_per_cell; ++face_no)
+          {
+            const typename DoFHandler<dim>::face_iterator face=cell->face(face_no);
+
+            // loop over all vertices of the cell and see if it is listed in the map
+            // given as first argument of the function
+            for (unsigned int vertex_no=0;
+                 vertex_no<GeometryInfo<dim>::vertices_per_face; ++vertex_no)
+              {
+                const unsigned int vertex_index=face->vertex_index(vertex_no);
+
+                const typename std::map<unsigned int,Point<dim> >::const_iterator map_iter
+                  = new_points.find(vertex_index);
+
+                if (map_iter!=map_end)
+                  for (unsigned int i=0; i<dim; ++i)
+                    m[i].insert(std::pair<unsigned int,double> (
+                                  face->vertex_dof_index(vertex_no, 0), map_iter->second(i)));
+              }
+          }
+      }
+
+    // solve the dim problems with
+    // different right hand sides.
+    Vector<double> us[dim];
+    for (unsigned int i=0; i<dim; ++i)
+      us[i].reinit (dof_handler.n_dofs());
+
+    // solve linear systems in parallel
+    Threads::TaskGroup<> tasks;
+    for (unsigned int i=0; i<dim; ++i)
+      tasks += Threads::new_task (&laplace_solve,
+                                  S, m[i], us[i]);
+    tasks.join_all ();
+
+    // change the coordinates of the
+    // points of the triangulation
+    // according to the computed values
+    for (cell=dof_handler.begin_active(); cell!=endc; ++cell)
+      for (unsigned int vertex_no=0;
+           vertex_no<GeometryInfo<dim>::vertices_per_cell; ++vertex_no)
+        {
+          Point<dim> &v=cell->vertex(vertex_no);
+          const unsigned int dof_index=cell->vertex_dof_index(vertex_no, 0);
+          for (unsigned int i=0; i<dim; ++i)
+            v(i)=us[i](dof_index);
+        }
   }
 
 
index 6749419546d6c172116842153d60b3f76eb2d569..d31b93d5eca25ed58964232ad4366aa400791616 100644 (file)
@@ -238,11 +238,13 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
        Triangulation<deal_II_dimension, deal_II_space_dimension>       &result);
 
 #if deal_II_dimension == deal_II_space_dimension
+#  if deal_II_dimension > 1
     template
     void
     laplace_transform (const std::map<unsigned int,Point<deal_II_dimension> > &new_points,
                        Triangulation<deal_II_dimension> &triangulation,
                        const Function<deal_II_dimension> *coefficient);
+#  endif
 
     template
       Triangulation<deal_II_dimension,deal_II_space_dimension>::DistortedCellList
index fc54388ab82d21e0836c5e0563945f86d9f5f7fb..ef8bbe12631211be9107c526204159a83ca7b91c 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_boundary_lib.h>
 #include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
 #include <deal.II/grid/grid_out.h>
 #include <deal.II/fe/mapping_q.h>
 
@@ -97,7 +98,7 @@ int main ()
           }
     }
 
-  GridGenerator::laplace_transformation (tria, new_points);
+  GridTools::laplace_transform (new_points, tria);
   HyperBallBoundary<dim> inner_ball(n_center, n_radius);
   tria.set_boundary(1, inner_ball);
 
index 47f1cb0e51595044982be8a94ef03e005e2cc573..ef9883f644e90ff870be3bdc123acd5540657411 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_boundary_lib.h>
 #include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
 #include <deal.II/grid/grid_out.h>
 
 #include <fstream>
@@ -86,7 +87,7 @@ int main ()
             }
       }
 
-  GridGenerator::laplace_transformation (tria, new_points);
+  GridTools::laplace_transform (new_points, tria);
 
 
   GridOut grid_out;
index 6c3d7c3a430bd8a2068f64bf5d6323386bdd8634..714d37d189701a5760346d6b1184929c629f9567 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_boundary_lib.h>
 #include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
 #include <deal.II/grid/grid_out.h>
 #include <deal.II/fe/mapping_q.h>
 
@@ -108,7 +109,7 @@ int main ()
     }
 
   Coefficient<dim> c;
-  GridGenerator::laplace_transformation (tria, new_points, &c);
+  GridTools::laplace_transform (new_points, tria, &c);
   HyperBallBoundary<dim> inner_ball(n_center, n_radius);
   tria.set_boundary(1, inner_ball);
 

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.