Fix of MappingQ<dim>::transform_real_to_unit_cell. Take MappingQ1::transform_real_to_...

git-svn-id: https://svn.dealii.org/trunk@5266 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/deal.II/source/fe/mapping_q.cc b/deal.II/deal.II/source/fe/mapping_q.cc
index 37bfe641b574403d3e7dd568ea8a9f660869885f..cf1db84ddf5e19197bace131afeace61499d6082 100644 (file)
--- a/deal.II/deal.II/source/fe/mapping_q.cc
+++ b/deal.II/deal.II/source/fe/mapping_q.cc
@@ -1232,38 +1232,29 @@ Point<dim> MappingQ<dim>::transform_real_to_unit_cell (
   const typename Triangulation<dim>::cell_iterator cell,
   const Point<dim> &p) const
 {
-                                  // Let the start value of the
-                                  // newton iteration be the center
-                                  // of the unit cell
-  Point<dim> p_unit;
-  for (unsigned int i=0; i<dim; ++i)
-    p_unit(i)=0.5;
-
-                                  // Use the get_data function to
-                                  // create an InternalData with data
-                                  // vectors of the right size and
-                                  // transformation shape values and
-                                  // derivatives already computed at
-                                  // point p_unit.
-  const Quadrature<dim> point_quadrature(p_unit);
-  InternalData *mdata=dynamic_cast<InternalData *> (
-    get_data(update_transformation_values | update_transformation_gradients,
-            point_quadrature));
-  Assert(mdata!=0, ExcInternalError());
-
-  mdata->use_mapping_q1_on_current_cell = !(use_mapping_q_on_all_cells
-                                            || cell->has_boundary_lines());
-
-  typename MappingQ1<dim>::InternalData *p_data=0;
-  if (mdata->use_mapping_q1_on_current_cell)
-    p_data=&mdata->mapping_q1_data;
-  else
-    p_data=mdata;
+                                  // first a Newton iteration based
+                                  // on a Q1 mapping
+  Point<dim> p_unit=MappingQ1<dim>::transform_real_to_unit_cell(cell, p);
+  
+  if (cell->has_boundary_lines() || use_mapping_q_on_all_cells)
+    {
+                                      // then a Newton iteration
+                                      // based on the full MappingQ
+      const Quadrature<dim> point_quadrature(p_unit);
+      InternalData *mdata=dynamic_cast<InternalData *> (
+       get_data(update_transformation_values | update_transformation_gradients,
+                point_quadrature));
+      Assert(mdata!=0, ExcInternalError());
+      mdata->use_mapping_q1_on_current_cell=false;
+
+      std::vector<Point<dim> > &points=mdata->mapping_support_points;
+      compute_mapping_support_points(cell, points);
+
+      transform_real_to_unit_cell_internal(cell, p, *mdata, p_unit);
   
-                                  // perform the newton iteration.
-  transform_real_to_unit_cell_internal(cell, p, *p_data, p_unit);
+      delete mdata;
+    }
   
-  delete mdata;
   return p_unit;
 }
 

diff --git a/deal.II/deal.II/source/fe/mapping_q1.cc b/deal.II/deal.II/source/fe/mapping_q1.cc
index 05c1eda9b7b3b5bd0f83cd7b7521fdd51479cce0..e0d0c50b98b6c7839961da6baa35027e9fdf6307 100644 (file)
--- a/deal.II/deal.II/source/fe/mapping_q1.cc
+++ b/deal.II/deal.II/source/fe/mapping_q1.cc
@@ -930,10 +930,14 @@ Point<dim> MappingQ1<dim>::transform_real_to_unit_cell (
                                   // point p_unit.
   const Quadrature<dim> point_quadrature(p_unit);
   InternalData *mdata=dynamic_cast<InternalData *> (
-    get_data(update_transformation_values | update_transformation_gradients,
-            point_quadrature));
+    MappingQ1<dim>::get_data(update_transformation_values
+                            | update_transformation_gradients,
+                            point_quadrature));
   Assert(mdata!=0, ExcInternalError());
 
+  MappingQ1<dim>::compute_mapping_support_points(cell, mdata->mapping_support_points);
+  Assert(mdata->mapping_support_points.size()==4, ExcInternalError());
+
                                   // perform the newton iteration.
   transform_real_to_unit_cell_internal(cell, p, *mdata, p_unit);
   
@@ -955,7 +959,6 @@ void MappingQ1<dim>::transform_real_to_unit_cell_internal (
   Assert(mdata.shape_derivatives.size()==n_shapes, ExcInternalError());
   
   std::vector<Point<dim> > &points=mdata.mapping_support_points;
-  compute_mapping_support_points(cell, points);
   Assert(points.size()==n_shapes, ExcInternalError());
   
                                   // Newton iteration to solve