/* @{ */
void initialize_update_flags();
+ /**
+ * Add FEValues UpdateFlags for
+ * integration on all objects
+ * (cells, boundary faces and
+ * all interior faces).
+ */
+ void add_update_flags_all (const UpdateFlags flags);
+
+ /**
+ * Add FEValues UpdateFlags for
+ * integration on cells.
+ */
+ void add_update_flags_cell(const UpdateFlags flags);
+
+ /**
+ * Add FEValues UpdateFlags for
+ * integration on boundary faces.
+ */
+ void add_update_flags_boundary(const UpdateFlags flags);
+
+ /**
+ * Add FEValues UpdateFlags for
+ * integration on interior faces.
+ */
+ void add_update_flags_face(const UpdateFlags flags);
+
/**
* Add additional update flags
* to the ones already set in
* boolean flags indicate
* wether the additional flags
* should be set for cell,
- * boundary, interelement face,
- * or neighbor integration, or
+ * boundary, interelement face
+ * for the cell itself
+ * or neighbor cell, or
* any combination thereof.
*/
void add_update_flags(const UpdateFlags flags,
{
cell_quadrature = QGauss<1>(cp);
}
+
+
+ template <int dim, int sdim>
+ inline
+ void
+ IntegrationInfoBox<dim,sdim>::add_update_flags_all (const UpdateFlags flags)
+ {
+ add_update_flags(flags, true, true, true, true);
+ }
+
+
+ template <int dim, int sdim>
+ inline
+ void
+ IntegrationInfoBox<dim,sdim>::add_update_flags_cell (const UpdateFlags flags)
+ {
+ add_update_flags(flags, true, false, false, false);
+ }
+
+
+ template <int dim, int sdim>
+ inline
+ void
+ IntegrationInfoBox<dim,sdim>::add_update_flags_boundary (const UpdateFlags flags)
+ {
+ add_update_flags(flags, false, true, false, false);
+ }
+
+
+ template <int dim, int sdim>
+ inline
+ void
+ IntegrationInfoBox<dim,sdim>::add_update_flags_face (const UpdateFlags flags)
+ {
+ add_update_flags(flags, false, false, true, true);
+ }
+
}
DEAL_II_NAMESPACE_CLOSE
const unsigned int n_gauss_points = dof_handler.get_fe().tensor_degree()+1;
info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points);
// which values to update in these
- // points. We call
- // <tt>initialize_update_flags</tt>
- // first in order to set default
- // values. Then, we add what we
- // need additionally.
- info_box.initialize_update_flags();
+ // points. Update flags are
+ // initialized to some default
+ // values to be able to
+ // integrate. Here, we add what we
+ // need additionally, namely the
+ // values and gradients of shape
+ // functions on all objects (cells,
+ // boundary and interior faces).
UpdateFlags update_flags = update_values | update_gradients;
- info_box.add_update_flags(update_flags, true, true, true, true);
+ info_box.add_update_flags_all(update_flags);
info_box.initialize(fe, mapping);
// This is the object into which we
MeshWorker::IntegrationInfoBox<dim> info_box;
const unsigned int n_gauss_points = mg_dof_handler.get_fe().tensor_degree()+1;
info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points);
- info_box.initialize_update_flags();
UpdateFlags update_flags = update_values | update_gradients;
- info_box.add_update_flags(update_flags, true, true, true, true);
+ info_box.add_update_flags_all(update_flags);
info_box.initialize(fe, mapping);
MeshWorker::DoFInfo<dim> dof_info(mg_dof_handler);
MeshWorker::IntegrationInfoBox<dim> info_box;
const unsigned int n_gauss_points = dof_handler.get_fe().tensor_degree()+1;
info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points+1, n_gauss_points);
- info_box.initialize_update_flags();
UpdateFlags update_flags = update_quadrature_points | update_values | update_gradients;
- info_box.add_update_flags(update_flags, true, true, true, true);
+ info_box.add_update_flags_all(update_flags);
info_box.initialize(fe, mapping);
MeshWorker::DoFInfo<dim> dof_info(dof_handler);
// Then, we tell the Meshworker::VectorSelector
// for cells, that we need the
// second derivatives of this
- // solution (to compute the Laplacian).
+ // solution (to compute the
+ // Laplacian). Therefore, the
+ // boolean arguments selecting
+ // function values and first
+ // derivatives a false, only the
+ // last one selecting second
+ // derivatives is true.
info_box.cell_selector.add("solution", false, false, true);
// On interior and boundary faces,
// we need the function values and
- // the first derivatives.
+ // the first derivatives, but not
+ // second derivatives.
info_box.boundary_selector.add("solution", true, true, false);
info_box.face_selector.add("solution", true, true, false);
// update flags are already
// adjusted to the values and
// derivatives we requested above.
- info_box.initialize_update_flags();
- info_box.add_update_flags(update_quadrature_points, false, true, false, false);
+ info_box.add_update_flags_boundary(update_quadrature_points);
info_box.initialize(fe, mapping, solution_data);
MeshWorker::DoFInfo<dim> dof_info(dof_handler);