]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update the way DataPostprocessor receives its input.
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 31 Oct 2016 21:37:19 +0000 (15:37 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 2 Nov 2016 22:27:33 +0000 (16:27 -0600)
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
include/deal.II/numerics/data_out_faces.h
include/deal.II/numerics/data_postprocessor.h
source/numerics/data_out.cc
source/numerics/data_out_dof_data.cc
source/numerics/data_out_faces.cc
source/numerics/data_out_rotation.cc
source/numerics/data_postprocessor.cc

index 1c0e7e6bdc067e682bd735fb55e09aa66b9f5631..f0c50264ae49e4030544f392eb9454c976b3d5ca 100644 (file)
@@ -345,12 +345,8 @@ namespace internal
       const unsigned int n_datasets;
       const unsigned int n_subdivisions;
 
-      std::vector<double>                                patch_values;
-      std::vector<dealii::Vector<double> >               patch_values_system;
-      std::vector<Tensor<1,spacedim> >                   patch_gradients;
-      std::vector<std::vector<Tensor<1,spacedim> > >     patch_gradients_system;
-      std::vector<Tensor<2,spacedim> >                   patch_hessians;
-      std::vector<std::vector<Tensor<2,spacedim> > >     patch_hessians_system;
+      DataPostprocessorInputs::Scalar<spacedim>          patch_values_scalar;
+      DataPostprocessorInputs::Vector<spacedim>          patch_values_system;
       std::vector<std::vector<dealii::Vector<double> > > postprocessed_values;
 
       const dealii::hp::MappingCollection<dim,spacedim> mapping_collection;
index a82feb84fea5cf7449aea595fc946f6618a46771..b89611b6aef48a1a9ba7d467e3665c86a2584cdb 100644 (file)
@@ -45,7 +45,6 @@ namespace internal
                     const std::vector<std_cxx11::shared_ptr<dealii::hp::FECollection<dim,spacedim> > > &finite_elements,
                     const UpdateFlags update_flags);
 
-      std::vector<Point<dim> > patch_normals;
       std::vector<Point<spacedim> > patch_evaluation_points;
     };
   }
index 0baa69632098882583c42ee9d23b4386b07ef91c..de652967afa8e329a7c7d476a55a27757ef552de 100644 (file)
@@ -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.
 //
 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 <int spacedim>
+  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<double>               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<Tensor<1, spacedim> > 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<Tensor<2, spacedim> > 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<Tensor<1, spacedim> > 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<Point<spacedim> >     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 <int spacedim>
+  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<dealii::Vector<double> >            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<std::vector<Tensor<1, spacedim> > > 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<std::vector<Tensor<2, spacedim> > > 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<Tensor<1, spacedim> >               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<Point<spacedim> >                   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 <int dim>
 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<double>         &uh,
-                                     const std::vector<Tensor<1,dim> > &duh,
-                                     const std::vector<Tensor<2,dim> > &dduh,
+  compute_derived_quantities_scalar (const DataPostprocessorInputs::Scalar<dim> &input_data,
+                                     std::vector<Vector<double> >               &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<double>         &solution_values,
+                                     const std::vector<Tensor<1,dim> > &solution_gradients,
+                                     const std::vector<Tensor<2,dim> > &solution_hessians,
                                      const std::vector<Point<dim> >    &normals,
                                      const std::vector<Point<dim> >    &evaluation_points,
-                                     std::vector<Vector<double> >      &computed_quantities) const;
+                                     std::vector<Vector<double> >      &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<Vector<double> >              &uh,
-                                     const std::vector<std::vector<Tensor<1,dim> > > &duh,
-                                     const std::vector<std::vector<Tensor<2,dim> > > &dduh,
+  compute_derived_quantities_vector (const DataPostprocessorInputs::Vector<dim> &input_data,
+                                     std::vector<Vector<double> >               &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<Vector<double> >              &solution_values,
+                                     const std::vector<std::vector<Tensor<1,dim> > > &solution_gradients,
+                                     const std::vector<std::vector<Tensor<2,dim> > > &solution_hessians,
                                      const std::vector<Point<dim> >                  &normals,
                                      const std::vector<Point<dim> >                  &evaluation_points,
-                                     std::vector<Vector<double> >                    &computed_quantities) const;
+                                     std::vector<Vector<double> >                    &computed_quantities) const DEAL_II_DEPRECATED;
 
   /**
    * Return the vector of strings describing the names of the computed
index c63f0ee5f8be2261033b7964d02502215cf02055..fdc4674f21e9cd2b761f3e3e34c76ce6e471f0b9 100644 (file)
@@ -141,7 +141,8 @@ build_one_patch
           const unsigned int n_components =
             this_fe_patch_values.get_fe().n_components();
 
-          const DataPostprocessor<DoFHandlerType::space_dimension> *postprocessor=this->dof_data[dataset]->postprocessor;
+          const DataPostprocessor<DoFHandlerType::space_dimension> *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<Point<DoFHandlerType::space_dimension> > 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<Point<DoFHandlerType::space_dimension> > 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; q<n_q_points; ++q)
-                  patch.data(offset,q) = scratch_data.patch_values[q];
+                  patch.data(offset,q) = scratch_data.patch_values_scalar.solution_values[q];
               }
             else
               {
                 scratch_data.resize_system_vectors(n_components);
                 this->dof_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; component<n_components;
                      ++component)
                   for (unsigned int q=0; q<n_q_points; ++q)
                     patch.data(offset+component,q) =
-                      scratch_data.patch_values_system[q](component);
+                      scratch_data.patch_values_system.solution_values[q](component);
               }
           // increment the counter for the actual data record
           offset+=this->dof_data[dataset]->n_output_variables;
index c28967aeac5c5dd73e6f46e027d9bfb5e9329dbe..29f50815087b5f6c5b8079397eae77eb295f1e10 100644 (file)
@@ -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<n_postprocessor_outputs.size(); ++dataset)
         if (n_postprocessor_outputs[dataset] != 0)
@@ -142,12 +142,8 @@ namespace internal
       :
       n_datasets (data.n_datasets),
       n_subdivisions (data.n_subdivisions),
-      patch_values (data.patch_values),
+      patch_values_scalar (data.patch_values_scalar),
       patch_values_system (data.patch_values_system),
-      patch_gradients (data.patch_gradients),
-      patch_gradients_system (data.patch_gradients_system),
-      patch_hessians (data.patch_hessians),
-      patch_hessians_system (data.patch_hessians_system),
       postprocessed_values (data.postprocessed_values),
       mapping_collection (data.mapping_collection),
       finite_elements (data.finite_elements),
@@ -270,18 +266,18 @@ namespace internal
     ParallelDataBase<dim,spacedim>::
     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; k<patch_values_system.size(); ++k)
+      for (unsigned int k=0; k<patch_values_system.solution_values.size(); ++k)
         {
-          patch_values_system[k].reinit(n_components);
-          patch_gradients_system[k].resize(n_components);
-          patch_hessians_system[k].resize(n_components);
+          patch_values_system.solution_values[k].reinit(n_components);
+          patch_values_system.solution_gradients[k].resize(n_components);
+          patch_values_system.solution_hessians[k].resize(n_components);
         }
     }
 
index 15284f736a7b97961e7bb8dec06f830166ec0b70..29e98f417315f16d4bbfba6755c43d5c4ea84bfe 100644 (file)
@@ -149,40 +149,28 @@ build_one_patch (const FaceDescriptor *cell_and_face,
               // have to be updated
               const UpdateFlags update_flags=postprocessor->get_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.n_quadrature_points; ++q)
-                    data.patch_normals[q] = Point<dim>(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; q<n_q_points; ++q)
-                  patch.data(offset,q) = data.patch_values[q];
+                  patch.data(offset,q) = data.patch_values_scalar.solution_values[q];
               }
             else
               {
                 data.resize_system_vectors(n_components);
                 this->dof_data[dataset]->get_function_values (this_fe_patch_values,
-                                                              data.patch_values_system);
+                                                              data.patch_values_system.solution_values);
                 for (unsigned int component=0; component<n_components;
                      ++component)
                   for (unsigned int q=0; q<n_q_points; ++q)
                     patch.data(offset+component,q) =
-                      data.patch_values_system[q](component);
+                      data.patch_values_system.solution_values[q](component);
               }
           // increment the counter for the actual data record
-          offset+=this->dof_data[dataset]->n_output_variables;
+          offset += this->dof_data[dataset]->n_output_variables;
         }
 
       // then do the cell data
index 653b6d27667a1c1ab0fcebb71a1f5c2ba0237565..6812cebc5a0ccc5845b93f26ff302fd65fe018c8 100644 (file)
@@ -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<Point<space_dimension> > 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<Point<space_dimension> > 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; y<n_points; ++y)
                           patches[angle].data(offset,
                                               x*n_points + y)
-                            = data.patch_values[x];
+                            = data.patch_values_scalar.solution_values[x];
                       break;
 
                     case 2:
@@ -307,7 +298,7 @@ build_one_patch (const cell_iterator *cell,
                                                 x*n_points*n_points +
                                                 y +
                                                 z*n_points)
-                              = data.patch_values[x*n_points+z];
+                              = data.patch_values_scalar.solution_values[x*n_points+z];
                       break;
 
                     default:
@@ -319,7 +310,7 @@ build_one_patch (const cell_iterator *cell,
                 {
                   data.resize_system_vectors(n_components);
                   this->dof_data[dataset]->get_function_values (fe_patch_values,
-                                                                data.patch_values_system);
+                                                                data.patch_values_system.solution_values);
 
                   for (unsigned int component=0; component<n_components;
                        ++component)
@@ -331,7 +322,7 @@ build_one_patch (const cell_iterator *cell,
                             for (unsigned int y=0; y<n_points; ++y)
                               patches[angle].data(offset+component,
                                                   x*n_points + y)
-                                = data.patch_values_system[x](component);
+                                = data.patch_values_system.solution_values[x](component);
                           break;
 
                         case 2:
@@ -342,7 +333,7 @@ build_one_patch (const cell_iterator *cell,
                                                     x*n_points*n_points +
                                                     y*n_points +
                                                     z)
-                                  = data.patch_values_system[x*n_points+z](component);
+                                  = data.patch_values_system.solution_values[x*n_points+z](component);
                           break;
 
                         default:
index 41c7cda412537ff00dd2397c7edaa3b081a4ecf8..a3ab1cb0928a1e15b0e0461478855bd91924a476 100644 (file)
@@ -30,9 +30,30 @@ DataPostprocessor<dim>::~DataPostprocessor()
 template <int dim>
 void
 DataPostprocessor<dim>::
-compute_derived_quantities_scalar (const std::vector<double>         &/*uh*/,
-                                   const std::vector<Tensor<1,dim> > &/*duh*/,
-                                   const std::vector<Tensor<2,dim> > &/*dduh*/,
+compute_derived_quantities_scalar (const DataPostprocessorInputs::Scalar<dim> &inputs,
+                                   std::vector<Vector<double> >               &computed_quantities) const
+{
+  // for backward compatibility, call the old function.
+  // this also requires converting the accidental use
+  // of Point<dim> for normal vectors
+  std::vector<Point<dim> > 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 <int dim>
+void
+DataPostprocessor<dim>::
+compute_derived_quantities_scalar (const std::vector<double>         &/*solution_values*/,
+                                   const std::vector<Tensor<1,dim> > &/*solution_gradients*/,
+                                   const std::vector<Tensor<2,dim> > &/*solution_hessians*/,
                                    const std::vector<Point<dim> >    &/*normals*/,
                                    const std::vector<Point<dim> >    &/*evaluation_points*/,
                                    std::vector<Vector<double> >      &computed_quantities) const
@@ -46,9 +67,30 @@ compute_derived_quantities_scalar (const std::vector<double>         &/*uh*/,
 template <int dim>
 void
 DataPostprocessor<dim>::
-compute_derived_quantities_vector (const std::vector<Vector<double> > &/*uh*/,
-                                   const std::vector<std::vector<Tensor<1,dim> > > &/*duh*/,
-                                   const std::vector<std::vector<Tensor<2,dim> > > &/*dduh*/,
+compute_derived_quantities_vector (const DataPostprocessorInputs::Vector<dim> &inputs,
+                                   std::vector<Vector<double> >               &computed_quantities) const
+{
+  // for backward compatibility, call the old function.
+  // this also requires converting the accidental use
+  // of Point<dim> for normal vectors
+  std::vector<Point<dim> > 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 <int dim>
+void
+DataPostprocessor<dim>::
+compute_derived_quantities_vector (const std::vector<Vector<double> > &/*solution_values*/,
+                                   const std::vector<std::vector<Tensor<1,dim> > > &/*solution_gradients*/,
+                                   const std::vector<std::vector<Tensor<2,dim> > > &/*solution_hessians*/,
                                    const std::vector<Point<dim> >                  &/*normals*/,
                                    const std::vector<Point<dim> >                  &/*evaluation_points*/,
                                    std::vector<Vector<double> >                    &computed_quantities) const

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.