From 69fa4b1a7ea2ed7edb370c83199372c126b973fa Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 4 Sep 2015 11:18:13 -0500 Subject: [PATCH] Some more preparatory work for separating MappingQ/Q1. Some documentation updates. Also make the Q1::InternalData object in MappingQ a pointer since that will fit better into the structure we will need later on. --- include/deal.II/fe/mapping_q.h | 10 ++--- include/deal.II/fe/mapping_q1.h | 24 ++++++----- source/fe/mapping_q.cc | 72 ++++++++++++++++++--------------- source/fe/mapping_q1.cc | 2 +- 4 files changed, 60 insertions(+), 48 deletions(-) diff --git a/include/deal.II/fe/mapping_q.h b/include/deal.II/fe/mapping_q.h index 7fb3cc0364..8697acdefd 100644 --- a/include/deal.II/fe/mapping_q.h +++ b/include/deal.II/fe/mapping_q.h @@ -47,10 +47,10 @@ template class TensorProductPolynomials; * the documentation of FiniteElement or the one of Triangulation. * * @note Since the boundary description is closely tied to the unit cell - * support points, new boundary descriptions need to explicitly use the Gauss- - * Lobatto points. + * support points, new boundary descriptions need to explicitly use the + * Gauss-Lobatto points. * - * @author Ralf Hartmann, 2000, 2001, 2005; Guido Kanschat 2000, 2001 + * @author Ralf Hartmann, 2000, 2001, 2005; Guido Kanschat 2000, 2001, Wolfgang Bangerth, 2015 */ template class MappingQ : public MappingQ1 @@ -272,10 +272,10 @@ protected: mutable bool use_mapping_q1_on_current_cell; /** - * A structure to store the corresponding information for the pure + * A pointer to a structure to store the information for the pure * $Q_1$ mapping that is, by default, used on all interior cells. */ - typename MappingQ1::InternalData mapping_q1_data; + std_cxx11::unique_ptr::InternalData> mapping_q1_data; }; protected: diff --git a/include/deal.II/fe/mapping_q1.h b/include/deal.II/fe/mapping_q1.h index bbdf4b4df8..92943025ac 100644 --- a/include/deal.II/fe/mapping_q1.h +++ b/include/deal.II/fe/mapping_q1.h @@ -37,19 +37,19 @@ DEAL_II_NAMESPACE_OPEN * Mapping of the reference to cell to a general * quadrilateral/hexahedra by $d$-linear shape functions. * - * This mapping implemented by this class maps the reference (unit) cell + * The mapping implemented by this class maps the reference (unit) cell * to a general grid cell with * straight lines in $d$ dimensions. (Note, however, that in 3D the * faces of a general, trilinearly mapped cell may be curved, even if the * edges are not). This is the standard mapping used for polyhedral domains. It - * is also the mapping used throughout deal.II for many functions that two - * variants, one that allows to pass a mapping argument explicitly and one + * is also the mapping used throughout deal.II for many functions that come in + * two variants, one that allows to pass a mapping argument explicitly and one * that simply falls back to the MappingQ1 class declared here. * * The shape functions for this mapping are the same as for the finite element FE_Q * of order 1. Therefore, coupling these two yields an isoparametric element. * - * @author Guido Kanschat, 2000, 2001; Ralf Hartmann, 2000, 2001, 2005 + * @author Guido Kanschat, 2000, 2001; Ralf Hartmann, 2000, 2001, 2005, Wolfgang Bangerth, 2015 */ template class MappingQ1 : public MappingQGeneric @@ -182,15 +182,19 @@ protected: /** * Computes the support points of the mapping. For @p MappingQ1 these are - * the vertices. However, other classes may override this function. In - * particular, the MappingQ1Eulerian class does exactly this by not - * computing the support points from the geometry of the current cell but + * the vertices, as obtained by calling Mapping::get_vertices(). + * + * By default, that function just computes the locations of the vertices as + * reported by the Triangulation. However, other classes may override this + * function. In particular, the MappingQ1Eulerian class does exactly this by + * not computing the support points from the geometry of the current cell but * instead evaluating an externally given displacement field in addition to * the geometry of the cell. */ - virtual void compute_mapping_support_points( - const typename Triangulation::cell_iterator &cell, - std::vector > &a) const; + virtual + void + compute_mapping_support_points(const typename Triangulation::cell_iterator &cell, + std::vector > &a) const; }; diff --git a/source/fe/mapping_q.cc b/source/fe/mapping_q.cc index dc81742075..207cec2ac8 100644 --- a/source/fe/mapping_q.cc +++ b/source/fe/mapping_q.cc @@ -37,8 +37,7 @@ template MappingQ::InternalData::InternalData (const unsigned int polynomial_degree) : MappingQ1::InternalData(polynomial_degree), - use_mapping_q1_on_current_cell(false), - mapping_q1_data(1) + use_mapping_q1_on_current_cell(false) {} @@ -116,14 +115,17 @@ MappingQ::get_data (const UpdateFlags update_flags, std_cxx11::cref(quadrature), quadrature.size()))); if (!use_mapping_q_on_all_cells) - tasks += Threads::new_task (std_cxx11::function - (std_cxx11::bind(&MappingQ1::InternalData::initialize, - &data->mapping_q1_data, - update_flags, - std_cxx11::cref(quadrature), - quadrature.size()))); - tasks.join_all (); + { + data->mapping_q1_data.reset (new typename MappingQ1::InternalData(1)); + tasks += Threads::new_task (std_cxx11::function + (std_cxx11::bind(&MappingQ1::InternalData::initialize, + data->mapping_q1_data.get(), + update_flags, + std_cxx11::cref(quadrature), + quadrature.size()))); + } + tasks.join_all (); return data; } @@ -146,14 +148,17 @@ MappingQ::get_face_data (const UpdateFlags update_flags, std_cxx11::cref(q), quadrature.size()))); if (!use_mapping_q_on_all_cells) - tasks += Threads::new_task (std_cxx11::function - (std_cxx11::bind(&MappingQ1::InternalData::initialize_face, - &data->mapping_q1_data, - update_flags, - std_cxx11::cref(q), - quadrature.size()))); - tasks.join_all (); + { + data->mapping_q1_data.reset (new typename MappingQ1::InternalData(1)); + tasks += Threads::new_task (std_cxx11::function + (std_cxx11::bind(&MappingQ1::InternalData::initialize_face, + data->mapping_q1_data.get(), + update_flags, + std_cxx11::cref(q), + quadrature.size()))); + } + tasks.join_all (); return data; } @@ -176,14 +181,17 @@ MappingQ::get_subface_data (const UpdateFlags update_flags, std_cxx11::cref(q), quadrature.size()))); if (!use_mapping_q_on_all_cells) - tasks += Threads::new_task (std_cxx11::function - (std_cxx11::bind(&MappingQ1::InternalData::initialize_face, - &data->mapping_q1_data, - update_flags, - std_cxx11::cref(q), - quadrature.size()))); - tasks.join_all (); + { + data->mapping_q1_data.reset (new typename MappingQ1::InternalData(1)); + tasks += Threads::new_task (std_cxx11::function + (std_cxx11::bind(&MappingQ1::InternalData::initialize_face, + data->mapping_q1_data.get(), + update_flags, + std_cxx11::cref(q), + quadrature.size()))); + } + tasks.join_all (); return data; } @@ -214,7 +222,7 @@ fill_fe_values (const typename Triangulation::cell_iterator &cell, const typename MappingQ1::InternalData * p_data = (data.use_mapping_q1_on_current_cell ? - &data.mapping_q1_data + data.mapping_q1_data.get() : &data); @@ -273,7 +281,7 @@ fill_fe_face_values (const typename Triangulation::cell_iterator & const typename MappingQ1::InternalData *p_data = (data.use_mapping_q1_on_current_cell ? - &data.mapping_q1_data + data.mapping_q1_data.get() : &data); MappingQ1::fill_fe_face_values(cell, @@ -315,7 +323,7 @@ fill_fe_subface_values (const typename Triangulation::cell_iterato const typename MappingQ1::InternalData *p_data = (data.use_mapping_q1_on_current_cell ? - &data.mapping_q1_data + data.mapping_q1_data.get() : &data); MappingQ1::fill_fe_subface_values(cell, @@ -876,7 +884,7 @@ transform (const VectorSlice > > input, { // If we only use the Q1-portion, we have to extract that data object if (data->use_mapping_q1_on_current_cell) - q1_data = &data->mapping_q1_data; + q1_data = data->mapping_q1_data.get(); } // Now, q1_data should have the right tensors in it and we call the base @@ -906,7 +914,7 @@ transform (const VectorSlice { // If we only use the Q1-portion, we have to extract that data object if (data->use_mapping_q1_on_current_cell) - q1_data = &data->mapping_q1_data; + q1_data = data->mapping_q1_data.get(); } // Now, q1_data should have the right tensors in it and we call the base @@ -934,7 +942,7 @@ transform (const VectorSlice > > input, { // If we only use the Q1-portion, we have to extract that data object if (data->use_mapping_q1_on_current_cell) - q1_data = &data->mapping_q1_data; + q1_data = data->mapping_q1_data.get(); } // Now, q1_data should have the right tensors in it and we call the base @@ -964,7 +972,7 @@ transform (const VectorSlice { // If we only use the Q1-portion, we have to extract that data object if (data->use_mapping_q1_on_current_cell) - q1_data = &data->mapping_q1_data; + q1_data = data->mapping_q1_data.get(); } // Now, q1_data should have the right tensors in it and we call the base @@ -992,7 +1000,7 @@ transform (const VectorSlice > > input, { // If we only use the Q1-portion, we have to extract that data object if (data->use_mapping_q1_on_current_cell) - q1_data = &data->mapping_q1_data; + q1_data = data->mapping_q1_data.get(); } // Now, q1_data should have the right tensors in it and we call the base @@ -1019,7 +1027,7 @@ transform_unit_to_real_cell (const typename Triangulation::cell_it typename MappingQ1::InternalData *p_data = (mdata->use_mapping_q1_on_current_cell ? - &mdata->mapping_q1_data : + mdata->mapping_q1_data.get() : &*mdata); compute_mapping_support_points(cell, p_data->mapping_support_points); diff --git a/source/fe/mapping_q1.cc b/source/fe/mapping_q1.cc index 2b522e6d0f..6bed941c46 100644 --- a/source/fe/mapping_q1.cc +++ b/source/fe/mapping_q1.cc @@ -182,8 +182,8 @@ MappingQ1::compute_mapping_support_points( { std_cxx11::array, GeometryInfo::vertices_per_cell> vertices = this->get_vertices(cell); - a.resize(GeometryInfo::vertices_per_cell); + a.resize(GeometryInfo::vertices_per_cell); for (unsigned int i=0; i::vertices_per_cell; ++i) a[i] = vertices[i]; } -- 2.39.5