From 3827c524a31d3ba90485802a31c4e80be9f86268 Mon Sep 17 00:00:00 2001
From: Wolfgang Bangerth <bangerth@math.tamu.edu>
Date: Tue, 8 Sep 2015 10:24:17 -0500
Subject: [PATCH] Override the Q1 mapping in MappingQ.

---
 source/fe/mapping_q_eulerian.cc      | 30 ++++++++++++++++++++--------
 source/fe/mapping_q_eulerian.inst.in |  3 ++-
 2 files changed, 24 insertions(+), 9 deletions(-)

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 <deal.II/fe/fe.h>
 #include <deal.II/fe/fe_tools.h>
 #include <deal.II/fe/mapping_q_eulerian.h>
+#include <deal.II/fe/mapping_q1_eulerian.h>
 
 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<dim,EulerVectorType,spacedim>(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<dim,EulerVectorType,spacedim>(euler_vector,
+                          euler_dof_handler));
+}
 
 
 
@@ -101,8 +114,8 @@ SupportQuadrature (const unsigned int map_degree)
   for (unsigned int i=1; i<dpo.size(); ++i)
     dpo[i]=dpo[i-1]*(map_degree-1);
 
-  FETools::lexicographic_to_hierarchic_numbering (
-    FiniteElementData<dim> (dpo, 1, map_degree), renumber);
+  FETools::lexicographic_to_hierarchic_numbering (FiniteElementData<dim> (dpo, 1, map_degree),
+                                                  renumber);
 
   // finally we assign the quadrature points in the required order.
   for (unsigned int q=0; q<n_q_points; ++q)
@@ -116,9 +129,8 @@ SupportQuadrature (const unsigned int map_degree)
 template <int dim, class EulerVectorType, int spacedim>
 void
 MappingQEulerian<dim, EulerVectorType, spacedim>::
-compute_mapping_support_points
-(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- std::vector<Point<spacedim> > &a) const
+compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                                std::vector<Point<spacedim> > &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<Vector<double> > shift_vector(n_support_pts,Vector<double>(n_components));
+  std::vector<Vector<typename EulerVectorType::value_type> >
+  shift_vector(n_support_pts,
+               Vector<typename EulerVectorType::value_type>(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_dimension, Vector<double>, deal_II_space_dimension>;
+    template class MappingQEulerian<deal_II_dimension, Vector<float>,  deal_II_space_dimension>;
 #  ifdef DEAL_II_WITH_PETSC
     template class MappingQEulerian<deal_II_dimension,
 				    PETScWrappers::Vector, deal_II_space_dimension>;
-- 
2.39.5