From: Wolfgang Bangerth Date: Tue, 8 Sep 2015 15:24:17 +0000 (-0500) Subject: Override the Q1 mapping in MappingQ. X-Git-Tag: v8.4.0-rc2~439^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3827c524a31d3ba90485802a31c4e80be9f86268;p=dealii.git Override the Q1 mapping in MappingQ. --- diff --git a/source/fe/mapping_q_eulerian.cc b/source/fe/mapping_q_eulerian.cc index 67f069b640..c7a91efb09 100644 --- a/source/fe/mapping_q_eulerian.cc +++ b/source/fe/mapping_q_eulerian.cc @@ -26,6 +26,7 @@ #include #include #include +#include DEAL_II_NAMESPACE_OPEN @@ -46,7 +47,13 @@ MappingQEulerian (const unsigned int degree, fe_values(euler_dof_handler.get_fe(), support_quadrature, update_values | update_q_points) -{ } +{ + // reset the q1 mapping we use for interior cells (and previously + // set by the MappingQ constructor) to a MappingQ1Eulerian with the + // current vector + this->q1_mapping.reset (new MappingQ1Eulerian(euler_vector, + euler_dof_handler)); +} @@ -63,7 +70,13 @@ MappingQEulerian (const unsigned int degree, fe_values(euler_dof_handler.get_fe(), support_quadrature, update_values | update_q_points) -{ } +{ + // reset the q1 mapping we use for interior cells (and previously + // set by the MappingQ constructor) to a MappingQ1Eulerian with the + // current vector + this->q1_mapping.reset (new MappingQ1Eulerian(euler_vector, + euler_dof_handler)); +} @@ -101,8 +114,8 @@ SupportQuadrature (const unsigned int map_degree) for (unsigned int i=1; i (dpo, 1, map_degree), renumber); + FETools::lexicographic_to_hierarchic_numbering (FiniteElementData (dpo, 1, map_degree), + renumber); // finally we assign the quadrature points in the required order. for (unsigned int q=0; q void MappingQEulerian:: -compute_mapping_support_points -(const typename Triangulation::cell_iterator &cell, - std::vector > &a) const +compute_mapping_support_points (const typename Triangulation::cell_iterator &cell, + std::vector > &a) const { // first, basic assertion with respect to vector size, @@ -153,7 +165,9 @@ compute_mapping_support_points Assert (n_components >= spacedim, ExcDimensionMismatch(n_components, spacedim) ); - std::vector > shift_vector(n_support_pts,Vector(n_components)); + std::vector > + shift_vector(n_support_pts, + Vector(n_components)); // fill shift vector for each support point using an fe_values object. make // sure that the fe_values variable isn't used simultaneously from different diff --git a/source/fe/mapping_q_eulerian.inst.in b/source/fe/mapping_q_eulerian.inst.in index caf91c5d5b..dd6397f73b 100644 --- a/source/fe/mapping_q_eulerian.inst.in +++ b/source/fe/mapping_q_eulerian.inst.in @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 1998 - 2014 by the deal.II authors +// Copyright (C) 1998 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -18,6 +18,7 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS { #if deal_II_dimension <= deal_II_space_dimension template class MappingQEulerian, deal_II_space_dimension>; + template class MappingQEulerian, deal_II_space_dimension>; # ifdef DEAL_II_WITH_PETSC template class MappingQEulerian;