]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix mapping cartesian 5313/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Tue, 24 Oct 2017 15:40:05 +0000 (09:40 -0600)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Tue, 24 Oct 2017 16:32:51 +0000 (10:32 -0600)
doc/news/changes/minor/20171024ReneGassmoeller [new file with mode: 0644]
source/fe/mapping_cartesian.cc
tests/mappings/mapping_cartesian_02.cc [new file with mode: 0644]
tests/mappings/mapping_cartesian_02.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20171024ReneGassmoeller b/doc/news/changes/minor/20171024ReneGassmoeller
new file mode 100644 (file)
index 0000000..18bb2c1
--- /dev/null
@@ -0,0 +1,6 @@
+Fixed: The computation of the inverse jacobians in the
+MappingCartesian class was broken. It would trigger an
+assertion at runtime complaining that only an assignment
+with zero is allowed. This is fixed now.
+<br>
+(Rene Gassmoeller, 2017/10/24)
index 9b42228ad2f4966313ac7e3a096bb87682e9c77c..b6e9ef78d87ec5c0a2af14c736b7eda5cd63541a 100644 (file)
@@ -420,7 +420,7 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
         {
           output_data.inverse_jacobians[i] = Tensor<2,dim>();
           for (unsigned int j=0; j<dim; ++j)
-            output_data.inverse_jacobians[j][j]=1./data.cell_extents[j];
+            output_data.inverse_jacobians[i][j][j]=1./data.cell_extents[j];
         }
 
   return cell_similarity;
diff --git a/tests/mappings/mapping_cartesian_02.cc b/tests/mappings/mapping_cartesian_02.cc
new file mode 100644 (file)
index 0000000..f720ce3
--- /dev/null
@@ -0,0 +1,96 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2003 - 2016 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// there used to be a bug in MappingCartesian, where the inverse jacobian
+// was not computed correctly. Check that all works as intended.
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping_cartesian.h>
+
+
+
+
+template <int dim>
+void check (const Triangulation<dim> &tria)
+{
+  MappingCartesian<dim> mapping;
+  FE_Q<dim> fe(1);
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+
+  QGauss<dim> quadrature(1);
+
+  UpdateFlags update_flags = update_quadrature_points |
+                             update_JxW_values |
+                             update_jacobians |
+                             update_jacobian_grads |
+                             update_jacobian_pushed_forward_grads |
+                             update_jacobian_2nd_derivatives |
+                             update_jacobian_pushed_forward_2nd_derivatives |
+                             update_jacobian_3rd_derivatives |
+                             update_jacobian_pushed_forward_3rd_derivatives |
+                             update_inverse_jacobians |
+                             // Transformation dependence
+                             update_covariant_transformation |
+                             update_contravariant_transformation |
+                             update_transformation_values |
+                             update_transformation_gradients |
+                             // Volume data
+                             update_volume_elements;
+
+  FEValues<dim>    fe_values (mapping, fe, quadrature,
+                              update_flags);
+
+  fe_values.reinit (dof_handler.begin_active());
+
+  for (unsigned int i=0; i<dim; ++i)
+    deallog << fe_values.inverse_jacobian(0)[i] << std::endl;
+
+  deallog << std::endl;
+
+  deallog << "OK" << std::endl;
+}
+
+
+int main ()
+{
+  initlog();
+
+  {
+    Triangulation<2> coarse_grid;
+    GridGenerator::hyper_cube (coarse_grid);
+    check (coarse_grid);
+  }
+
+  {
+    Triangulation<3> coarse_grid;
+    GridGenerator::hyper_cube (coarse_grid);
+    check (coarse_grid);
+  }
+}
+
+
+
diff --git a/tests/mappings/mapping_cartesian_02.output b/tests/mappings/mapping_cartesian_02.output
new file mode 100644 (file)
index 0000000..68b7984
--- /dev/null
@@ -0,0 +1,10 @@
+
+DEAL::1.00000 0.00000
+DEAL::0.00000 1.00000
+DEAL::
+DEAL::OK
+DEAL::1.00000 0.00000 0.00000
+DEAL::0.00000 1.00000 0.00000
+DEAL::0.00000 0.00000 1.00000
+DEAL::
+DEAL::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.