From: bangerth Date: Thu, 11 Dec 2008 21:34:23 +0000 (+0000) Subject: Change the way we store data to keep thing accessed at the same time together in... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cb800e3d63cb556070ec0b7740d2a4def43f2e5f;p=dealii-svn.git Change the way we store data to keep thing accessed at the same time together in memory. Also reduce the number of indirections through this scheme. git-svn-id: https://svn.dealii.org/trunk@17918 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index fa4a4b5eae..4fc65ba2c1 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -257,34 +257,52 @@ namespace FEValuesViews const unsigned int component; /** - * 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). + * A structure where for each shape + * function we pre-compute a bunch of + * data that will make later accesses + * much cheaper. */ - Table<1,bool> is_nonzero_shape_function_component; + 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; + }; /** - * 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. + * Store the data about shape + * functions. */ - Table<1,unsigned int> row_index; + std::vector shape_function_data; }; @@ -467,53 +485,74 @@ namespace FEValuesViews const unsigned int first_vector_component; /** - * 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). - */ - Table<2,bool> is_nonzero_shape_function_component; - - /** - * 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. + * A structure where for each shape + * function we pre-compute a bunch of + * data that will make later accesses + * much cheaper. */ - Table<2,unsigned int> row_index; + 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[dim]; + + /** + * 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[dim]; + + /** + * 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; + }; /** - * 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. + * Store the data about shape + * functions. */ - Table<1,int> single_nonzero_component; - Table<1,unsigned int> single_nonzero_component_index; + std::vector shape_function_data; }; } @@ -3418,8 +3457,10 @@ namespace FEValuesViews // component as fixed and we have // pre-computed and cached a bunch of // information. see the comments there - if (is_nonzero_shape_function_component[shape_function]) - return fe_values.shape_values(row_index[shape_function],q_point); + if (shape_function_data[shape_function].is_nonzero_shape_function_component) + return fe_values.shape_values(shape_function_data[shape_function] + .row_index, + q_point); else return 0; } @@ -3445,8 +3486,9 @@ namespace FEValuesViews // component as fixed and we have // pre-computed and cached a bunch of // information. see the comments there - if (is_nonzero_shape_function_component[shape_function]) - return fe_values.shape_gradients[row_index[shape_function]][q_point]; + if (shape_function_data[shape_function].is_nonzero_shape_function_component) + return fe_values.shape_gradients[shape_function_data[shape_function] + .row_index][q_point]; else return gradient_type(); } @@ -3471,8 +3513,8 @@ namespace FEValuesViews // component as fixed and we have // pre-computed and cached a bunch of // information. see the comments there - if (is_nonzero_shape_function_component[shape_function]) - return fe_values.shape_hessians[row_index[shape_function]][q_point]; + if (shape_function_data[shape_function].is_nonzero_shape_function_component) + return fe_values.shape_hessians[shape_function_data[shape_function].row_index][q_point]; else return hessian_type(); } @@ -3493,13 +3535,13 @@ namespace FEValuesViews // same as for the scalar case except // that we have one more index - const int snc = single_nonzero_component(shape_function); + const int snc = shape_function_data[shape_function].single_nonzero_component; if (snc == -2) return value_type(); else if (snc != -1) { value_type return_value; - return_value[single_nonzero_component_index[shape_function]] + return_value[shape_function_data[shape_function].single_nonzero_component_index] = fe_values.shape_values(snc,q_point); return return_value; } @@ -3507,9 +3549,9 @@ namespace FEValuesViews { value_type return_value; for (unsigned int d=0; ddofs_per_cell), - row_index (fe_values.fe->dofs_per_cell) + shape_function_data (fe_values.fe->dofs_per_cell) { Assert (component < fe_values.fe->n_components(), ExcIndexRange(component, 0, fe_values.fe->n_components())); @@ -76,28 +75,29 @@ namespace FEValuesViews fe_values.fe->is_primitive(i)); if (is_primitive == true) - is_nonzero_shape_function_component[i] + shape_function_data[i].is_nonzero_shape_function_component = (component == fe_values.fe->system_to_component_index(i).first); else - is_nonzero_shape_function_component[i] + shape_function_data[i].is_nonzero_shape_function_component = (fe_values.fe->get_nonzero_components(i)[component] == true); - if (is_nonzero_shape_function_component[i] == true) + if (shape_function_data[i].is_nonzero_shape_function_component == true) { if (is_primitive == true) - row_index[i] = shape_function_to_row_table[i]; + shape_function_data[i].row_index = shape_function_to_row_table[i]; else - row_index[i] = (shape_function_to_row_table[i] - + - std::count (fe_values.fe->get_nonzero_components(i).begin(), - fe_values.fe->get_nonzero_components(i).begin()+ - component, - true)); + shape_function_data[i].row_index + = (shape_function_to_row_table[i] + + + std::count (fe_values.fe->get_nonzero_components(i).begin(), + fe_values.fe->get_nonzero_components(i).begin()+ + component, + true)); } else - row_index[i] = numbers::invalid_unsigned_int; + shape_function_data[i].row_index = numbers::invalid_unsigned_int; } } @@ -128,12 +128,7 @@ namespace FEValuesViews : fe_values (fe_values), first_vector_component (first_vector_component), - is_nonzero_shape_function_component (fe_values.fe->dofs_per_cell, - dim), - row_index (fe_values.fe->dofs_per_cell, - dim), - single_nonzero_component (fe_values.fe->dofs_per_cell), - single_nonzero_component_index (fe_values.fe->dofs_per_cell) + shape_function_data (fe_values.fe->dofs_per_cell) { Assert (first_vector_component+dim-1 < fe_values.fe->n_components(), ExcIndexRange(first_vector_component+dim-1, 0, @@ -152,28 +147,32 @@ namespace FEValuesViews fe_values.fe->is_primitive(i)); if (is_primitive == true) - is_nonzero_shape_function_component[i][d] + shape_function_data[i].is_nonzero_shape_function_component[d] = (component == fe_values.fe->system_to_component_index(i).first); else - is_nonzero_shape_function_component[i][d] + shape_function_data[i].is_nonzero_shape_function_component[d] = (fe_values.fe->get_nonzero_components(i)[component] == true); - if (is_nonzero_shape_function_component[i][d] == true) + if (shape_function_data[i].is_nonzero_shape_function_component[d] + == true) { if (is_primitive == true) - row_index[i][d] = shape_function_to_row_table[i]; + shape_function_data[i].row_index[d] + = shape_function_to_row_table[i]; else - row_index[i][d] = (shape_function_to_row_table[i] - + - std::count (fe_values.fe->get_nonzero_components(i).begin(), - fe_values.fe->get_nonzero_components(i).begin()+ - component, - true)); + shape_function_data[i].row_index[d] + = (shape_function_to_row_table[i] + + + std::count (fe_values.fe->get_nonzero_components(i).begin(), + fe_values.fe->get_nonzero_components(i).begin()+ + component, + true)); } else - row_index[i][d] = numbers::invalid_unsigned_int; + shape_function_data[i].row_index[d] + = numbers::invalid_unsigned_int; } } @@ -181,20 +180,24 @@ namespace FEValuesViews { unsigned int n_nonzero_components = 0; for (unsigned int d=0; d 1) - single_nonzero_component(i) = -1; + shape_function_data[i].single_nonzero_component = -1; else { for (unsigned int d=0; d