From 2cad3b8509790e9d537f1ad70cc14d1452936781 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 21 Dec 2016 08:28:52 -0700 Subject: [PATCH] Make cells available in data-postprocessors. --- include/deal.II/numerics/data_postprocessor.h | 312 ++++++++++++++---- 1 file changed, 240 insertions(+), 72 deletions(-) diff --git a/include/deal.II/numerics/data_postprocessor.h b/include/deal.II/numerics/data_postprocessor.h index 4348172a75..230874d9ef 100644 --- a/include/deal.II/numerics/data_postprocessor.h +++ b/include/deal.II/numerics/data_postprocessor.h @@ -25,6 +25,8 @@ #include #include +#include + #include #include @@ -37,16 +39,190 @@ DEAL_II_NAMESPACE_OPEN */ namespace DataPostprocessorInputs { + /** + * A base class containing common elements for the Scalar and Vector classes + * that are passed as arguments to + * DataPostprocessor::evaluate_scalar_field() and + * DataPostprocessor::evaluate_vector_field(). This common base class + * provides access to the points at which the solution is being evaluated, + * and a few other fields, as described in the following. + * + *

Normal vector access

+ * + * If appropriate, i.e., if the object that is currently being processed + * is a face of a cell and the DataPostprocessor object is called from + * DataOutFaces or a similar class, then the current object also + * stores the normal vectors to the geometry on which output + * is generated, at these evaluation points. + * + * On the other hand, if the solution is being evaluated on a cell, + * then the @p normal_vectors member variable does not contain anything + * useful. + * + *

Cell access

+ * + * DataPostprocessor is typically called from classes such as DataOut + * or DataOutFaces that evaluate solution fields on a cell-by-cell + * basis. As a consequence, classes derived from DataPostprocessor + * (or DataPostprocessorScalar or DataPostprocessorVector) sometimes + * need to use which cell is currently under investigation. Consequently, + * DataOut and similar classes pass the cell they are currently working + * on to DataPostprocessor via the classes in this namespace (and + * specifically the current base class). + * + * However, the situation is not so simple. This is because the current + * class (and those derived from it) only knows the space dimension in + * which the output lives. But this can come from many sources. First, + * the cell may be a cell in a DoFHandler or hp::DoFHandler object. + * Second, if we are in 3d, this may be because we are working on a + * DoFHandler<3>, or a DoFHandler<2,3> (i.e., either a 3d mesh, or a + * 2d meshes of a 2d surface embedded in 3d space). Finally, if one + * considers classes such as DataOutRotation or DataOutStack, then + * @p spacedim being equal to 3 might mean that we are actually + * working on a DoFHandler<2> or hp::DoFHandler<2>. + * + * In other words, just because we know the value of the @p spacedim + * template argument of the current class does not mean that the + * data type of the cell iterator that is currently being worked on + * is obvious. + * + * To make the cell iterator accessible nevertheless, this class uses + * an object of type boost::any to store the cell iterator. You can + * think of this as being a void pointer that can point to anything. + * To use what is being used therefore requires the user to know the + * data type of the thing being pointed to. + * + * To make this work, the DataOut and related classes store in objects + * of the current type a representation of the cell. To get it back out, + * you would use the get_cell() function that requires you to say, + * as a template parameter, the DoFHandler type to which the cell that + * is currently being processed belongs. This is knowledge you typically + * have in an application: for example, if your application runs in + * @p dim space dimensions, uses a hp::DoFHandler, and you are currently + * using the DataOut class, then the cells that are worked on have data + * type DataOut::cell_iterator. Consequently, in a + * postprocessor, you can call + * inputs.get_cell@ @> . For technical + * reasons, however, C++ will typically require you to write this as + * inputs.template get_cell@ @> + * because the member function we call here requires that we explicitly + * provide the template argument. + * + * Let us consider a complete example of a postprocessor that computes + * the fluid norm of the stress $\|\sigma\| = \|\eta \nabla u\|$ from the + * viscosity $\eta$ and the gradient of the fluid velocity, $\nabla u$, + * assuming that the viscosity is something that depends on the cell's material + * id. This can be done using a class we derive from DataPostprocessorScalar + * where we overload the DataPostprocessor::evaluate_vector_field() function + * that receives the values and gradients of the velocity (plus of + * other solution variables such as the pressure, but let's ignore those + * for the moment). Then we could use code such as this, assuming that we + * use a hp::DoFHandler: + * @code + * template + * class ComputeStress : public DataPostprocessorScalar + * { + * public: + * ... // overload other necessary member variables + * virtual + * void + * evaluate_vector_field (const DataPostprocessorInputs::Vector &input_data, + * std::vector > &computed_quantities) const + * { + * const typename hp::DoFHandler::cell_iterator + * current_cell = input_data.template get_cell >(); + * const viscosity = look_up_viscosity (current_cell->material_id()); + * + * for (unsigned int q=0; q + struct CommonInputs + { + /** + * An array of vectors normal to the faces of cells, evaluated at the points + * at which we are generating graphical output. This array is only used by + * the DataOutFaces class, and is left empty by all other classes for + * which the DataPostprocessor framework can be used. In the case of + * DataOutFaces, the array contains the outward normal vectors to the + * face, seen from the interior of the cell. + * + * This array is only filled if a user-derived class overloads the + * DataPostprocessor::get_needed_update_flags(), and the function + * returns (possibly among other flags) + * UpdateFlags::update_normal_vectors. Alternatively, a class + * derived from DataPostprocessorScalar or DataPostprocessorVector + * may pass this flag to the constructor of + * DataPostprocessorScalar or DataPostprocessorVector. + */ + std::vector > normals; + + /** + * An array of coordinates corresponding to the locations at which + * we are generating graphical output on one cell. + * + * This array is only filled if a user-derived class overloads the + * DataPostprocessor::get_needed_update_flags(), and the function + * returns (possibly among other flags) + * UpdateFlags::update_quadrature_points. Alternatively, a class + * derived from DataPostprocessorScalar or DataPostprocessorVector + * may pass this flag to the constructor of + * DataPostprocessorScalar or DataPostprocessorVector. + */ + std::vector > evaluation_points; + + /** + * Set the cell that is currently being used in evaluating the data + * for which the DataPostprocessor object is being called. + * + * This function is not usually called from user space, but is instead + * called by DataOut and similar classes when creating the object that + * is then passed to DataPostprocessor. + */ + template + void + set_cell (const typename DoFHandlerType::cell_iterator &cell); + + /** + * Query the cell on which we currently produce graphical output. + * See the documentation of the current class for an example on how + * to use this function. + */ + template + typename DoFHandlerType::cell_iterator + get_cell () const; + + private: + /** + * The place where set_cell() stores the cell. Since the actual data + * type of the cell iterator can be many different things, the + * interface uses boost::any here. This makes assignment in set_cell() + * simple, but requires knowing the data type of the stored object in + * get_cell(). + */ + boost::any cell; + }; + /** * A structure that is used to pass information to * DataPostprocessor::evaluate_scalar_field(). It contains * the values and (if requested) derivatives of a scalar solution - * variable at the evaluation points on a cell or face. If appropriate, - * it also contains the normal vectors to the geometry on which output - * is generated, at these evaluation points. + * variable at the evaluation points on a cell or face. + * + * Through the fields in the CommonInputs base class, this class also + * makes available access to the locations of evaluations points, + * normal vectors (if appropriate), and which cell data is currently + * being evaluated on (also if appropriate). + * + * @author Wolfgang Bangerth, 2016 */ template - struct Scalar + struct Scalar : public CommonInputs { /** * An array of values of the (scalar) solution at each of the evaluation @@ -84,38 +260,6 @@ namespace DataPostprocessorInputs * DataPostprocessorVector. */ std::vector > solution_hessians; - - /** - * An array of vectors normal to the faces of cells, evaluated at the points - * at which we are generating graphical output. This array is only used by - * the DataOutFaces class, and is left empty by all other classes for - * which the DataPostprocessor framework can be used. In the case of - * DataOutFaces, the array contains the outward normal vectors to the - * face, seen from the interior of the cell. - * - * This array is only filled if a user-derived class overloads the - * DataPostprocessor::get_needed_update_flags(), and the function - * returns (possibly among other flags) - * UpdateFlags::update_normal_vectors. Alternatively, a class - * derived from DataPostprocessorScalar or DataPostprocessorVector - * may pass this flag to the constructor of - * DataPostprocessorScalar or DataPostprocessorVector. - */ - std::vector > normals; - - /** - * An array of coordinates corresponding to the locations at which - * we are generating graphical output on one cell. - * - * This array is only filled if a user-derived class overloads the - * DataPostprocessor::get_needed_update_flags(), and the function - * returns (possibly among other flags) - * UpdateFlags::update_quadrature_points. Alternatively, a class - * derived from DataPostprocessorScalar or DataPostprocessorVector - * may pass this flag to the constructor of - * DataPostprocessorScalar or DataPostprocessorVector. - */ - std::vector > evaluation_points; }; @@ -124,12 +268,17 @@ namespace DataPostprocessorInputs * A structure that is used to pass information to * DataPostprocessor::evaluate_vector_field(). It contains * the values and (if requested) derivatives of a vector-valued solution - * variable at the evaluation points on a cell or face. If appropriate, - * it also contains the normal vectors to the geometry on which output - * is generated, at these evaluation points. + * variable at the evaluation points on a cell or face. + * + * Through the fields in the CommonInputs base class, this class also + * makes available access to the locations of evaluations points, + * normal vectors (if appropriate), and which cell data is currently + * being evaluated on (also if appropriate). + * + * @author Wolfgang Bangerth, 2016 */ template - struct Vector + struct Vector : public CommonInputs { /** * An array of values of a vector-valued solution at each of the evaluation @@ -179,39 +328,8 @@ namespace DataPostprocessorInputs * DataPostprocessorVector. */ std::vector > > solution_hessians; - - /** - * An array of vectors normal to the faces of cells, evaluated at the points - * at which we are generating graphical output. This array is only used by - * the DataOutFaces class, and is left empty by all other classes for - * which the DataPostprocessor framework can be used. In the case of - * DataOutFaces, the array contains the outward normal vectors to the - * face, seen from the interior of the cell. - * - * This array is only filled if a user-derived class overloads the - * DataPostprocessor::get_needed_update_flags(), and the function - * returns (possibly among other flags) - * UpdateFlags::update_normal_vectors. Alternatively, a class - * derived from DataPostprocessorScalar or DataPostprocessorVector - * may pass this flag to the constructor of - * DataPostprocessorScalar or DataPostprocessorVector. - */ - std::vector > normals; - - /** - * An array of coordinates corresponding to the locations at which - * we are generating graphical output on one cell. - * - * This array is only filled if a user-derived class overloads the - * DataPostprocessor::get_needed_update_flags(), and the function - * returns (possibly among other flags) - * UpdateFlags::update_quadrature_points. Alternatively, a class - * derived from DataPostprocessorScalar or DataPostprocessorVector - * may pass this flag to the constructor of - * DataPostprocessorScalar or DataPostprocessorVector. - */ - std::vector > evaluation_points; }; + } @@ -587,6 +705,56 @@ private: }; + +#ifndef DOXYGEN +// -------------------- template functions ---------------------- + +namespace DataPostprocessorInputs +{ + + template + template + void + CommonInputs::set_cell (const typename DoFHandlerType::cell_iterator &new_cell) + { + // see if we had previously already stored a cell that has the same + // data type; if so, reuse the memory location and avoid calling 'new' + // inside boost::any + if (typename DoFHandlerType::cell_iterator * storage_location + = boost::any_cast(&cell)) + *storage_location = new_cell; + else + // if we had nothing stored before, or if we had stored a different + // data type, just let boost::any replace things + cell = new_cell; + } + + + + template + template + typename DoFHandlerType::cell_iterator + CommonInputs::get_cell () const + { + Assert(cell.empty() == false, + ExcMessage ("You are trying to access the cell associated with a " + "DataPostprocessorInputs::Scalar object for which no cell has " + "been set.")); + Assert(boost::any_cast(&cell) != 0, + ExcMessage ("You are trying to access the cell associated with a " + "DataPostprocessorInputs::Scalar with a DoFHandler type that " + "is different from the type with which it has been set. For " + "example, if the cell for which output is currently being " + "generated belongs to a hp::DoFHandler<2,3> object, then you can " + "only call the current function with a template argument " + "equal to hp::DoFHandler<2,3>, but not with any other class " + "type or dimension template argument.")); + return boost::any_cast(cell); + } +} + +#endif + DEAL_II_NAMESPACE_CLOSE #endif -- 2.39.5