From: Martin Kronbichler Date: Wed, 12 Dec 2012 13:38:38 +0000 (+0000) Subject: Avoid reinterpret_cast by making the declaration of struct ShapeFunctionData in FEVal... X-Git-Tag: v8.0.0~1710 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cd5f7d2544f56bc85ca340b9652d7bf2ad00ed09;p=dealii.git Avoid reinterpret_cast by making the declaration of struct ShapeFunctionData in FEValuesViews::Scalar/Vector public. git-svn-id: https://svn.dealii.org/trunk@27796 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/fe/fe_values.h b/deal.II/include/deal.II/fe/fe_values.h index 9cb34a481f..62d3d24295 100644 --- a/deal.II/include/deal.II/fe/fe_values.h +++ b/deal.II/include/deal.II/fe/fe_values.h @@ -179,6 +179,48 @@ namespace FEValuesViews */ typedef Tensor<2,spacedim> hessian_type; + /** + * A structure where for each shape + * function we pre-compute a bunch of + * data that will make later accesses + * much cheaper. + */ + struct ShapeFunctionData + { + /** + * For each shape function, store + * whether the selected vector + * component may be nonzero. For + * primitive shape functions we + * know for sure whether a certain + * scalar component of a given + * shape function is nonzero, + * whereas for non-primitive shape + * functions this may not be + * entirely clear (e.g. for RT + * elements it depends on the shape + * of a cell). + */ + bool is_nonzero_shape_function_component; + + /** + * For each shape function, store + * the row index within the + * shape_values, shape_gradients, + * and shape_hessians tables (the + * column index is the quadrature + * point index). If the shape + * function is primitive, then we + * can get this information from + * the shape_function_to_row_table + * of the FEValues object; + * otherwise, we have to work a bit + * harder to compute this + * information. + */ + unsigned int row_index; + }; + /** * Default constructor. Creates an * invalid object. @@ -341,7 +383,6 @@ namespace FEValuesViews void get_function_laplacians (const InputVector &fe_function, std::vector &laplacians) const; - private: /** * A reference to the FEValuesBase object @@ -356,48 +397,6 @@ namespace FEValuesViews */ const unsigned int component; - /** - * A structure where for each shape - * function we pre-compute a bunch of - * data that will make later accesses - * much cheaper. - */ - struct ShapeFunctionData - { - /** - * For each shape function, store - * whether the selected vector - * component may be nonzero. For - * primitive shape functions we - * know for sure whether a certain - * scalar component of a given - * shape function is nonzero, - * whereas for non-primitive shape - * functions this may not be - * entirely clear (e.g. for RT - * elements it depends on the shape - * of a cell). - */ - bool is_nonzero_shape_function_component; - - /** - * For each shape function, store - * the row index within the - * shape_values, shape_gradients, - * and shape_hessians tables (the - * column index is the quadrature - * point index). If the shape - * function is primitive, then we - * can get this information from - * the shape_function_to_row_table - * of the FEValues object; - * otherwise, we have to work a bit - * harder to compute this - * information. - */ - unsigned int row_index; - }; - /** * Store the data about shape * functions. @@ -516,6 +515,70 @@ namespace FEValuesViews */ typedef Tensor<3,spacedim> hessian_type; + /** + * A structure where for each shape + * function we pre-compute a bunch of + * data that will make later accesses + * much cheaper. + */ + struct ShapeFunctionData + { + /** + * For each pair (shape + * function,component within + * vector), store whether the + * selected vector component may be + * nonzero. For primitive shape + * functions we know for sure + * whether a certain scalar + * component of a given shape + * function is nonzero, whereas for + * non-primitive shape functions + * this may not be entirely clear + * (e.g. for RT elements it depends + * on the shape of a cell). + */ + bool is_nonzero_shape_function_component[spacedim]; + + /** + * For each pair (shape function, + * component within vector), store + * the row index within the + * shape_values, shape_gradients, + * and shape_hessians tables (the + * column index is the quadrature + * point index). If the shape + * function is primitive, then we + * can get this information from + * the shape_function_to_row_table + * of the FEValues object; + * otherwise, we have to work a bit + * harder to compute this + * information. + */ + unsigned int row_index[spacedim]; + + /** + * For each shape function say the + * following: if only a single + * entry in + * is_nonzero_shape_function_component + * for this shape function is + * nonzero, then store the + * corresponding value of row_index + * and + * single_nonzero_component_index + * represents the index between 0 + * and dim for which it is + * attained. If multiple components + * are nonzero, then store -1. If + * no components are nonzero then + * store -2. + */ + int single_nonzero_component; + unsigned int single_nonzero_component_index; + }; + /** * Default constructor. Creates an * invalid object. @@ -849,70 +912,6 @@ namespace FEValuesViews */ const unsigned int first_vector_component; - /** - * A structure where for each shape - * function we pre-compute a bunch of - * data that will make later accesses - * much cheaper. - */ - struct ShapeFunctionData - { - /** - * For each pair (shape - * function,component within - * vector), store whether the - * selected vector component may be - * nonzero. For primitive shape - * functions we know for sure - * whether a certain scalar - * component of a given shape - * function is nonzero, whereas for - * non-primitive shape functions - * this may not be entirely clear - * (e.g. for RT elements it depends - * on the shape of a cell). - */ - bool is_nonzero_shape_function_component[spacedim]; - - /** - * For each pair (shape function, - * component within vector), store - * the row index within the - * shape_values, shape_gradients, - * and shape_hessians tables (the - * column index is the quadrature - * point index). If the shape - * function is primitive, then we - * can get this information from - * the shape_function_to_row_table - * of the FEValues object; - * otherwise, we have to work a bit - * harder to compute this - * information. - */ - unsigned int row_index[spacedim]; - - /** - * For each shape function say the - * following: if only a single - * entry in - * is_nonzero_shape_function_component - * for this shape function is - * nonzero, then store the - * corresponding value of row_index - * and - * single_nonzero_component_index - * represents the index between 0 - * and dim for which it is - * attained. If multiple components - * are nonzero, then store -1. If - * no components are nonzero then - * store -2. - */ - int single_nonzero_component; - unsigned int single_nonzero_component_index; - }; - /** * Store the data about shape * functions. @@ -989,6 +988,70 @@ namespace FEValuesViews */ typedef Tensor<1, spacedim> divergence_type; + /** + * A structure where for each shape + * function we pre-compute a bunch of + * data that will make later accesses + * much cheaper. + */ + struct ShapeFunctionData + { + /** + * For each pair (shape + * function,component within + * vector), store whether the + * selected vector component may be + * nonzero. For primitive shape + * functions we know for sure + * whether a certain scalar + * component of a given shape + * function is nonzero, whereas for + * non-primitive shape functions + * this may not be entirely clear + * (e.g. for RT elements it depends + * on the shape of a cell). + */ + bool is_nonzero_shape_function_component[value_type::n_independent_components]; + + /** + * For each pair (shape function, + * component within vector), store + * the row index within the + * shape_values, shape_gradients, + * and shape_hessians tables (the + * column index is the quadrature + * point index). If the shape + * function is primitive, then we + * can get this information from + * the shape_function_to_row_table + * of the FEValues object; + * otherwise, we have to work a bit + * harder to compute this + * information. + */ + unsigned int row_index[value_type::n_independent_components]; + + /** + * For each shape function say the + * following: if only a single + * entry in + * is_nonzero_shape_function_component + * for this shape function is + * nonzero, then store the + * corresponding value of row_index + * and + * single_nonzero_component_index + * represents the index between 0 + * and (dim^2 + dim)/2 for which it is + * attained. If multiple components + * are nonzero, then store -1. If + * no components are nonzero then + * store -2. + */ + int single_nonzero_component; + unsigned int single_nonzero_component_index; + }; + /** * Default constructor. Creates an * invalid object. @@ -1130,70 +1193,6 @@ namespace FEValuesViews */ const unsigned int first_tensor_component; - /** - * A structure where for each shape - * function we pre-compute a bunch of - * data that will make later accesses - * much cheaper. - */ - struct ShapeFunctionData - { - /** - * For each pair (shape - * function,component within - * vector), store whether the - * selected vector component may be - * nonzero. For primitive shape - * functions we know for sure - * whether a certain scalar - * component of a given shape - * function is nonzero, whereas for - * non-primitive shape functions - * this may not be entirely clear - * (e.g. for RT elements it depends - * on the shape of a cell). - */ - bool is_nonzero_shape_function_component[value_type::n_independent_components]; - - /** - * For each pair (shape function, - * component within vector), store - * the row index within the - * shape_values, shape_gradients, - * and shape_hessians tables (the - * column index is the quadrature - * point index). If the shape - * function is primitive, then we - * can get this information from - * the shape_function_to_row_table - * of the FEValues object; - * otherwise, we have to work a bit - * harder to compute this - * information. - */ - unsigned int row_index[value_type::n_independent_components]; - - /** - * For each shape function say the - * following: if only a single - * entry in - * is_nonzero_shape_function_component - * for this shape function is - * nonzero, then store the - * corresponding value of row_index - * and - * single_nonzero_component_index - * represents the index between 0 - * and (dim^2 + dim)/2 for which it is - * attained. If multiple components - * are nonzero, then store -1. If - * no components are nonzero then - * store -2. - */ - int single_nonzero_component; - unsigned int single_nonzero_component_index; - }; - /** * Store the data about shape * functions. diff --git a/deal.II/source/fe/fe_values.cc b/deal.II/source/fe/fe_values.cc index 9e4305b0a7..5a63aa75f2 100644 --- a/deal.II/source/fe/fe_values.cc +++ b/deal.II/source/fe/fe_values.cc @@ -355,17 +355,11 @@ namespace FEValuesViews // consistently throughout the code since revision 17903 at least. // ------------------------- scalar functions -------------------------- - - struct ShapeFunctionDataScalar - { - bool is_nonzero_shape_function_component; - unsigned int row_index; - }; - + template void do_function_values (const ::dealii::Vector &dof_values, const Table<2,double> &shape_values, - const std::vector &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector &values) { const unsigned int dofs_per_cell = dof_values.size(); @@ -394,11 +388,11 @@ namespace FEValuesViews // same code for gradient and Hessian, template argument 'order' to give // the order of the derivative (= rank of gradient/Hessian tensor) - template + template void do_function_derivatives (const ::dealii::Vector &dof_values, const std::vector > > &shape_derivatives, - const std::vector &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &derivatives) { const unsigned int dofs_per_cell = dof_values.size(); @@ -426,11 +420,11 @@ namespace FEValuesViews - template + template void do_function_laplacians (const ::dealii::Vector &dof_values, const std::vector > > &shape_hessians, - const std::vector &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector &laplacians) { const unsigned int dofs_per_cell = dof_values.size(); @@ -459,21 +453,10 @@ namespace FEValuesViews // ----------------------------- vector part --------------------------- - template - struct ShapeFunctionDataVector - { - bool is_nonzero_shape_function_component[spacedim]; - unsigned int row_index[spacedim]; - int single_nonzero_component; - unsigned int single_nonzero_component_index; - }; - - - - template + template void do_function_values (const ::dealii::Vector &dof_values, const Table<2,double> &shape_values, - const std::vector > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &values) { const unsigned int dofs_per_cell = dof_values.size(); @@ -518,11 +501,11 @@ namespace FEValuesViews - template + template void do_function_derivatives (const ::dealii::Vector &dof_values, const std::vector > > &shape_derivatives, - const std::vector > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &derivatives) { const unsigned int dofs_per_cell = dof_values.size(); @@ -570,11 +553,11 @@ namespace FEValuesViews - template + template void do_function_symmetric_gradients (const ::dealii::Vector &dof_values, const std::vector > > &shape_gradients, - const std::vector > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &symmetric_gradients) { const unsigned int dofs_per_cell = dof_values.size(); @@ -623,11 +606,11 @@ namespace FEValuesViews - template + template void do_function_divergences (const ::dealii::Vector &dof_values, const std::vector > > &shape_gradients, - const std::vector > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector &divergences) { const unsigned int dofs_per_cell = dof_values.size(); @@ -673,11 +656,11 @@ namespace FEValuesViews - template + template void do_function_curls (const ::dealii::Vector &dof_values, const std::vector > > &shape_gradients, - const std::vector > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector::type> &curls) { const unsigned int dofs_per_cell = dof_values.size(); @@ -858,11 +841,11 @@ namespace FEValuesViews - template + template void do_function_laplacians (const ::dealii::Vector &dof_values, const std::vector > > &shape_hessians, - const std::vector > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &laplacians) { const unsigned int dofs_per_cell = dof_values.size(); @@ -912,11 +895,11 @@ namespace FEValuesViews // ---------------------- symmetric tensor part ------------------------ - template + template void do_function_values (const ::dealii::Vector &dof_values, const Table<2,double> &shape_values, - const std::vector::n_independent_components> > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &values) { const unsigned int dofs_per_cell = dof_values.size(); @@ -966,11 +949,11 @@ namespace FEValuesViews - template + template void do_function_divergences (const ::dealii::Vector &dof_values, const std::vector > > &shape_gradients, - const std::vector::n_independent_components> > &shape_function_data, + const std::vector::ShapeFunctionData> &shape_function_data, std::vector > &divergences) { const unsigned int dofs_per_cell = dof_values.size(); @@ -1050,6 +1033,7 @@ namespace FEValuesViews } } } + } // end of namespace internal @@ -1072,9 +1056,8 @@ namespace FEValuesViews // get function values of dofs on this cell and call internal worker function dealii::Vector dof_values(fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_values (dof_values, fe_values.shape_values, - reinterpret_cast&>(shape_function_data), - values); + internal::do_function_values + (dof_values, fe_values.shape_values, shape_function_data, values); } @@ -1097,9 +1080,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_derivatives (dof_values, fe_values.shape_gradients, - reinterpret_cast&>(shape_function_data), - gradients); + internal::do_function_derivatives<1,dim,spacedim> + (dof_values, fe_values.shape_gradients, shape_function_data, gradients); } @@ -1122,9 +1104,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_derivatives (dof_values, fe_values.shape_hessians, - reinterpret_cast&>(shape_function_data), - hessians); + internal::do_function_derivatives<2,dim,spacedim> + (dof_values, fe_values.shape_hessians, shape_function_data, hessians); } @@ -1147,9 +1128,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_laplacians (dof_values, fe_values.shape_hessians, - reinterpret_cast&>(shape_function_data), - laplacians); + internal::do_function_laplacians + (dof_values, fe_values.shape_hessians, shape_function_data, laplacians); } @@ -1172,9 +1152,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_values (dof_values, fe_values.shape_values, - reinterpret_cast >&>(shape_function_data), - values); + internal::do_function_values + (dof_values, fe_values.shape_values, shape_function_data, values); } @@ -1198,9 +1177,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_derivatives (dof_values, fe_values.shape_gradients, - reinterpret_cast >&>(shape_function_data), - gradients); + internal::do_function_derivatives<1,dim,spacedim> + (dof_values, fe_values.shape_gradients, shape_function_data, gradients); } @@ -1223,10 +1201,9 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_symmetric_gradients (dof_values, - fe_values.shape_gradients, - reinterpret_cast >&>(shape_function_data), - symmetric_gradients); + internal::do_function_symmetric_gradients + (dof_values, fe_values.shape_gradients, shape_function_data, + symmetric_gradients); } @@ -1250,10 +1227,8 @@ namespace FEValuesViews // on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_divergences (dof_values, - fe_values.shape_gradients, - reinterpret_cast >&>(shape_function_data), - divergences); + internal::do_function_divergences + (dof_values, fe_values.shape_gradients, shape_function_data, divergences); } template @@ -1275,9 +1250,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values (fe_function, dof_values); - internal::do_function_curls (dof_values, fe_values.shape_gradients, - reinterpret_cast >&>(shape_function_data), - curls); + internal::do_function_curls + (dof_values, fe_values.shape_gradients, shape_function_data, curls); } @@ -1299,9 +1273,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_derivatives (dof_values, fe_values.shape_hessians, - reinterpret_cast >&>(shape_function_data), - hessians); + internal::do_function_derivatives<2,dim,spacedim> + (dof_values, fe_values.shape_hessians, shape_function_data, hessians); } @@ -1327,9 +1300,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values (fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_laplacians (dof_values, fe_values.shape_hessians, - reinterpret_cast >&>(shape_function_data), - laplacians); + internal::do_function_laplacians + (dof_values, fe_values.shape_hessians, shape_function_data, laplacians); } @@ -1352,9 +1324,8 @@ namespace FEValuesViews // get function values of dofs on this cell dealii::Vector dof_values(fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_values (dof_values, fe_values.shape_values, - reinterpret_cast::n_independent_components> >&>(shape_function_data), - values); + internal::do_function_values + (dof_values, fe_values.shape_values, shape_function_data, values); } @@ -1378,9 +1349,8 @@ namespace FEValuesViews // on this cell dealii::Vector dof_values(fe_values.dofs_per_cell); fe_values.present_cell->get_interpolated_dof_values(fe_function, dof_values); - internal::do_function_divergences (dof_values, fe_values.shape_gradients, - reinterpret_cast::n_independent_components> >&>(shape_function_data), - divergences); + internal::do_function_divergences + (dof_values, fe_values.shape_gradients, shape_function_data, divergences); } }