From fc022ee5768cc94d5dd934648402a2523a0129cc Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 31 Oct 2016 15:37:19 -0600 Subject: [PATCH] Update the way DataPostprocessor receives its input. In the past, the two functions we use to postprocess data got lists of arguments for the solution, its derivatives, the evaluation points, and normal vectors. This is not flexible enough: We can not easily add other information that we have needed in the past or that users have requested, such as a pointer to the cell we're currently on, or the material-id of the cell. Rather than adding each possible argument anyone may want to use individually to the list of the postprocessor function arguments, we use a technique here that has proven useful in ASPECT: the functions simply take a reference to a structure that contains these individual pieces of information. We can extend the members of these structures without backward compatibility issues because the functions still get a reference to the same structure, we just grow the structure itself. Functions that never used the new members of the structure will continue to work as they always did. --- include/deal.II/numerics/data_out_dof_data.h | 8 +- include/deal.II/numerics/data_out_faces.h | 1 - include/deal.II/numerics/data_postprocessor.h | 277 ++++++++++++++++-- source/numerics/data_out.cc | 40 +-- source/numerics/data_out_dof_data.cc | 38 ++- source/numerics/data_out_faces.cc | 53 ++-- source/numerics/data_out_rotation.cc | 39 +-- source/numerics/data_postprocessor.cc | 54 +++- 8 files changed, 375 insertions(+), 135 deletions(-) diff --git a/include/deal.II/numerics/data_out_dof_data.h b/include/deal.II/numerics/data_out_dof_data.h index 1c0e7e6bdc..f0c50264ae 100644 --- a/include/deal.II/numerics/data_out_dof_data.h +++ b/include/deal.II/numerics/data_out_dof_data.h @@ -345,12 +345,8 @@ namespace internal const unsigned int n_datasets; const unsigned int n_subdivisions; - std::vector patch_values; - std::vector > patch_values_system; - std::vector > patch_gradients; - std::vector > > patch_gradients_system; - std::vector > patch_hessians; - std::vector > > patch_hessians_system; + DataPostprocessorInputs::Scalar patch_values_scalar; + DataPostprocessorInputs::Vector patch_values_system; std::vector > > postprocessed_values; const dealii::hp::MappingCollection mapping_collection; diff --git a/include/deal.II/numerics/data_out_faces.h b/include/deal.II/numerics/data_out_faces.h index a82feb84fe..b89611b6ae 100644 --- a/include/deal.II/numerics/data_out_faces.h +++ b/include/deal.II/numerics/data_out_faces.h @@ -45,7 +45,6 @@ namespace internal const std::vector > > &finite_elements, const UpdateFlags update_flags); - std::vector > patch_normals; std::vector > patch_evaluation_points; }; } diff --git a/include/deal.II/numerics/data_postprocessor.h b/include/deal.II/numerics/data_postprocessor.h index 0baa696320..de652967af 100644 --- a/include/deal.II/numerics/data_postprocessor.h +++ b/include/deal.II/numerics/data_postprocessor.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2007 - 2015 by the deal.II authors +// Copyright (C) 2007 - 2016 by the deal.II authors // // This file is part of the deal.II library. // @@ -31,6 +31,190 @@ DEAL_II_NAMESPACE_OPEN +/** + * A namespace for data structures that are going to be passed from + * DataOut to the member functions of DataPostprocessor. + */ +namespace DataPostprocessorInputs +{ + /** + * A structure that is used to pass information to + * DataPostprocessor::compute_derived_quantities_scalar(). 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. + */ + template + struct Scalar + { + /** + * An array of values of the (scalar) solution at each of the evaluation + * points used to create graphical output from one cell, face, or other + * object. + */ + std::vector solution_values; + + /** + * An array of gradients of the (scalar) solution at each of the evaluation + * points used to create graphical output from one cell, face, or other + * object. + * + * 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_gradients. Alternatively, a class derived + * from DataPostprocessorScalar or DataPostprocessorVector may + * pass this flag to the constructor of DataPostprocessorScalar or + * DataPostprocessorVector. + */ + std::vector > solution_gradients; + + /** + * An array of second derivatives of the (scalar) solution at each of the evaluation + * points used to create graphical output from one cell, face, or other + * object. + * + * 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_hessians. Alternatively, a class derived + * from DataPostprocessorScalar or DataPostprocessorVector may + * pass this flag to the constructor of DataPostprocessorScalar or + * 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; + }; + + + + /** + * A structure that is used to pass information to + * DataPostprocessor::compute_derived_quantities_vector(). 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. + */ + template + struct Vector + { + /** + * An array of values of a vector-valued solution at each of the evaluation + * points used to create graphical output from one cell, face, or other + * object. + * + * The outer vector runs over the evaluation points, whereas the inner + * vector runs over the components of the finite element field for which + * output will be generated. + */ + std::vector > solution_values; + + /** + * An array of gradients of a vector-valued solution at each of the evaluation + * points used to create graphical output from one cell, face, or other + * object. + * + * The outer vector runs over the evaluation points, whereas the inner + * vector runs over the components of the finite element field for which + * output will be generated. + * + * 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_gradients. Alternatively, a class derived + * from DataPostprocessorScalar or DataPostprocessorVector may + * pass this flag to the constructor of DataPostprocessorScalar or + * DataPostprocessorVector. + */ + std::vector > > solution_gradients; + + /** + * An array of second derivatives of a vector-valued solution at each of the evaluation + * points used to create graphical output from one cell, face, or other + * object. + * + * The outer vector runs over the evaluation points, whereas the inner + * vector runs over the components of the finite element field for which + * output will be generated. + * + * 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_hessians. Alternatively, a class derived + * from DataPostprocessorScalar or DataPostprocessorVector may + * pass this flag to the constructor of DataPostprocessorScalar or + * 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; + }; +} + + /** * This class provides an interface to compute derived quantities from a * solution that can then be output in graphical formats for visualization, @@ -108,12 +292,13 @@ DEAL_II_NAMESPACE_OPEN * actually computes the results. * * @ingroup output - * @author Tobias Leicht, 2007 + * @author Tobias Leicht, 2007, Wolfgang Bangerth, 2016 */ template class DataPostprocessor: public Subscriptor { public: + /** * Destructor. This function doesn't actually do anything but is marked as * virtual to ensure that data postprocessors can be destroyed through @@ -123,26 +308,55 @@ public: /** * This is the main function which actually performs the postprocessing. The - * last argument is a reference to the postprocessed data which has correct - * size already and must be filled by this function. @p uh is a reference to - * a vector of data values at all points, @p duh the same for gradients, @p - * dduh for second derivatives and @p normals is a reference to the normal - * vectors. Note, that the last four references will only contain valid - * data, if the respective flags are returned by @p get_needed_update_flags, - * otherwise those vectors will be in an unspecified state. @p normals will - * always be an empty vector when working on cells, not on faces. + * second argument is a reference to the postprocessed data which already has + * correct size and must be filled by this function. + * + * The function takes the values, gradients, and higher derivatives of the + * solution at all evaluation points, as well as other data such as the + * cell, via the first argument. Not all of the member vectors of this + * argument will be filled with data -- in fact, derivatives and other + * quantities will only be contain valid data if the corresponding flags + * are returned by by an overloaded version of the get_needed_update_flags() + * function (implemented in a user's derived class). + * Otherwise those vectors will be in an unspecified state. * - * This function is called when the original data vector represents scalar + * This function is called when the finite element field that is being + * converted into graphical data by DataOut or similar classes represents scalar * data, i.e. the finite element in use has only a single vector component. */ virtual void - compute_derived_quantities_scalar (const std::vector &uh, - const std::vector > &duh, - const std::vector > &dduh, + compute_derived_quantities_scalar (const DataPostprocessorInputs::Scalar &input_data, + std::vector > &computed_quantities) const; + + /** + * @deprecated This function is deprecated. It has been superseded by + * function of same name above that receives a superset of the + * information provided to the current function through the members + * of the structure it receives as the first argument. + * + * If a user class derived from the current class (or from + * DataPostprocessorScalar) does not overload the function above, + * but instead overloads the current (legacy) form of the function, + * then the default implementation of the function above will simply + * call the current function. However, not all elements of the + * DataPostprocessorInputs::Scalar argument the function above + * receives have corresponding function arguments in the current + * function, and consequently not all information that function has + * available is passed on to the current one. In other words, there + * are pieces of information you may need in an implementation of a + * postprocess that are available if you overload the new form of this + * function above, but that are not available if you overload the old + * form of the function here. + */ + virtual + void + compute_derived_quantities_scalar (const std::vector &solution_values, + const std::vector > &solution_gradients, + const std::vector > &solution_hessians, const std::vector > &normals, const std::vector > &evaluation_points, - std::vector > &computed_quantities) const; + std::vector > &computed_quantities) const DEAL_II_DEPRECATED; /** * Same as the compute_derived_quantities_scalar() function, but this @@ -151,12 +365,37 @@ public: */ virtual void - compute_derived_quantities_vector (const std::vector > &uh, - const std::vector > > &duh, - const std::vector > > &dduh, + compute_derived_quantities_vector (const DataPostprocessorInputs::Vector &input_data, + std::vector > &computed_quantities) const; + + /** + * @deprecated This function is deprecated. It has been superseded by + * function of same name above that receives a superset of the + * information provided to the current function through the members + * of the structure it receives as the first argument. + * + * If a user class derived from the current class (or from + * DataPostprocessorVector) does not overload the function above, + * but instead overloads the current (legacy) form of the function, + * then the default implementation of the function above will simply + * call the current function. However, not all elements of the + * DataPostprocessorInputs::Vector argument the function above + * receives have corresponding function arguments in the current + * function, and consequently not all information that function has + * available is passed on to the current one. In other words, there + * are pieces of information you may need in an implementation of a + * postprocess that are available if you overload the new form of this + * function above, but that are not available if you overload the old + * form of the function here. + */ + virtual + void + compute_derived_quantities_vector (const std::vector > &solution_values, + const std::vector > > &solution_gradients, + const std::vector > > &solution_hessians, const std::vector > &normals, const std::vector > &evaluation_points, - std::vector > &computed_quantities) const; + std::vector > &computed_quantities) const DEAL_II_DEPRECATED; /** * Return the vector of strings describing the names of the computed diff --git a/source/numerics/data_out.cc b/source/numerics/data_out.cc index c63f0ee5f8..fdc4674f21 100644 --- a/source/numerics/data_out.cc +++ b/source/numerics/data_out.cc @@ -141,7 +141,8 @@ build_one_patch const unsigned int n_components = this_fe_patch_values.get_fe().n_components(); - const DataPostprocessor *postprocessor=this->dof_data[dataset]->postprocessor; + const DataPostprocessor *postprocessor + = this->dof_data[dataset]->postprocessor; if (postprocessor != 0) { @@ -154,25 +155,20 @@ build_one_patch // gradient etc. if (update_flags & update_values) this->dof_data[dataset]->get_function_values (this_fe_patch_values, - scratch_data.patch_values); + scratch_data.patch_values_scalar.solution_values); if (update_flags & update_gradients) this->dof_data[dataset]->get_function_gradients (this_fe_patch_values, - scratch_data.patch_gradients); + scratch_data.patch_values_scalar.solution_gradients); if (update_flags & update_hessians) this->dof_data[dataset]->get_function_hessians (this_fe_patch_values, - scratch_data.patch_hessians); + scratch_data.patch_values_scalar.solution_hessians); if (update_flags & update_quadrature_points) - scratch_data.patch_evaluation_points = this_fe_patch_values.get_quadrature_points(); + scratch_data.patch_values_scalar.evaluation_points = this_fe_patch_values.get_quadrature_points(); - std::vector > dummy_normals; postprocessor-> - compute_derived_quantities_scalar(scratch_data.patch_values, - scratch_data.patch_gradients, - scratch_data.patch_hessians, - dummy_normals, - scratch_data.patch_evaluation_points, + compute_derived_quantities_scalar(scratch_data.patch_values_scalar, scratch_data.postprocessed_values[dataset]); } else @@ -183,25 +179,19 @@ build_one_patch // derivative... if (update_flags & update_values) this->dof_data[dataset]->get_function_values (this_fe_patch_values, - scratch_data.patch_values_system); + scratch_data.patch_values_system.solution_values); if (update_flags & update_gradients) this->dof_data[dataset]->get_function_gradients (this_fe_patch_values, - scratch_data.patch_gradients_system); + scratch_data.patch_values_system.solution_gradients); if (update_flags & update_hessians) this->dof_data[dataset]->get_function_hessians (this_fe_patch_values, - scratch_data.patch_hessians_system); + scratch_data.patch_values_system.solution_hessians); if (update_flags & update_quadrature_points) - scratch_data.patch_evaluation_points = this_fe_patch_values.get_quadrature_points(); - - std::vector > dummy_normals; + scratch_data.patch_values_system.evaluation_points = this_fe_patch_values.get_quadrature_points(); postprocessor-> compute_derived_quantities_vector(scratch_data.patch_values_system, - scratch_data.patch_gradients_system, - scratch_data.patch_hessians_system, - dummy_normals, - scratch_data.patch_evaluation_points, scratch_data.postprocessed_values[dataset]); } @@ -219,20 +209,20 @@ build_one_patch if (n_components == 1) { this->dof_data[dataset]->get_function_values (this_fe_patch_values, - scratch_data.patch_values); + scratch_data.patch_values_scalar.solution_values); for (unsigned int q=0; qdof_data[dataset]->get_function_values (this_fe_patch_values, - scratch_data.patch_values_system); + scratch_data.patch_values_system.solution_values); for (unsigned int component=0; componentdof_data[dataset]->n_output_variables; diff --git a/source/numerics/data_out_dof_data.cc b/source/numerics/data_out_dof_data.cc index c28967aeac..29f5081508 100644 --- a/source/numerics/data_out_dof_data.cc +++ b/source/numerics/data_out_dof_data.cc @@ -116,12 +116,12 @@ namespace internal } } - patch_values.resize (n_q_points); - patch_values_system.resize (n_q_points); - patch_gradients.resize (n_q_points); - patch_gradients_system.resize (n_q_points); - patch_hessians.resize (n_q_points); - patch_hessians_system.resize (n_q_points); + patch_values_scalar.solution_values.resize (n_q_points); + patch_values_scalar.solution_gradients.resize (n_q_points); + patch_values_scalar.solution_hessians.resize (n_q_points); + patch_values_system.solution_values.resize (n_q_points); + patch_values_system.solution_gradients.resize (n_q_points); + patch_values_system.solution_hessians.resize (n_q_points); for (unsigned int dataset=0; dataset:: resize_system_vectors(const unsigned int n_components) { - Assert(patch_values_system.size() > 0, ExcInternalError()); - AssertDimension(patch_values_system.size(), - patch_gradients_system.size()); - AssertDimension(patch_values_system.size(), - patch_hessians_system.size()); - if (patch_values_system[0].size() == n_components) + Assert(patch_values_system.solution_values.size() > 0, ExcInternalError()); + AssertDimension(patch_values_system.solution_values.size(), + patch_values_system.solution_gradients.size()); + AssertDimension(patch_values_system.solution_values.size(), + patch_values_system.solution_hessians.size()); + if (patch_values_system.solution_values[0].size() == n_components) return; - for (unsigned int k=0; kget_needed_update_flags(); - // get normals, if needed. this is a geometrical information and - // thus does not depend on the number of components of the data - // vector - if (update_flags & update_normal_vectors) - { -//TODO: undo this copying when we can change the data type of -// data.patch_normals to Tensor<1,spacedim> as well - for (unsigned int q=0; q(this_fe_patch_values.get_all_normal_vectors()[q]); - } - if (n_components == 1) { // at each point there is only one component of value, // gradient etc. if (update_flags & update_values) this->dof_data[dataset]->get_function_values (this_fe_patch_values, - data.patch_values); + data.patch_values_scalar.solution_values); if (update_flags & update_gradients) this->dof_data[dataset]->get_function_gradients (this_fe_patch_values, - data.patch_gradients); + data.patch_values_scalar.solution_gradients); if (update_flags & update_hessians) this->dof_data[dataset]->get_function_hessians (this_fe_patch_values, - data.patch_hessians); + data.patch_values_scalar.solution_hessians); if (update_flags & update_quadrature_points) - data.patch_evaluation_points = this_fe_patch_values.get_quadrature_points(); + data.patch_values_scalar.evaluation_points = this_fe_patch_values.get_quadrature_points(); + + if (update_flags & update_normal_vectors) + data.patch_values_scalar.normals = this_fe_patch_values.get_all_normal_vectors(); postprocessor-> - compute_derived_quantities_scalar(data.patch_values, - data.patch_gradients, - data.patch_hessians, - data.patch_normals, - data.patch_evaluation_points, + compute_derived_quantities_scalar(data.patch_values_scalar, data.postprocessed_values[dataset]); } else @@ -192,23 +180,22 @@ build_one_patch (const FaceDescriptor *cell_and_face, data.resize_system_vectors(n_components); if (update_flags & update_values) this->dof_data[dataset]->get_function_values (this_fe_patch_values, - data.patch_values_system); + data.patch_values_system.solution_values); if (update_flags & update_gradients) this->dof_data[dataset]->get_function_gradients (this_fe_patch_values, - data.patch_gradients_system); + data.patch_values_system.solution_gradients); if (update_flags & update_hessians) this->dof_data[dataset]->get_function_hessians (this_fe_patch_values, - data.patch_hessians_system); + data.patch_values_system.solution_hessians); if (update_flags & update_quadrature_points) - data.patch_evaluation_points = this_fe_patch_values.get_quadrature_points(); + data.patch_values_system.evaluation_points = this_fe_patch_values.get_quadrature_points(); + + if (update_flags & update_normal_vectors) + data.patch_values_system.normals = this_fe_patch_values.get_all_normal_vectors(); postprocessor-> compute_derived_quantities_vector(data.patch_values_system, - data.patch_gradients_system, - data.patch_hessians_system, - data.patch_normals, - data.patch_evaluation_points, data.postprocessed_values[dataset]); } @@ -225,23 +212,23 @@ build_one_patch (const FaceDescriptor *cell_and_face, if (n_components == 1) { this->dof_data[dataset]->get_function_values (this_fe_patch_values, - data.patch_values); + data.patch_values_scalar.solution_values); for (unsigned int q=0; qdof_data[dataset]->get_function_values (this_fe_patch_values, - data.patch_values_system); + data.patch_values_system.solution_values); for (unsigned int component=0; componentdof_data[dataset]->n_output_variables; + offset += this->dof_data[dataset]->n_output_variables; } // then do the cell data diff --git a/source/numerics/data_out_rotation.cc b/source/numerics/data_out_rotation.cc index 653b6d2766..6812cebc5a 100644 --- a/source/numerics/data_out_rotation.cc +++ b/source/numerics/data_out_rotation.cc @@ -205,24 +205,19 @@ build_one_patch (const cell_iterator *cell, // value, gradient etc. if (update_flags & update_values) this->dof_data[dataset]->get_function_values (fe_patch_values, - data.patch_values); + data.patch_values_scalar.solution_values); if (update_flags & update_gradients) this->dof_data[dataset]->get_function_gradients (fe_patch_values, - data.patch_gradients); + data.patch_values_scalar.solution_gradients); if (update_flags & update_hessians) this->dof_data[dataset]->get_function_hessians (fe_patch_values, - data.patch_hessians); + data.patch_values_scalar.solution_hessians); if (update_flags & update_quadrature_points) - data.patch_evaluation_points = fe_patch_values.get_quadrature_points(); + data.patch_values_scalar.evaluation_points = fe_patch_values.get_quadrature_points(); - std::vector > dummy_normals; postprocessor-> - compute_derived_quantities_scalar(data.patch_values, - data.patch_gradients, - data.patch_hessians, - dummy_normals, - data.patch_evaluation_points, + compute_derived_quantities_scalar(data.patch_values_scalar, data.postprocessed_values[dataset]); } else @@ -233,24 +228,20 @@ build_one_patch (const cell_iterator *cell, // its derivative... if (update_flags & update_values) this->dof_data[dataset]->get_function_values (fe_patch_values, - data.patch_values_system); + data.patch_values_system.solution_values); if (update_flags & update_gradients) this->dof_data[dataset]->get_function_gradients (fe_patch_values, - data.patch_gradients_system); + data.patch_values_system.solution_gradients); if (update_flags & update_hessians) this->dof_data[dataset]->get_function_hessians (fe_patch_values, - data.patch_hessians_system); + data.patch_values_system.solution_hessians); if (update_flags & update_quadrature_points) - data.patch_evaluation_points = fe_patch_values.get_quadrature_points(); + data.patch_values_system.evaluation_points = fe_patch_values.get_quadrature_points(); std::vector > dummy_normals; postprocessor-> compute_derived_quantities_vector(data.patch_values_system, - data.patch_gradients_system, - data.patch_hessians_system, - dummy_normals, - data.patch_evaluation_points, data.postprocessed_values[dataset]); } @@ -287,7 +278,7 @@ build_one_patch (const cell_iterator *cell, else if (n_components == 1) { this->dof_data[dataset]->get_function_values (fe_patch_values, - data.patch_values); + data.patch_values_scalar.solution_values); switch (dimension) { @@ -296,7 +287,7 @@ build_one_patch (const cell_iterator *cell, for (unsigned int y=0; ydof_data[dataset]->get_function_values (fe_patch_values, - data.patch_values_system); + data.patch_values_system.solution_values); for (unsigned int component=0; component::~DataPostprocessor() template void DataPostprocessor:: -compute_derived_quantities_scalar (const std::vector &/*uh*/, - const std::vector > &/*duh*/, - const std::vector > &/*dduh*/, +compute_derived_quantities_scalar (const DataPostprocessorInputs::Scalar &inputs, + std::vector > &computed_quantities) const +{ + // for backward compatibility, call the old function. + // this also requires converting the accidental use + // of Point for normal vectors + std::vector > normals (inputs.normals.begin(), + inputs.normals.end()); + compute_derived_quantities_scalar(inputs.solution_values, + inputs.solution_gradients, + inputs.solution_hessians, + normals, + inputs.evaluation_points, + computed_quantities); +} + + + +template +void +DataPostprocessor:: +compute_derived_quantities_scalar (const std::vector &/*solution_values*/, + const std::vector > &/*solution_gradients*/, + const std::vector > &/*solution_hessians*/, const std::vector > &/*normals*/, const std::vector > &/*evaluation_points*/, std::vector > &computed_quantities) const @@ -46,9 +67,30 @@ compute_derived_quantities_scalar (const std::vector &/*uh*/, template void DataPostprocessor:: -compute_derived_quantities_vector (const std::vector > &/*uh*/, - const std::vector > > &/*duh*/, - const std::vector > > &/*dduh*/, +compute_derived_quantities_vector (const DataPostprocessorInputs::Vector &inputs, + std::vector > &computed_quantities) const +{ + // for backward compatibility, call the old function. + // this also requires converting the accidental use + // of Point for normal vectors + std::vector > normals (inputs.normals.begin(), + inputs.normals.end()); + compute_derived_quantities_vector(inputs.solution_values, + inputs.solution_gradients, + inputs.solution_hessians, + normals, + inputs.evaluation_points, + computed_quantities); +} + + + +template +void +DataPostprocessor:: +compute_derived_quantities_vector (const std::vector > &/*solution_values*/, + const std::vector > > &/*solution_gradients*/, + const std::vector > > &/*solution_hessians*/, const std::vector > &/*normals*/, const std::vector > &/*evaluation_points*/, std::vector > &computed_quantities) const -- 2.39.5