/**
* For the (graphical) output of a FE solution one frequently wants to include
- * derived quantities, which are calculated from the values of the solution and
- * possibly the first and second derivates of the solution. This class offers
- * the interface to perform such a postprocessing. Given the values and
- * derivatives of provided dataon given points of a cell, new quantities can be
- * calculated.
+ * derived quantities, which are calculated from the values of the solution
+ * and possibly the first and second derivates of the solution. Examples are
+ * the calculation Mach numbers from velocity and density in supersonic flow
+ * computations, or the computation of the magnitude of a complex-valued
+ * solution as demonstrated in @ref step_29 "step-29". This class offers the
+ * interface to perform such postprocessing. Given the values and derivatives
+ * of the solution at those points where we want to generated output, the
+ * functions of this class can be overloaded to compute new quantities.
*
* A data vector and an object of a derived class can be given to the
* <tt>DataOut::add_data_vector</tt> function, which will write the derived
/**
* This is the main function which actually
- * performs the postprocessing. The first
+ * performs the postprocessing. The last
* argument is a reference to the
* postprocessed data which has correct
* size already and must be filled by this
*/
virtual
void
- compute_derived_quantities_scalar (std::vector<Vector<double> > &computed_quantities,
- const std::vector<double> &uh,
+ compute_derived_quantities_scalar (const std::vector<double> &uh,
const std::vector<Tensor<1,dim> > &duh,
const std::vector<Tensor<2,dim> > &dduh,
- const std::vector<Point<dim> > &normals) const;
+ const std::vector<Point<dim> > &normals,
+ std::vector<Vector<double> > &computed_quantities) const;
/**
* Same as above function, but
*/
virtual
void
- compute_derived_quantities_vector (std::vector<Vector<double> > &computed_quantities,
- const std::vector<Vector<double> > &uh,
+ 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,
- const std::vector<Point<dim> > &normals) const;
+ const std::vector<Point<dim> > &normals,
+ std::vector<Vector<double> > &computed_quantities) const;
/**
* Return the vector of strings describing
this->dof_data[dataset]->get_function_second_derivatives (fe_patch_values,
data.patch_second_derivatives);
postprocessor->
- compute_derived_quantities_scalar(data.postprocessed_values[dataset],
- data.patch_values,
+ compute_derived_quantities_scalar(data.patch_values,
data.patch_gradients,
data.patch_second_derivatives,
- data.dummy_normals);
+ data.dummy_normals,
+ data.postprocessed_values[dataset]);
}
else
{
this->dof_data[dataset]->get_function_second_derivatives (fe_patch_values,
data.patch_second_derivatives_system);
postprocessor->
- compute_derived_quantities_vector(data.postprocessed_values[dataset],
- data.patch_values_system,
+ compute_derived_quantities_vector(data.patch_values_system,
data.patch_gradients_system,
data.patch_second_derivatives_system,
- data.dummy_normals);
+ data.dummy_normals,
+ data.postprocessed_values[dataset]);
}
for (unsigned int q=0; q<n_q_points; ++q)
this->dof_data[dataset]->get_function_second_derivatives (fe_patch_values,
data.patch_second_derivatives);
postprocessor->
- compute_derived_quantities_scalar(data.postprocessed_values[dataset],
- data.patch_values,
+ compute_derived_quantities_scalar(data.patch_values,
data.patch_gradients,
data.patch_second_derivatives,
- data.patch_normals);
+ data.patch_normals,
+ data.postprocessed_values[dataset]);
}
else
{
this->dof_data[dataset]->get_function_second_derivatives (fe_patch_values,
data.patch_second_derivatives_system);
postprocessor->
- compute_derived_quantities_vector(data.postprocessed_values[dataset],
- data.patch_values_system,
+ compute_derived_quantities_vector(data.patch_values_system,
data.patch_gradients_system,
data.patch_second_derivatives_system,
- data.patch_normals);
+ data.patch_normals,
+ data.postprocessed_values[dataset]);
}
for (unsigned int q=0; q<n_q_points; ++q)
this->dof_data[dataset]->get_function_second_derivatives (fe_patch_values,
data.patch_second_derivatives);
postprocessor->
- compute_derived_quantities_scalar(data.postprocessed_values[dataset],
- data.patch_values,
+ compute_derived_quantities_scalar(data.patch_values,
data.patch_gradients,
data.patch_second_derivatives,
- data.dummy_normals);
+ data.dummy_normals,
+ data.postprocessed_values[dataset]);
}
else
{
this->dof_data[dataset]->get_function_second_derivatives (fe_patch_values,
data.patch_second_derivatives_system);
postprocessor->
- compute_derived_quantities_vector(data.postprocessed_values[dataset],
- data.patch_values_system,
+ compute_derived_quantities_vector(data.patch_values_system,
data.patch_gradients_system,
data.patch_second_derivatives_system,
- data.dummy_normals);
+ data.dummy_normals,
+ data.postprocessed_values[dataset]);
}
for (unsigned int component=0;
template <int dim>
void
DataPostprocessor<dim>::
-compute_derived_quantities_scalar (std::vector<Vector<double> > &computed_quantities,
- const std::vector<double> &/*uh*/,
+compute_derived_quantities_scalar (const std::vector<double> &/*uh*/,
const std::vector<Tensor<1,dim> > &/*duh*/,
const std::vector<Tensor<2,dim> > &/*dduh*/,
- const std::vector<Point<dim> > &/*normals*/) const
+ const std::vector<Point<dim> > &/*normals*/,
+ std::vector<Vector<double> > &computed_quantities) const
{
computed_quantities.clear();
AssertThrow(false,ExcPureFunctionCalled());
template <int dim>
void
DataPostprocessor<dim>::
-compute_derived_quantities_vector (std::vector<Vector<double> > &computed_quantities,
- const std::vector<Vector<double> > &/*uh*/,
+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*/,
- const std::vector<Point<dim> > &/*normals*/) const
+ const std::vector<Point<dim> > &/*normals*/,
+ std::vector<Vector<double> > &computed_quantities) const
{
computed_quantities.clear();
AssertThrow(false,ExcPureFunctionCalled());
public:
void compute_derived_quantities_vector (
- std::vector< Vector< double > > &,
const std::vector< Vector< double > > &,
const std::vector< std::vector< Tensor< 1, dim > > > &,
const std::vector< std::vector< Tensor< 2, dim > > > &,
- const std::vector< Point< dim > > &
+ const std::vector< Point< dim > > &,
+ std::vector< Vector< double > > &
) const;
std::vector<std::string> get_names () const;
template <int dim>
void
Postprocessor<dim>::compute_derived_quantities_vector (
- std::vector< Vector< double > > &computed_quantities,
const std::vector< Vector< double > > &uh,
const std::vector< std::vector< Tensor< 1, dim > > > &/*duh*/,
const std::vector< std::vector< Tensor< 2, dim > > > &/*dduh*/,
- const std::vector< Point< dim > > &/*normals*/
+ const std::vector< Point< dim > > &/*normals*/,
+ std::vector< Vector< double > > &computed_quantities
) const
{
Assert(computed_quantities.size() == uh.size(),