From: Wolfgang Bangerth Date: Tue, 29 Aug 2017 23:17:47 +0000 (-0600) Subject: Provide an example for DataPostprocessorVector. X-Git-Tag: v9.0.0-rc1~1142^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=45c25714473c9e702f408e7a5b223f6ce86de14c;p=dealii.git Provide an example for DataPostprocessorVector. --- diff --git a/doc/doxygen/images/data_postprocessor_vector_0.png b/doc/doxygen/images/data_postprocessor_vector_0.png new file mode 100644 index 0000000000..f3011a70fb Binary files /dev/null and b/doc/doxygen/images/data_postprocessor_vector_0.png differ diff --git a/doc/doxygen/images/data_postprocessor_vector_1.png b/doc/doxygen/images/data_postprocessor_vector_1.png new file mode 100644 index 0000000000..aa8e776c5a Binary files /dev/null and b/doc/doxygen/images/data_postprocessor_vector_1.png differ diff --git a/doc/doxygen/images/data_postprocessor_vector_2.png b/doc/doxygen/images/data_postprocessor_vector_2.png new file mode 100644 index 0000000000..bfda9a383c Binary files /dev/null and b/doc/doxygen/images/data_postprocessor_vector_2.png differ diff --git a/include/deal.II/numerics/data_postprocessor.h b/include/deal.II/numerics/data_postprocessor.h index 243dc62596..8244dd9271 100644 --- a/include/deal.II/numerics/data_postprocessor.h +++ b/include/deal.II/numerics/data_postprocessor.h @@ -591,8 +591,172 @@ private: * An example of how the closely related class DataPostprocessorScalar is used * can be found in step-29. * + * + *

An example

+ * + * A common example of what one wants to do with postprocessors is to visualize + * not just the value of the solution, but the gradient. Let's, for simplicity, + * assume that you have only a scalar solution. In fact, because it's readily + * available, let us simply take the step-6 solver to produce such a scalar + * solution. The gradient is a vector (with exactly @p dim components), so the + * current class fits the bill to produce the gradient through postprocessing. + * Then, the following code snippet implements everything you need to have + * to visualize the gradient: + * @code + * template + * class GradientPostprocessor : public DataPostprocessorVector + * { + * public: + * GradientPostprocessor () + * : + * // call the constructor of the base class. call the variable to + * // be output "grad_u" and make sure that DataOut provides us + * // with the gradients: + * DataPostprocessorVector ("grad_u", + * update_gradients) + * {} + * + * virtual + * void + * evaluate_scalar_field (const DataPostprocessorInputs::Scalar &input_data, + * std::vector > &computed_quantities) const + * { + * // ensure that there really are as many output slots + * // as there are points at which DataOut provides the + * // gradients: + * AssertDimension (input_data.solution_gradients.size(), + * computed_quantities.size()); + * + * // then loop over all of these inputs: + * for (unsigned int p=0; p gradient_postprocessor; + * + * DataOut data_out; + * data_out.attach_dof_handler (dof_handler); + * data_out.add_data_vector (solution, "solution"); + * data_out.add_data_vector (solution, gradient_postprocessor); + * data_out.build_patches (); + * + * std::ofstream output ("solution.vtu"); + * data_out.write_vtu (output); + * @endcode + * + * This leads to the following output for the solution and the + * gradients (you may want to compare with the solution shown in the results + * section of step-6; the current data is generated on a coarser mesh for + * simplicity): + * + * @image html data_postprocessor_vector_0.png + * @image html data_postprocessor_vector_1.png + * + * In the second image, the background color corresponds to the magnitude of the + * gradient vector and the vector glyphs to the gradient itself. It may be surprising + * at first to see that from each vertex, multiple vectors originate, going in + * different directions. But that is because the solution is only continuous: + * in general, the gradient is discontinuous across edges, and so the multiple + * vectors originating from each vertex simply represent the differing + * gradients of the solution at each adjacent cell. + * + * The output above -- namely, the gradient $\nabla u$ of the solution -- + * corresponds to the temperature gradient if one interpreted step-6 as solving + * a steady-state heat transfer problem. + * It is very small in the central part of the domain because in step-6 we are + * solving an equation that has a coefficient $a(\mathbf x)$ that is large in + * the central part and small on the outside. This can be thought as a material + * that conducts heat well, and consequently the temperature gradient is small. + * On the other hand, the "heat flux" corresponds to the quantity + * $a(\mathbf x) \nabla u(\mathbf x)$. For the + * solution of that equation, the flux should be continuous across the interface. + * This is easily verified by the following modification of the postprocessor: + * @code + * template + * class HeatFluxPostprocessor : public DataPostprocessorVector + * { + * public: + * HeatFluxPostprocessor () + * : + * // like above, but now also make sure that DataOut provides + * // us with coordinates of the evaluation points: + * DataPostprocessorVector ("heat_flux", + * update_gradients | update_quadrature_points) + * {} + * + * virtual + * void + * evaluate_scalar_field (const DataPostprocessorInputs::Scalar &input_data, + * std::vector > &computed_quantities) const + * { + * AssertDimension (input_data.solution_gradients.size(), + * computed_quantities.size()); + * + * for (unsigned int p=0; p Extension to the gradients of vector-valued problems + * + * The example above uses a scalar solution and its gradient as an example. + * On the other hand, one may want to do something similar for the gradient + * of a vector-valued displacement field (such as the strain or + * stress of a displacement field, like those computed in step-8, step-17, + * step-18, or step-44). In that case, the solution is already vector + * valued and the stress is a (symmetric) tensor. + * + * deal.II does not currently support outputting tensor-valued quantities, but + * they can of course be output as a collection of scalar-valued components of + * the tensor. This means that the current class is not applicable any more, + * but it is not very difficult to derive a true "stress" postprocessor directly + * from the DataPostprocessor class that simply outputs dim*dim + * components of the stress vector as scalars and that is structured in a + * similar way to the postprocessors above. (It has to overload the + * @p evaluate_vector_field() function, however, given that the solution is + * vector valued already.) + * + * * @ingroup output - * @author Wolfgang Bangerth, 2011 + * @author Wolfgang Bangerth, 2011, 2017 */ template class DataPostprocessorVector : public DataPostprocessor