From 1a4fb99a015a094bd6b1695c67276abefe296902 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Thu, 19 Feb 2015 10:50:55 +0100 Subject: [PATCH] Remove a few 1D specializations for mapping. The old implementation of MappingQEulerian only moved the vertices in 1D, not interior points for higher order mappings. While it might be debatable how useful non-linear mappings are in 1D, this was not documented and lead to surprising behavior. By making the class more general we can actually make our code base simpler because we can remove a few template specialization and use the same code as in higher dimensions. --- doc/news/changes.h | 7 ++ include/deal.II/fe/mapping_q.h | 7 -- include/deal.II/grid/tria_boundary.h | 6 -- source/fe/mapping_q.cc | 48 ------------- source/grid/tria_boundary.cc | 36 ---------- tests/bits/mapping_q_eulerian_06.cc | 91 +++++++++++++++++++++++++ tests/bits/mapping_q_eulerian_06.output | 13 ++++ 7 files changed, 111 insertions(+), 97 deletions(-) create mode 100644 tests/bits/mapping_q_eulerian_06.cc create mode 100644 tests/bits/mapping_q_eulerian_06.output diff --git a/doc/news/changes.h b/doc/news/changes.h index 23518e1e70..77ca21aa83 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -350,6 +350,13 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: MappingQEulerian would previously not move interior points + in 1D for higher order mappings. This has been fixed by removing a few + specializations of MappingQ for 1D that are no longer necessary. +
    + (Martin Kronbichler, 2015/02/19) +
    2. +
  2. Fixed: The implementation of the class GrowingVectorMemory has been moved from source/lac/vector_memory.cc to the new file include/deal.II/lac/vector_memory.templates.h. This allows users to diff --git a/include/deal.II/fe/mapping_q.h b/include/deal.II/fe/mapping_q.h index 87a8e628ee..9f8957068b 100644 --- a/include/deal.II/fe/mapping_q.h +++ b/include/deal.II/fe/mapping_q.h @@ -487,13 +487,6 @@ private: #ifndef DOXYGEN -template<> MappingQ<1>::MappingQ (const unsigned int, - const bool); -template<> MappingQ<1>::~MappingQ (); - -template<> -void MappingQ<1>::compute_shapes_virtual (const std::vector > &unit_points, - MappingQ1<1>::InternalData &data) const; template <> void MappingQ<1>::set_laplace_on_quad_vector(Table<2,double> &) const; diff --git a/include/deal.II/grid/tria_boundary.h b/include/deal.II/grid/tria_boundary.h index c5d59a4ddf..145d69461d 100644 --- a/include/deal.II/grid/tria_boundary.h +++ b/include/deal.II/grid/tria_boundary.h @@ -465,12 +465,6 @@ Point<3> StraightBoundary<3,3>:: get_new_point_on_quad (const Triangulation<3,3>::quad_iterator &quad) const; -template <> -void -StraightBoundary<1,1>:: -get_intermediate_points_on_line (const Triangulation<1,1>::line_iterator &, - std::vector > &) const; - template <> void StraightBoundary<3,3>:: diff --git a/source/fe/mapping_q.cc b/source/fe/mapping_q.cc index 43bc3f7814..79d53e46d9 100644 --- a/source/fe/mapping_q.cc +++ b/source/fe/mapping_q.cc @@ -57,44 +57,6 @@ MappingQ::InternalData::memory_consumption () const -// in 1d, it is irrelevant which polynomial degree to use, since all -// cells are scaled linearly. Unless codimension is equal to two -template<> -MappingQ<1>::MappingQ (const unsigned int, - const bool /*use_mapping_q_on_all_cells*/) - : - degree(1), - n_inner(0), - n_outer(0), - tensor_pols(0), - n_shape_functions(2), - renumber(0), - use_mapping_q_on_all_cells (false), - feq(degree), - line_support_points(degree+1) -{} - - -template<> -MappingQ<1>::MappingQ (const MappingQ<1> &m): - MappingQ1<1> (), - degree(1), - n_inner(0), - n_outer(0), - tensor_pols(0), - n_shape_functions(2), - renumber(0), - use_mapping_q_on_all_cells (m.use_mapping_q_on_all_cells), - feq(degree), - line_support_points(degree+1) -{} - -template<> -MappingQ<1>::~MappingQ () -{} - - - namespace { template @@ -179,16 +141,6 @@ MappingQ::~MappingQ () -template<> -void -MappingQ<1>::compute_shapes_virtual (const std::vector > &unit_points, - MappingQ1<1>::InternalData &data) const -{ - MappingQ1<1>::compute_shapes_virtual(unit_points, data); -} - - - template void MappingQ::compute_shapes_virtual (const std::vector > &unit_points, diff --git a/source/grid/tria_boundary.cc b/source/grid/tria_boundary.cc index e32916554f..e68d058cd5 100644 --- a/source/grid/tria_boundary.cc +++ b/source/grid/tria_boundary.cc @@ -335,42 +335,6 @@ get_new_point_on_quad (const Triangulation<3>::quad_iterator &quad) const -template <> -void -StraightBoundary<1>:: -get_intermediate_points_on_line (const Triangulation<1>::line_iterator &, - std::vector > &) const -{ - Assert(false, ExcImpossibleInDim(1)); -} - -template <> -void -StraightBoundary<1, 2>:: -get_intermediate_points_on_line (const Triangulation<1, 2>::line_iterator &line, - std::vector > &points) const -{ - const unsigned int spacedim = 2; - const unsigned int n=points.size(); - Assert(n>0, ExcInternalError()); - - // Use interior points of QGaussLobatto quadrature formula support points - // for consistency with MappingQ - const std::vector > &line_points = this->get_line_support_points(n); - const Point vertices[2] = { line->vertex(0), - line->vertex(1) - }; - - for (unsigned int i=0; i void StraightBoundary:: diff --git a/tests/bits/mapping_q_eulerian_06.cc b/tests/bits/mapping_q_eulerian_06.cc new file mode 100644 index 0000000000..4d627e153f --- /dev/null +++ b/tests/bits/mapping_q_eulerian_06.cc @@ -0,0 +1,91 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2005 - 2015 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. +// +// --------------------------------------------------------------------- + + + +// computes points in real space for 1D Eulerian mapping where middle points +// are moved + +#include "../tests.h" +#include +#include + +// all include files you need here + +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +std::ofstream logfile("output"); + +void test(unsigned int degree) +{ + const unsigned int dim = 1; + const unsigned int spacedim = 1; + Triangulation tria; + GridGenerator::hyper_cube(tria, 0, 1); + FE_Q fe(degree); + + DoFHandler shift_dh(tria); + + shift_dh.distribute_dofs(fe); + + // Shift just the interior points but not the boundary points + Vector shift(shift_dh.n_dofs()); + for (unsigned int i=2; i<=degree; ++i) + shift(i) = 0.1; + + QGauss quad(degree+1); + MappingQEulerian,spacedim> mapping(degree,shift, shift_dh); + + Triangulation::active_cell_iterator cell=tria.begin_active(), + endc=tria.end() ; + Point real; + Point unit; + double eps = 1e-10; + for (; cell!=endc; ++cell) + { + deallog< " << real << std::endl; + if ( (unit-quad.point(q)).norm()>eps) + deallog<