const dealii::hp::QCollection<dim> q_collection;
const dealii::hp::MappingCollection<dim,spacedim> mapping_collection;
dealii::hp::FEValues<dim,spacedim> x_fe_values;
+
+ std::vector<Point<spacedim> > patch_evaluation_points;
const std::vector<std::vector<unsigned int> > *cell_to_patch_index_map;
};
dealii::hp::FEFaceValues<dim> x_fe_values;
std::vector<Point<dim> > patch_normals;
+ std::vector<Point<spacedim> > patch_evaluation_points;
};
}
}
const dealii::hp::QCollection<dim> q_collection;
dealii::hp::FEValues<dim,spacedim> x_fe_values;
+
+ std::vector<Point<spacedim> > patch_evaluation_points;
};
}
}
// $Id$
// Version: $Name$
//
-// Copyright (C) 2007, 2008 by the deal.II authors
+// Copyright (C) 2007, 2008, 2010 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
virtual ~DataPostprocessor ();
+ /**
+ * @deprecated
+ *
+ * This function only exists for backward
+ * compatibility as this is the interface
+ * provided by previous versions of the
+ * library. The default implementation of
+ * the other function of same name below
+ * calls this function by simply dropping
+ * the argument that denotes the
+ * evaluation points. Since this function
+ * might at one point go away, you should
+ * overload the other function instead.
+ */
+ 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,
+ const std::vector<Point<dim> > &normals,
+ std::vector<Vector<double> > &computed_quantities) const;
+
/**
* This is the main function which actually
* performs the postprocessing. The last
const std::vector<Tensor<1,dim> > &duh,
const std::vector<Tensor<2,dim> > &dduh,
const std::vector<Point<dim> > &normals,
- std::vector<Vector<double> > &computed_quantities) const;
+ const std::vector<Point<dim> > &evaluation_points,
+ std::vector<Vector<double> > &computed_quantities) const;
/**
- * Same as above function, but
- * this function is called when
- * the original data vector
- * represents vector data,
- * i.e. the finite element in use
- * has multiple vector
- * components.
+ * @deprecated
+ *
+ * This function only exists for backward
+ * compatibility as this is the interface
+ * provided by previous versions of the
+ * library. The default implementation of
+ * the other function of same name below
+ * calls this function by simply dropping
+ * the argument that denotes the
+ * evaluation points. Since this function
+ * might at one point go away, you should
+ * overload the other function instead.
*/
virtual
void
const std::vector<std::vector<Tensor<2,dim> > > &dduh,
const std::vector<Point<dim> > &normals,
std::vector<Vector<double> > &computed_quantities) const;
+
+ /**
+ * Same as the
+ * compute_derived_quantities_scalar()
+ * function, but this function is called
+ * when the original data vector
+ * represents vector data, i.e. the
+ * finite element in use has multiple
+ * vector components.
+ */
+ 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,
+ const std::vector<Point<dim> > &normals,
+ const std::vector<Point<dim> > &evaluation_points,
+ std::vector<Vector<double> > &computed_quantities) const;
/**
* Return the vector of strings describing
if (update_flags & update_hessians)
this->dof_data[dataset]->get_function_hessians (fe_patch_values,
data.patch_hessians);
+
+ if (update_flags & update_quadrature_points)
+ data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+
+
std::vector<Point<DH::space_dimension> > dummy_normals;
postprocessor->
compute_derived_quantities_scalar(data.patch_values,
data.patch_gradients,
data.patch_hessians,
dummy_normals,
+ data.patch_evaluation_points,
data.postprocessed_values[dataset]);
}
else
if (update_flags & update_hessians)
this->dof_data[dataset]->get_function_hessians (fe_patch_values,
data.patch_hessians_system);
+
+ if (update_flags & update_quadrature_points)
+ data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+
std::vector<Point<DH::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]);
}
if (update_flags & update_hessians)
this->dof_data[dataset]->get_function_hessians (fe_patch_values,
data.patch_hessians);
+
+ if (update_flags & update_quadrature_points)
+ data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+
postprocessor->
compute_derived_quantities_scalar(data.patch_values,
data.patch_gradients,
data.patch_hessians,
data.patch_normals,
+ data.patch_evaluation_points,
data.postprocessed_values[dataset]);
}
else
if (update_flags & update_hessians)
this->dof_data[dataset]->get_function_hessians (fe_patch_values,
data.patch_hessians_system);
+
+ if (update_flags & update_quadrature_points)
+ data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+
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]);
}
if (update_flags & update_hessians)
this->dof_data[dataset]->get_function_hessians (fe_patch_values,
data.patch_hessians);
+
+ if (update_flags & update_quadrature_points)
+ data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
std::vector<Point<DH::space_dimension> > dummy_normals;
postprocessor->
data.patch_gradients,
data.patch_hessians,
dummy_normals,
+ data.patch_evaluation_points,
data.postprocessed_values[dataset]);
}
else
if (update_flags & update_hessians)
this->dof_data[dataset]->get_function_hessians (fe_patch_values,
data.patch_hessians_system);
+
+ if (update_flags & update_quadrature_points)
+ data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
+
std::vector<Point<DH::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]);
}
// $Id$
// Version: $Name$
//
-// Copyright (C) 2007 by the deal.II authors
+// Copyright (C) 2007, 2010 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
}
+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,
+ const std::vector<Point<dim> > &normals,
+ const std::vector<Point<dim> > &/*evaluation_points*/,
+ std::vector<Vector<double> > &computed_quantities) const
+{
+ compute_derived_quantities_scalar(uh, duh, dduh, normals, computed_quantities);
+}
+
+
template <int dim>
void
DataPostprocessor<dim>::
-compute_derived_quantities_vector (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*/,
+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,
+ const std::vector<Point<dim> > &normals,
+ const std::vector<Point<dim> > &/*evaluation_points*/,
+ std::vector<Vector<double> > &computed_quantities) const
+{
+ compute_derived_quantities_vector(uh, duh, dduh, normals, computed_quantities);
+}
+
+
+
template <int dim>
std::vector<DataComponentInterpretation::DataComponentInterpretation>
DataPostprocessor<dim>::get_data_component_interpretation () const
<h3>deal.II</h3>
<ol>
+ <li><p>Changed: The DataPostprocessor functions now take an additional
+ argument that indicates the location of an evaluation point. For backward
+ compatibility, the old function signature still exists so that applications
+ that overload one of the existing functions continue to work.
+ The old signature has been deprecated, however, and will be removed in a
+ future version.
+ <br>
+ (Scott Miller 2010/10/08)
+ </p></li>
<li><p>Changed: FETools is now a namespace rather than a class with only
static member functions.
/* Author: Moritz Allmaras, Texas A&M University, 2007 */
/* */
-/* Copyright (C) 2007, 2008 by the deal.II authors and M. Allmaras */
+/* Copyright (C) 2007, 2008, 2010 by the deal.II authors and M. Allmaras */
/* */
/* This file is subject to QPL and may not be distributed */
/* without copyright and license information. Please refer */
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;
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> > & /*evaluation_points*/,
std::vector< Vector< double > > & computed_quantities
) const
{
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> > &evaluation_points,
std::vector<Vector<double> > &computed_quantities) const;
virtual std::vector<std::string> get_names () const;
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> > &/*evaluation_points*/,
std::vector<Vector<double> > &computed_quantities) const
{
// At the beginning of the function, let us