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<ShapeFunctionData> shape_function_data;
};
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<ShapeFunctionData> shape_function_data;
};
}
// 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;
}
// 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();
}
// 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();
}
// 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;
}
{
value_type return_value;
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(shape_function,d))
+ if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
return_value[d]
- = fe_values.shape_values(row_index(shape_function,d),q_point);
+ = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
return return_value;
}
// 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 gradient_type();
else if (snc != -1)
{
gradient_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_gradients[snc][q_point];
return return_value;
}
{
gradient_type return_value;
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(shape_function,d))
+ if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
return_value[d]
- = fe_values.shape_gradients[row_index(shape_function,d)][q_point];
+ = fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
return return_value;
}
// 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 divergence_type();
else if (snc != -1)
return
- fe_values.shape_gradients[snc][q_point][single_nonzero_component_index[shape_function]];
+ fe_values.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
else
{
divergence_type return_value = 0;
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(shape_function,d))
+ if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
return_value
- += fe_values.shape_gradients[row_index(shape_function,d)][q_point][d];
+ += fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
return return_value;
}
// 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 hessian_type();
else if (snc != -1)
{
hessian_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_hessians[snc][q_point];
return return_value;
}
{
hessian_type return_value;
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(shape_function,d))
+ if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
return_value[d]
- = fe_values.shape_hessians[row_index(shape_function,d)][q_point];
+ = fe_values.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
return return_value;
}
// 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 symmetric_gradient_type();
else if (snc != -1)
{
gradient_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_gradients[snc][q_point];
return symmetrize(return_value);
}
{
gradient_type return_value;
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(shape_function,d))
+ if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
return_value[d]
- = fe_values.shape_gradients[row_index(shape_function,d)][q_point];
+ = fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
return symmetrize(return_value);
}
:
fe_values (fe_values),
component (component),
- is_nonzero_shape_function_component (fe_values.fe->dofs_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()));
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;
}
}
:
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,
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;
}
}
{
unsigned int n_nonzero_components = 0;
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(i,d) == true)
+ if (shape_function_data[i].is_nonzero_shape_function_component[d]
+ == true)
++n_nonzero_components;
if (n_nonzero_components == 0)
- single_nonzero_component(i) = -2;
+ shape_function_data[i].single_nonzero_component = -2;
else if (n_nonzero_components > 1)
- single_nonzero_component(i) = -1;
+ shape_function_data[i].single_nonzero_component = -1;
else
{
for (unsigned int d=0; d<dim; ++d)
- if (is_nonzero_shape_function_component(i,d) == true)
+ if (shape_function_data[i].is_nonzero_shape_function_component[d]
+ == true)
{
- single_nonzero_component(i) = row_index(i,d);
- single_nonzero_component_index(i) = d;
+ shape_function_data[i].single_nonzero_component
+ = shape_function_data[i].row_index[d];
+ shape_function_data[i].single_nonzero_component_index
+ = d;
break;
}
}