]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Removed smart pointer to dh.get_fe() 5720/head
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 11 Jan 2018 09:18:09 +0000 (10:18 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 11 Jan 2018 09:18:44 +0000 (10:18 +0100)
include/deal.II/fe/mapping_fe_field.h
source/fe/mapping_fe_field.cc
tests/mappings/mapping_fe_field_02.cc [new file with mode: 0644]
tests/mappings/mapping_fe_field_02.output [new file with mode: 0644]

index 95fb0938ecff01da66c1ecc540e0d948638f03f9..6067d80f65732a0228487b1d726ddb17dcce6fc1 100644 (file)
@@ -496,16 +496,6 @@ private:
    */
   SmartPointer<const VectorType, MappingFEField<dim,spacedim,VectorType,DoFHandlerType> > euler_vector;
 
-  /**
-   * A FiniteElement object which is only needed in 3D, since it knows how to
-   * reorder shape functions/DoFs on non-standard faces. This is used to
-   * reorder support points in the same way. We could make this a pointer to
-   * prevent construction in 1D and 2D, but since memory and time requirements
-   * are not particularly high this seems unnecessary at the moment.
-   */
-  SmartPointer<const FiniteElement<dim,spacedim>, MappingFEField<dim,spacedim,VectorType,DoFHandlerType> > fe;
-
-
   /**
    * Pointer to the DoFHandler to which the mapping vector is associated.
    */
index 1f2e35e624a6fcdcdf628ad96e3cd0e8e9f6247c..a06dcf75c66e0da2f8e6798395f855ca89d6b86e 100644 (file)
@@ -207,10 +207,9 @@ MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::MappingFEField
  const ComponentMask  &mask)
   :
   euler_vector(&euler_vector),
-  fe(&euler_dof_handler.get_fe()),
   euler_dof_handler(&euler_dof_handler),
   fe_mask(mask.size() ? mask :
-          ComponentMask(fe->get_nonzero_components(0).size(), true)),
+          ComponentMask(euler_dof_handler.get_fe().get_nonzero_components(0).size(), true)),
   fe_to_real(fe_mask.size(), numbers::invalid_unsigned_int)
 {
   unsigned int size = 0;
@@ -228,7 +227,6 @@ MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::MappingFEField
 (const MappingFEField<dim,spacedim,VectorType,DoFHandlerType> &mapping)
   :
   euler_vector(mapping.euler_vector),
-  fe(mapping.fe),
   euler_dof_handler(mapping.euler_dof_handler),
   fe_mask(mapping.fe_mask),
   fe_to_real(mapping.fe_to_real)
@@ -266,6 +264,7 @@ MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::
 compute_shapes_virtual (const std::vector<Point<dim> >                    &unit_points,
                         typename MappingFEField<dim, spacedim, VectorType, DoFHandlerType>::InternalData &data) const
 {
+  const auto fe = &euler_dof_handler->get_fe();
   const unsigned int n_points=unit_points.size();
 
   for (unsigned int point=0; point<n_points; ++point)
@@ -480,7 +479,7 @@ MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::InternalData *
 MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::get_data (const UpdateFlags      update_flags,
     const Quadrature<dim> &quadrature) const
 {
-  InternalData *data = new InternalData(*fe, fe_mask);
+  InternalData *data = new InternalData(euler_dof_handler->get_fe(), fe_mask);
   this->compute_data (update_flags, quadrature,
                       quadrature.size(), *data);
   return data;
@@ -494,7 +493,7 @@ MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::get_face_data
 (const UpdateFlags        update_flags,
  const Quadrature<dim-1> &quadrature) const
 {
-  InternalData *data = new InternalData(*fe, fe_mask);
+  InternalData *data = new InternalData(euler_dof_handler->get_fe(), fe_mask);
   const Quadrature<dim> q (QProjector<dim>::project_to_all_faces(quadrature));
   this->compute_face_data (update_flags, q,
                            quadrature.size(), *data);
@@ -509,7 +508,7 @@ MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::get_subface_data
 (const UpdateFlags        update_flags,
  const Quadrature<dim-1> &quadrature) const
 {
-  InternalData *data = new InternalData(*fe, fe_mask);
+  InternalData *data = new InternalData(euler_dof_handler->get_fe(), fe_mask);
   const Quadrature<dim> q (QProjector<dim>::project_to_all_subfaces(quadrature));
   this->compute_face_data (update_flags, q,
                            quadrature.size(), *data);
@@ -1324,13 +1323,13 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
 
   internal::MappingFEField::maybe_compute_q_points<dim,spacedim,VectorType,DoFHandlerType>
   (QProjector<dim>::DataSetDescriptor::cell (),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.quadrature_points);
 
   internal::MappingFEField::maybe_update_Jacobians<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell (),
-   data, *fe, fe_mask, fe_to_real);
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real);
 
   const UpdateFlags update_flags = data.update_each;
   const std::vector<double> &weights=quadrature.get_weights();
@@ -1434,35 +1433,35 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
   internal::MappingFEField::maybe_update_jacobian_grads<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell(),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.jacobian_grads);
 
   // calculate derivatives of the Jacobians pushed forward to real cell coordinates
   internal::MappingFEField::maybe_update_jacobian_pushed_forward_grads<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell(),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.jacobian_pushed_forward_grads);
 
   // calculate hessians of the Jacobians
   internal::MappingFEField::maybe_update_jacobian_2nd_derivatives<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell(),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.jacobian_2nd_derivatives);
 
   // calculate hessians of the Jacobians pushed forward to real cell coordinates
   internal::MappingFEField::maybe_update_jacobian_pushed_forward_2nd_derivatives<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell(),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.jacobian_pushed_forward_2nd_derivatives);
 
   // calculate gradients of the hessians of the Jacobians
   internal::MappingFEField::maybe_update_jacobian_3rd_derivatives<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell(),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.jacobian_3rd_derivatives);
 
   // calculate gradients of the hessians of the Jacobians pushed forward to real
@@ -1470,7 +1469,7 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
   internal::MappingFEField::maybe_update_jacobian_pushed_forward_3rd_derivatives<dim,spacedim,VectorType,DoFHandlerType>
   (cell_similarity,
    QProjector<dim>::DataSetDescriptor::cell(),
-   data, *fe, fe_mask, fe_to_real,
+   data, euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data.jacobian_pushed_forward_3rd_derivatives);
 
   return updated_cell_similarity;
@@ -1506,7 +1505,7 @@ fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &
          quadrature.size()),
    quadrature,
    data,
-   *fe, fe_mask, fe_to_real,
+   euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data);
 }
 
@@ -1541,7 +1540,7 @@ fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterato
             cell->subface_case(face_no)),
    quadrature,
    data,
-   *fe, fe_mask, fe_to_real,
+   euler_dof_handler->get_fe(), fe_mask, fe_to_real,
    output_data);
 }
 
@@ -1797,7 +1796,7 @@ do_transform_unit_to_real_cell (const InternalData &data) const
 
   for (unsigned int i=0; i<data.n_shape_functions; ++i)
     {
-      unsigned int comp_i = fe->system_to_component_index(i).first;
+      unsigned int comp_i = euler_dof_handler->get_fe().system_to_component_index(i).first;
       if (fe_mask[comp_i])
         p_real[fe_to_real[comp_i]] += data.local_dof_values[i] * data.shape(0,i);
     }
@@ -1890,7 +1889,7 @@ do_transform_real_to_unit_cell
       for (unsigned int k=0; k<mdata.n_shape_functions; ++k)
         {
           const Tensor<1,dim> &grad_k = mdata.derivative(0,k);
-          unsigned int comp_k = fe->system_to_component_index(k).first;
+          unsigned int comp_k = euler_dof_handler->get_fe().system_to_component_index(k).first;
           if (fe_mask[comp_k])
             for (unsigned int j=0; j<dim; ++j)
               DF[j][fe_to_real[comp_k]] += mdata.local_dof_values[k] * grad_k[j];
@@ -1962,7 +1961,7 @@ template <int dim, int spacedim, typename VectorType, typename DoFHandlerType>
 unsigned int
 MappingFEField<dim,spacedim,VectorType,DoFHandlerType>::get_degree() const
 {
-  return fe->degree;
+  return euler_dof_handler->get_fe().degree;
 }
 
 
diff --git a/tests/mappings/mapping_fe_field_02.cc b/tests/mappings/mapping_fe_field_02.cc
new file mode 100644 (file)
index 0000000..b7a1599
--- /dev/null
@@ -0,0 +1,69 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 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.
+//
+// ---------------------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/geometry_info.h>
+#include<deal.II/dofs/dof_handler.h>
+#include<deal.II/dofs/dof_accessor.h>
+#include<deal.II/dofs/dof_tools.h>
+#include<deal.II/fe/fe_q.h>
+#include<deal.II/fe/fe_values.h>
+#include<deal.II/fe/fe_system.h>
+#include<deal.II/fe/mapping_q1_eulerian.h>
+#include<deal.II/fe/mapping_q1.h>
+#include <deal.II/fe/mapping_fe_field.h>
+#include <deal.II/grid/grid_generator.h>
+#include<deal.II/numerics/vector_tools.h>
+
+template<int dim, int spacedim>
+void test()
+{
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_cube (tria);
+
+  FESystem<dim,spacedim> fe(FE_Q<dim,spacedim>(1), spacedim);
+  DoFHandler<dim,spacedim> dh(tria);
+
+  dh.distribute_dofs(fe);
+
+  deallog << "dim, spacedim: " << dim << ", " << spacedim << std::endl
+          << "cells: " << tria.n_active_cells() << ", dofs: "
+          << dh.n_dofs() <<std::endl;
+
+  // Create a Mapping
+  Vector<double> map_vector(dh.n_dofs());
+  VectorTools::get_position_vector(dh, map_vector);
+  MappingFEField<dim,spacedim> mapping(dh, map_vector);
+
+  tria.refine_global (1);
+
+  dh.distribute_dofs(fe);
+
+  deallog << "After refine:" << std::endl
+          << "cells: " << tria.n_active_cells() << ", dofs: "
+          << dh.n_dofs() <<std::endl;
+}
+
+using namespace dealii;
+int main()
+{
+  initlog();
+  test<1,1>();
+  test<1,2>();
+  test<2,2>();
+  test<2,3>();
+  test<3,3>();
+}
diff --git a/tests/mappings/mapping_fe_field_02.output b/tests/mappings/mapping_fe_field_02.output
new file mode 100644 (file)
index 0000000..562627c
--- /dev/null
@@ -0,0 +1,21 @@
+
+DEAL::dim, spacedim: 1, 1
+DEAL::cells: 1, dofs: 2
+DEAL::After refine:
+DEAL::cells: 2, dofs: 3
+DEAL::dim, spacedim: 1, 2
+DEAL::cells: 1, dofs: 4
+DEAL::After refine:
+DEAL::cells: 2, dofs: 6
+DEAL::dim, spacedim: 2, 2
+DEAL::cells: 1, dofs: 8
+DEAL::After refine:
+DEAL::cells: 4, dofs: 18
+DEAL::dim, spacedim: 2, 3
+DEAL::cells: 1, dofs: 12
+DEAL::After refine:
+DEAL::cells: 4, dofs: 27
+DEAL::dim, spacedim: 3, 3
+DEAL::cells: 1, dofs: 24
+DEAL::After refine:
+DEAL::cells: 8, dofs: 81

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.