From: Wolfgang Bangerth Date: Tue, 30 Dec 2014 14:18:05 +0000 (-0600) Subject: Remove deprecated function GridGenerator::laplace_transformation. Move the implementa... X-Git-Tag: v8.3.0-rc1~566^2~15 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cff0b2aa1c8b45a80da52285c7a7026d01d34020;p=dealii.git Remove deprecated function GridGenerator::laplace_transformation. Move the implementation to its new home in GridTools. --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 7ab87361e9..147bcaf230 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -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.
  • Removed: The config.h file no longer exports HAVE_* definitions. diff --git a/include/deal.II/grid/grid_generator.h b/include/deal.II/grid/grid_generator.h index a4e7fd64a4..46745406ab 100644 --- a/include/deal.II/grid/grid_generator.h +++ b/include/deal.II/grid/grid_generator.h @@ -42,11 +42,6 @@ template 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 - void laplace_transformation (Triangulation &tria, - const std::map > &new_points, - const Function *coefficient = 0) DEAL_II_DEPRECATED; - - /* - * @} - */ - /** * @name Exceptions * @{ diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc index 147d5d6357..e27fe2e1fc 100644 --- a/source/grid/grid_generator.cc +++ b/source/grid/grid_generator.cc @@ -13,16 +13,6 @@ // // --------------------------------------------------------------------- -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include #include #include #include @@ -31,12 +21,6 @@ #include #include #include -#include -#include -#include -#include -#include -#include #include @@ -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 &S, - const std::map &m, - Vector &u) - { - const unsigned int n_dofs=S.n(); - FilteredMatrix > SF (S); - PreconditionJacobi > prec; - prec.initialize(S, 1.2); - FilteredMatrix > PF (prec); - - SolverControl control (n_dofs, 1.e-10, false, false); - GrowingVectorMemory > mem; - SolverCG > solver (control, mem); - - Vector 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 > &, - const Function<1> *) - { - Assert(false, ExcNotImplemented()); - } - - -// Implementation for dimensions except 1 - template - void laplace_transformation (Triangulation &tria, - const std::map > &new_points, - const Function *coefficient) - { - // first provide everything that is - // needed for solving a Laplace - // equation. - MappingQ1 mapping_q1; - FE_Q q1(1); - - DoFHandler 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 S(sparsity_pattern); - - QGauss quadrature(4); - - MatrixCreator::create_laplace_matrix(mapping_q1, dof_handler, quadrature, S,coefficient); - - // set up the boundary values for - // the laplace problem - std::vector > m(dim); - typename std::map >::const_iterator map_end=new_points.end(); - - // fill these maps using the data - // given by new_points - typename DoFHandler::cell_iterator cell=dof_handler.begin_active(), - endc=dof_handler.end(); - for (; cell!=endc; ++cell) - { - for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) - { - const typename DoFHandler::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::vertices_per_face; ++vertex_no) - { - const unsigned int vertex_index=face->vertex_index(vertex_no); - - const typename std::map >::const_iterator map_iter - = new_points.find(vertex_index); - - if (map_iter!=map_end) - for (unsigned int i=0; i ( - face->vertex_dof_index(vertex_no, 0), map_iter->second(i))); - } - } - } - - // solve the dim problems with - // different right hand sides. - Vector us[dim]; - for (unsigned int i=0; i tasks; - for (unsigned int i=0; i::vertices_per_cell; ++vertex_no) - { - Point &v=cell->vertex(vertex_no); - const unsigned int dof_index=cell->vertex_dof_index(vertex_no, 0); - for (unsigned int i=0; i void hyper_cube_with_cylindrical_hole (Triangulation<1> &, const double, diff --git a/source/grid/grid_generator.inst.in b/source/grid/grid_generator.inst.in index c74fb017ea..017dd2e5f2 100644 --- a/source/grid/grid_generator.inst.in +++ b/source/grid/grid_generator.inst.in @@ -127,13 +127,6 @@ namespace GridGenerator \{ const bool); #endif -#if deal_II_dimension > 1 - template void - laplace_transformation (Triangulation &, - const std::map > &, - const Function *); -#endif - \} } diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 95ee495eee..0e7b739693 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -15,6 +15,18 @@ #include #include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include #include #include #include @@ -22,9 +34,6 @@ #include #include #include -#include -#include -#include #include #include #include @@ -34,6 +43,7 @@ #include #include #include +#include #include #include @@ -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 &S, + const std::map &m, + Vector &u) + { + const unsigned int n_dofs=S.n(); + FilteredMatrix > SF (S); + PreconditionJacobi > prec; + prec.initialize(S, 1.2); + FilteredMatrix > PF (prec); + + SolverControl control (n_dofs, 1.e-10, false, false); + GrowingVectorMemory > mem; + SolverCG > solver (control, mem); + + Vector 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 > &, + Triangulation<1> &, + const Function<1> *) + { + Assert(false, ExcNotImplemented()); + } + + + // Implementation for dimensions except 1 template void laplace_transform (const std::map > &new_points, Triangulation &triangulation, const Function *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 mapping_q1; + FE_Q q1(1); + + DoFHandler 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 S(sparsity_pattern); + + QGauss quadrature(4); + + MatrixCreator::create_laplace_matrix(mapping_q1, dof_handler, quadrature, S,coefficient); + + // set up the boundary values for + // the laplace problem + std::vector > m(dim); + typename std::map >::const_iterator map_end=new_points.end(); + + // fill these maps using the data + // given by new_points + typename DoFHandler::cell_iterator cell=dof_handler.begin_active(), + endc=dof_handler.end(); + for (; cell!=endc; ++cell) + { + for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) + { + const typename DoFHandler::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::vertices_per_face; ++vertex_no) + { + const unsigned int vertex_index=face->vertex_index(vertex_no); + + const typename std::map >::const_iterator map_iter + = new_points.find(vertex_index); + + if (map_iter!=map_end) + for (unsigned int i=0; i ( + face->vertex_dof_index(vertex_no, 0), map_iter->second(i))); + } + } + } + + // solve the dim problems with + // different right hand sides. + Vector us[dim]; + for (unsigned int i=0; i tasks; + for (unsigned int i=0; i::vertices_per_cell; ++vertex_no) + { + Point &v=cell->vertex(vertex_no); + const unsigned int dof_index=cell->vertex_dof_index(vertex_no, 0); + for (unsigned int i=0; i &result); #if deal_II_dimension == deal_II_space_dimension +# if deal_II_dimension > 1 template void laplace_transform (const std::map > &new_points, Triangulation &triangulation, const Function *coefficient); +# endif template Triangulation::DistortedCellList diff --git a/tests/deal.II/grid_transform.cc b/tests/deal.II/grid_transform.cc index fc54388ab8..ef8bbe1263 100644 --- a/tests/deal.II/grid_transform.cc +++ b/tests/deal.II/grid_transform.cc @@ -24,6 +24,7 @@ #include #include #include +#include #include #include @@ -97,7 +98,7 @@ int main () } } - GridGenerator::laplace_transformation (tria, new_points); + GridTools::laplace_transform (new_points, tria); HyperBallBoundary inner_ball(n_center, n_radius); tria.set_boundary(1, inner_ball); diff --git a/tests/deal.II/grid_transform_3d.cc b/tests/deal.II/grid_transform_3d.cc index 47f1cb0e51..ef9883f644 100644 --- a/tests/deal.II/grid_transform_3d.cc +++ b/tests/deal.II/grid_transform_3d.cc @@ -24,6 +24,7 @@ #include #include #include +#include #include #include @@ -86,7 +87,7 @@ int main () } } - GridGenerator::laplace_transformation (tria, new_points); + GridTools::laplace_transform (new_points, tria); GridOut grid_out; diff --git a/tests/deal.II/grid_transform_coefficient.cc b/tests/deal.II/grid_transform_coefficient.cc index 6c3d7c3a43..714d37d189 100644 --- a/tests/deal.II/grid_transform_coefficient.cc +++ b/tests/deal.II/grid_transform_coefficient.cc @@ -22,6 +22,7 @@ #include #include #include +#include #include #include @@ -108,7 +109,7 @@ int main () } Coefficient c; - GridGenerator::laplace_transformation (tria, new_points, &c); + GridTools::laplace_transform (new_points, tria, &c); HyperBallBoundary inner_ball(n_center, n_radius); tria.set_boundary(1, inner_ball);