From 98c88143f8cc41f58369020c48382f13e6387b59 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Sun, 26 Jul 2015 12:04:43 +0200 Subject: [PATCH] Moved dof values and indices inside InternalData. Closes #941. --- include/deal.II/fe/mapping_fe_field.h | 29 +++++++++--------- source/fe/mapping_fe_field.cc | 42 ++++++++++++--------------- 2 files changed, 32 insertions(+), 39 deletions(-) diff --git a/include/deal.II/fe/mapping_fe_field.h b/include/deal.II/fe/mapping_fe_field.h index 23fbcdf404..31438418dc 100644 --- a/include/deal.II/fe/mapping_fe_field.h +++ b/include/deal.II/fe/mapping_fe_field.h @@ -329,6 +329,18 @@ public: * used for the euler vector and the euler dof handler. */ ComponentMask mask; + + + /** + * Storage for the indices of the local degrees of freedom. + */ + std::vector local_dof_indices; + + + /** + * Storage for local degrees of freedom. + */ + std::vector local_dof_values; }; @@ -508,21 +520,8 @@ private: /** * Update internal degrees of freedom. */ - void update_internal_dofs(const typename Triangulation::cell_iterator &cell) const; - - /** - * It stores the local degrees of freedom of the DH for each cell (i.e. - * euler_vector * dof_indices, see method update_internal_dofs for more - * clarifications.). - */ - mutable Threads::ThreadLocalStorage > local_dof_values; - - /** - * Store the degrees of freedom of the DH for each cell (i.e. - * cell->get_dof_indices(dof_indices), see method update_internal_dofs for - * more clarifications.). Thread safe. - */ - mutable Threads::ThreadLocalStorage > local_dof_indices; + void update_internal_dofs(const typename Triangulation::cell_iterator &cell, + typename MappingFEField::InternalData &data) const; /** * Reimplemented from Mapping. See the documentation of the base class for diff --git a/source/fe/mapping_fe_field.cc b/source/fe/mapping_fe_field.cc index a218480ab7..90fcbc780c 100644 --- a/source/fe/mapping_fe_field.cc +++ b/source/fe/mapping_fe_field.cc @@ -49,7 +49,9 @@ MappingFEField::InternalData::InternalData (const Finite const ComponentMask mask) : n_shape_functions (fe.dofs_per_cell), - mask (mask) + mask (mask), + local_dof_indices(fe.dofs_per_cell), + local_dof_values(fe.dofs_per_cell) {} @@ -70,8 +72,6 @@ MappingFEField::MappingFEField (const DH &euler_dof euler_vector(&euler_vector), fe(&euler_dof_handler.get_fe()), euler_dof_handler(&euler_dof_handler), - local_dof_values(std::vector(fe->dofs_per_cell)), - local_dof_indices(std::vector(fe->dofs_per_cell)), fe_mask(mask.size() ? mask : ComponentMask(fe->get_nonzero_components(0).size(), true)), fe_to_real(fe_mask.size(), numbers::invalid_unsigned_int) @@ -92,8 +92,6 @@ MappingFEField::MappingFEField (const MappingFEField(fe->dofs_per_cell)), - local_dof_indices(std::vector(fe->dofs_per_cell)), fe_mask(mapping.fe_mask), fe_to_real(mapping.fe_to_real) {} @@ -386,9 +384,6 @@ MappingFEField::fill_fe_values ( std::vector > &normal_vectors, CellSimilarity::Similarity &cell_similarity) const { - AssertDimension(fe->dofs_per_cell, local_dof_values.get().size()); - Assert(local_dof_values.get().size()>0, ExcMessage("Cannot do anything with zero degrees of freedom")); - // convert data object to internal data for this class. fails with an // exception if that is not possible Assert (dynamic_cast (&mapping_data) != 0, ExcInternalError()); @@ -520,7 +515,7 @@ MappingFEField::fill_fe_values ( for (unsigned int j=0; j::cell_it // Use the get_data function to create an InternalData with data vectors of // the right size and transformation shape values already computed at point // p. - update_internal_dofs(cell); - const Quadrature point_quadrature(p); std_cxx11::unique_ptr mdata (dynamic_cast ( get_data(update_transformation_values, point_quadrature))); + update_internal_dofs(cell, *mdata); + return this->transform_unit_to_real_cell_internal(*mdata); } @@ -801,7 +796,7 @@ transform_unit_to_real_cell_internal (const InternalData &data) const { unsigned int comp_i = fe->system_to_component_index(i).first; if (fe_mask[comp_i]) - p_real[fe_to_real[comp_i]] += local_dof_values.get()[i] * data.shape(0,i); + p_real[fe_to_real[comp_i]] += data.local_dof_values[i] * data.shape(0,i); } return p_real; @@ -815,8 +810,6 @@ MappingFEField:: transform_real_to_unit_cell (const typename Triangulation::cell_iterator &cell, const Point &p) const { - update_internal_dofs(cell); - // first a Newton iteration based on the real mapping. It uses the center // point of the cell as a starting point Point initial_p_unit; @@ -848,6 +841,8 @@ transform_real_to_unit_cell (const typename Triangulation::cell_it mdata (dynamic_cast ( get_data(update_flags,point_quadrature))); + update_internal_dofs(cell, *mdata); + return this->transform_real_to_unit_cell_internal(cell, p, initial_p_unit, *mdata); } @@ -896,7 +891,7 @@ transform_real_to_unit_cell_internal unsigned int comp_k = fe->system_to_component_index(k).first; if (fe_mask[comp_k]) for (unsigned int j=0; j::compute_fill ( std::vector > &quadrature_points) const { const UpdateFlags update_flags(data.current_update_flags()); - { - update_internal_dofs(cell); - } + update_internal_dofs(cell, data); // first compute quadrature points if (update_flags & update_quadrature_points) @@ -990,7 +983,7 @@ MappingFEField::compute_fill ( { unsigned int comp_k = fe->system_to_component_index(k).first; if (fe_mask[comp_k]) - result[fe_to_real[comp_k]] += local_dof_values.get()[k] * shape[k]; + result[fe_to_real[comp_k]] += data.local_dof_values[k] * shape[k]; } quadrature_points[point] = result; @@ -1024,7 +1017,7 @@ MappingFEField::compute_fill ( { unsigned int comp_k = fe->system_to_component_index(k).first; if (fe_mask[comp_k]) - result[fe_to_real[comp_k]] += local_dof_values.get()[k] * data_derv[k]; + result[fe_to_real[comp_k]] += data.local_dof_values[k] * data_derv[k]; } // write result into contravariant data. for @@ -1224,17 +1217,18 @@ MappingFEField::clone () const template void MappingFEField::update_internal_dofs ( - const typename Triangulation::cell_iterator &cell) const + const typename Triangulation::cell_iterator &cell, + typename MappingFEField::InternalData &data) const { Assert(euler_dof_handler != 0, ExcMessage("euler_dof_handler is empty")); typename DH::cell_iterator dof_cell(*cell, euler_dof_handler); Assert (dof_cell->active() == true, ExcInactiveCell()); - dof_cell->get_dof_indices(local_dof_indices.get()); + dof_cell->get_dof_indices(data.local_dof_indices); - for (unsigned int i=0; i