]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use ReferenceCell instead of StaticMappingQ1.
authorDavid Wells <drwells@email.unc.edu>
Tue, 16 Aug 2022 21:16:45 +0000 (17:16 -0400)
committerDavid Wells <drwells@email.unc.edu>
Thu, 18 Aug 2022 14:39:11 +0000 (10:39 -0400)
source/grid/grid_tools.cc
source/grid/tria_accessor.cc
source/numerics/data_out_rotation.cc
source/numerics/derivative_approximation.cc
source/numerics/error_estimator_1d.cc
source/numerics/point_value_history.cc

index de63f709d97d7d419ba2b7ec8394c57f202c9fb2..601c9c9bbe171f9b6f158aea8ce9027a44f1810f 100644 (file)
@@ -32,7 +32,6 @@
 #include <deal.II/fe/fe_q.h>
 #include <deal.II/fe/fe_values.h>
 #include <deal.II/fe/mapping_q.h>
-#include <deal.II/fe/mapping_q1.h>
 
 #include <deal.II/grid/filtered_iterator.h>
 #include <deal.II/grid/grid_reordering.h>
@@ -2155,8 +2154,15 @@ namespace GridTools
 
     QGauss<dim> quadrature(4);
 
+    Assert(triangulation.all_reference_cells_are_hyper_cube(),
+           ExcNotImplemented());
+    const auto reference_cell = ReferenceCells::get_hypercube<dim>();
     MatrixCreator::create_laplace_matrix(
-      StaticMappingQ1<dim>::mapping, dof_handler, quadrature, S, coefficient);
+      reference_cell.template get_default_linear_mapping<dim, dim>(),
+      dof_handler,
+      quadrature,
+      S,
+      coefficient);
 
     // set up the boundary values for the laplace problem
     std::array<AffineConstraints<double>, dim>                  constraints;
index 83f0d9e8e482f3fd4e97d2716a68c1b6e94c6be3..845a6d489b33bab916eb8390ff731fd7a73a9ac3 100644 (file)
@@ -2157,9 +2157,10 @@ CellAccessor<dim, spacedim>::point_inside_codim(const Point<spacedim_> &p) const
   Assert(this->reference_cell().is_hyper_cube(), ExcNotImplemented());
 
   const TriaRawIterator<CellAccessor<dim_, spacedim_>> cell_iterator(*this);
-  const Point<dim_>                                    p_unit =
-    StaticMappingQ1<dim_, spacedim_>::mapping.transform_real_to_unit_cell(
-      cell_iterator, p);
+
+  const Point<dim_> p_unit =
+    this->reference_cell().template get_default_linear_mapping<dim_, spacedim_>()
+      .transform_real_to_unit_cell(cell_iterator, p);
 
   return GeometryInfo<dim_>::is_inside_unit_cell(p_unit);
 }
index 16b34d3fc2c150c7abd90f655406bfbfdaeeca03..5577903b26e40acfd6fd6a475823b7a529fc8f63 100644 (file)
@@ -511,14 +511,18 @@ DataOutRotation<dim, spacedim>::build_patches(
     else
       n_postprocessor_outputs[dataset] = 0;
 
+  Assert(!this->triangulation->is_mixed_mesh(), ExcNotImplemented());
+  const auto reference_cell =
+    this->triangulation->get_reference_cells()[0];
   internal::DataOutRotationImplementation::ParallelData<dim, spacedim>
-    thread_data(n_datasets,
-                n_subdivisions,
-                n_patches_per_circle,
-                n_postprocessor_outputs,
-                StaticMappingQ1<dim, spacedim>::mapping,
-                this->get_fes(),
-                update_flags);
+    thread_data(
+      n_datasets,
+      n_subdivisions,
+      n_patches_per_circle,
+      n_postprocessor_outputs,
+      reference_cell.template get_default_linear_mapping<dim, spacedim>(),
+      this->get_fes(),
+      update_flags);
   std::vector<DataOutBase::Patch<patch_dim, patch_spacedim>> new_patches(
     n_patches_per_circle);
   for (unsigned int i = 0; i < new_patches.size(); ++i)
index d7d01a75a92c523081a6f8aa4a64e53f0184dc76..30930b0313da146ea9055c224732b95ff8e54c60 100644 (file)
@@ -21,7 +21,6 @@
 
 #include <deal.II/fe/fe.h>
 #include <deal.II/fe/fe_values.h>
-#include <deal.II/fe/mapping_q1.h>
 
 #include <deal.II/grid/filtered_iterator.h>
 #include <deal.II/grid/grid_tools.h>
@@ -1030,8 +1029,12 @@ namespace DerivativeApproximation
                        Vector<float> &                  derivative_norm,
                        const unsigned int               component)
   {
+    Assert(!dof_handler.get_triangulation().is_mixed_mesh(),
+           ExcNotImplemented());
+    const auto reference_cell =
+      dof_handler.get_triangulation().get_reference_cells()[0];
     internal::approximate_derivative<internal::Gradient<dim>, dim>(
-      StaticMappingQ1<dim>::mapping,
+      reference_cell.template get_default_linear_mapping<dim, spacedim>(),
       dof_handler,
       solution,
       component,
@@ -1059,8 +1062,12 @@ namespace DerivativeApproximation
                                 Vector<float> &    derivative_norm,
                                 const unsigned int component)
   {
+    Assert(!dof_handler.get_triangulation().is_mixed_mesh(),
+           ExcNotImplemented());
+    const auto reference_cell =
+      dof_handler.get_triangulation().get_reference_cells()[0];
     internal::approximate_derivative<internal::SecondDerivative<dim>, dim>(
-      StaticMappingQ1<dim>::mapping,
+      reference_cell.template get_default_linear_mapping<dim, spacedim>(),
       dof_handler,
       solution,
       component,
index 413ce5d18b2147f2af237111ac9695d303f83aaf..0d085e152c697a9ef25b8b77db3de834a0c2012b 100644 (file)
@@ -111,7 +111,8 @@ KellyErrorEstimator<1, spacedim>::estimate(
   const types::material_id  material_id,
   const Strategy            strategy)
 {
-  estimate(StaticMappingQ1<1, spacedim>::mapping,
+  const auto reference_cell = ReferenceCells::Line;
+  estimate(reference_cell.template get_default_linear_mapping<1, spacedim>(),
            dof_handler,
            quadrature,
            neumann_bc,
@@ -145,7 +146,8 @@ KellyErrorEstimator<1, spacedim>::estimate(
   const types::material_id                material_id,
   const Strategy                          strategy)
 {
-  estimate(StaticMappingQ1<1, spacedim>::mapping,
+  const auto reference_cell = ReferenceCells::Line;
+  estimate(reference_cell.template get_default_linear_mapping<1, spacedim>(),
            dof_handler,
            quadrature,
            neumann_bc,
@@ -217,7 +219,8 @@ KellyErrorEstimator<1, spacedim>::estimate(
   const types::material_id  material_id,
   const Strategy            strategy)
 {
-  estimate(StaticMappingQ1<1, spacedim>::mapping,
+  const auto reference_cell = ReferenceCells::Line;
+  estimate(reference_cell.template get_default_linear_mapping<1, spacedim>(),
            dof_handler,
            quadrature,
            neumann_bc,
@@ -251,7 +254,8 @@ KellyErrorEstimator<1, spacedim>::estimate(
   const types::material_id                material_id,
   const Strategy                          strategy)
 {
-  estimate(StaticMappingQ1<1, spacedim>::mapping,
+  const auto reference_cell = ReferenceCells::Line;
+  estimate(reference_cell.template get_default_linear_mapping<1, spacedim>(),
            dof_handler,
            quadrature,
            neumann_bc,
index 98f175f8eb4f1533e37ce22f3f8d129564b1b4b6..fb95e223344617555a7d0fd0fed8b1b589307d40 100644 (file)
@@ -710,15 +710,19 @@ PointValueHistory<dim>::evaluate_field(
   typename std::vector<
     internal::PointValueHistoryImplementation::PointGeometryData<dim>>::iterator
     point = point_geometry_data.begin();
+  Assert(!dof_handler->get_triangulation().is_mixed_mesh(), ExcNotImplemented());
+  const auto reference_cell =
+    dof_handler->get_triangulation().get_reference_cells()[0];
   for (unsigned int data_store_index = 0; point != point_geometry_data.end();
        ++point, ++data_store_index)
     {
       // we now have a point to query, need to know what cell it is in
       const Point<dim> requested_location = point->requested_location;
       const typename DoFHandler<dim>::active_cell_iterator cell =
-        GridTools::find_active_cell_around_point(StaticMappingQ1<dim>::mapping,
-                                                 *dof_handler,
-                                                 requested_location)
+        GridTools::find_active_cell_around_point(
+          reference_cell.template get_default_linear_mapping<dim, dim>(),
+          *dof_handler,
+          requested_location)
           .first;
 
 
@@ -1271,6 +1275,9 @@ PointValueHistory<dim>::get_postprocessor_locations(
   typename std::vector<
     internal::PointValueHistoryImplementation::PointGeometryData<dim>>::iterator
     point = point_geometry_data.begin();
+  Assert(!dof_handler->get_triangulation().is_mixed_mesh(), ExcNotImplemented());
+  const auto reference_cell =
+    dof_handler->get_triangulation().get_reference_cells()[0];
   for (unsigned int data_store_index = 0; point != point_geometry_data.end();
        ++point, ++data_store_index)
     {
@@ -1278,9 +1285,10 @@ PointValueHistory<dim>::get_postprocessor_locations(
       // need to know what cell it is in
       Point<dim> requested_location = point->requested_location;
       typename DoFHandler<dim>::active_cell_iterator cell =
-        GridTools::find_active_cell_around_point(StaticMappingQ1<dim>::mapping,
-                                                 *dof_handler,
-                                                 requested_location)
+        GridTools::find_active_cell_around_point(
+          reference_cell.template get_default_linear_mapping<dim, dim>(),
+          *dof_handler,
+          requested_location)
           .first;
       fe_values.reinit(cell);
 

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.