From: David Wells Date: Sat, 18 Jul 2015 03:54:21 +0000 (-0400) Subject: Use a separate function to extract vertices. X-Git-Tag: v8.4.0-rc2~667^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a13893459431d81de0a5b4fc6f49699351fa70e2;p=dealii.git Use a separate function to extract vertices. --- diff --git a/include/deal.II/fe/mapping.h b/include/deal.II/fe/mapping.h index 6fc969415f..144e769594 100644 --- a/include/deal.II/fe/mapping.h +++ b/include/deal.II/fe/mapping.h @@ -19,6 +19,7 @@ #include #include +#include #include #include #include @@ -230,6 +231,16 @@ public: */ virtual ~Mapping (); + /** + * Return the mapped vertices of a cell. These values are not equal to the + * vertex coordinates stored by the triangulation for MappingQEulerian and + * MappingQ1Eulerian. + */ + virtual + std_cxx11::array, GeometryInfo::vertices_per_cell> + get_vertices ( + const typename Triangulation::cell_iterator &cell) const; + /** * Transforms the point @p p on the unit cell to the point @p p_real on the * real cell @p cell and returns @p p_real. diff --git a/include/deal.II/fe/mapping_q1_eulerian.h b/include/deal.II/fe/mapping_q1_eulerian.h index 8b31afdc44..06fdb8c7be 100644 --- a/include/deal.II/fe/mapping_q1_eulerian.h +++ b/include/deal.II/fe/mapping_q1_eulerian.h @@ -17,6 +17,7 @@ #define dealii__mapping_q1_eulerian_h #include +#include #include #include @@ -99,6 +100,14 @@ public: MappingQ1Eulerian (const VECTOR &euler_transform_vectors, const DoFHandler &shiftmap_dof_handler); + /** + * Return the mapped vertices of the cell. + */ + virtual + std_cxx11::array, GeometryInfo::vertices_per_cell> + get_vertices + (const typename Triangulation::cell_iterator &cell) const; + /** * Return a pointer to a copy of the present object. The caller of this copy * then assumes ownership of it. @@ -146,17 +155,6 @@ protected: * Pointer to the DoFHandler to which the mapping vector is associated. */ SmartPointer,MappingQ1Eulerian > shiftmap_dof_handler; - - -private: - /** - * Computes the support points of the mapping. For @p MappingQ1Eulerian - * these are the vertices. - */ - virtual void compute_mapping_support_points( - const typename Triangulation::cell_iterator &cell, - std::vector > &a) const; - }; /*@}*/ diff --git a/source/fe/mapping.cc b/source/fe/mapping.cc index af4a2c4359..6160a21dba 100644 --- a/source/fe/mapping.cc +++ b/source/fe/mapping.cc @@ -25,6 +25,21 @@ Mapping::~Mapping () {} +template +std_cxx11::array, GeometryInfo::vertices_per_cell> +Mapping::get_vertices ( + const typename Triangulation::cell_iterator &cell) const +{ + std_cxx11::array, GeometryInfo::vertices_per_cell> vertices; + for (unsigned int i = 0; i < GeometryInfo::vertices_per_cell; ++i) + { + vertices[i] = cell->vertex(i); + } + return vertices; +} + + + /*------------------------------ InternalDataBase ------------------------------*/ diff --git a/source/fe/mapping_q1.cc b/source/fe/mapping_q1.cc index 638d46d588..210f01965f 100644 --- a/source/fe/mapping_q1.cc +++ b/source/fe/mapping_q1.cc @@ -18,6 +18,7 @@ #include #include #include +#include #include #include #include @@ -725,10 +726,12 @@ MappingQ1::compute_mapping_support_points( const typename Triangulation::cell_iterator &cell, std::vector > &a) const { + std_cxx11::array, GeometryInfo::vertices_per_cell> + vertices = this->get_vertices(cell); a.resize(GeometryInfo::vertices_per_cell); for (unsigned int i=0; i::vertices_per_cell; ++i) - a[i] = cell->vertex(i); + a[i] = vertices[i]; } diff --git a/source/fe/mapping_q1_eulerian.cc b/source/fe/mapping_q1_eulerian.cc index 0903bd44da..77101747d1 100644 --- a/source/fe/mapping_q1_eulerian.cc +++ b/source/fe/mapping_q1_eulerian.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2001 - 2014 by the deal.II authors +// Copyright (C) 2001 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -13,6 +13,7 @@ // // --------------------------------------------------------------------- +#include #include #include #include @@ -36,67 +37,50 @@ MappingQ1Eulerian (const EulerVectorType &euler_transform_vectors, template -void +std_cxx11::array, GeometryInfo::vertices_per_cell> MappingQ1Eulerian:: -compute_mapping_support_points(const typename Triangulation::cell_iterator &cell, - std::vector > &a) const +get_vertices +(const typename Triangulation::cell_iterator &cell) const { + std_cxx11::array, GeometryInfo::vertices_per_cell> vertices; + // The assertions can not be in the constructor, since this would + // require to call dof_handler.distribute_dofs(fe) *before* the mapping + // object is constructed, which is not necessarily what we want. - // The assertions can not be in the - // constructor, since this would - // require to call - // dof_handler.distribute_dofs(fe) - // *before* the mapping object is - // constructed, which is not - // necessarily what we want. - -//TODO: Only one of these two assertions should be relevant + //TODO: Only one of these two assertions should be relevant AssertDimension (spacedim, shiftmap_dof_handler->get_fe().n_dofs_per_vertex()); - AssertDimension(shiftmap_dof_handler->get_fe().n_components(), spacedim); + AssertDimension (shiftmap_dof_handler->get_fe().n_components(), spacedim); AssertDimension (shiftmap_dof_handler->n_dofs(), euler_transform_vectors->size()); - // cast the - // Triangulation::cell_iterator - // into a - // DoFHandler::cell_iterator - // which is necessary for access to + // cast the Triangulation::cell_iterator into a + // DoFHandler::cell_iterator which is necessary for access to // DoFCellAccessor::get_dof_values() typename DoFHandler::cell_iterator dof_cell (*cell, shiftmap_dof_handler); - // We require the cell to be active - // since we can only then get nodal + // We require the cell to be active since we can only then get nodal // values for the shifts Assert (dof_cell->active() == true, ExcInactiveCell()); - // for Q1 elements, the number of - // support points should equal the - // number of vertices - a.resize(GeometryInfo::vertices_per_cell); - - // now get the values of the shift - // vectors at the vertices + // now get the values of the shift vectors at the vertices Vector mapping_values (shiftmap_dof_handler->get_fe().dofs_per_cell); dof_cell->get_dof_values (*euler_transform_vectors, mapping_values); - for (unsigned int i=0; i::vertices_per_cell; ++i) { Point shift_vector; - // pick out the value of the - // shift vector at the present - // vertex. since vertex dofs - // are always numbered first, - // we can access them easily + // pick out the value of the shift vector at the present + // vertex. since vertex dofs are always numbered first, we can + // access them easily for (unsigned int j=0; jvertex(i) + shift_vector; + // compute new support point by old (reference) value and added + // shift + vertices[i] = cell->vertex(i) + shift_vector; } + return vertices; }