17 #ifndef dealii__matrix_free_mapping_data_on_the_fly_h 18 #define dealii__matrix_free_mapping_data_on_the_fly_h 21 #include <deal.II/base/config.h> 22 #include <deal.II/base/exceptions.h> 23 #include <deal.II/base/subscriptor.h> 24 #include <deal.II/base/vectorization.h> 25 #include <deal.II/base/aligned_vector.h> 26 #include <deal.II/matrix_free/shape_info.h> 27 #include <deal.II/matrix_free/mapping_info.h> 28 #include <deal.II/fe/fe_values.h> 29 #include <deal.II/fe/fe_nothing.h> 32 DEAL_II_NAMESPACE_OPEN
37 namespace MatrixFreeFunctions
57 template <
int dim,
typename Number=
double>
82 void reinit(typename ::Triangulation<dim>::cell_iterator cell);
93 typename ::Triangulation<dim>::cell_iterator
get_cell ()
const;
188 template <
int dim,
typename Number>
195 internal::MatrixFreeFunctions::
MappingInfo<dim,Number>::compute_update_flags(update_flags)),
208 template <
int dim,
typename Number>
214 internal::MatrixFreeFunctions::
MappingInfo<dim,Number>::compute_update_flags(update_flags)),
227 template <
int dim,
typename Number>
239 for (
unsigned int d=0; d<dim; ++d)
240 for (
unsigned int e=0; e<dim; ++e)
243 for (
unsigned int d=0; d<dim; ++d)
246 for (
unsigned int d=0; d<dim; ++d)
255 template <
int dim,
typename Number>
260 return present_cell != typename ::Triangulation<dim>::cell_iterator();
265 template <
int dim,
typename Number>
267 typename ::Triangulation<dim>::cell_iterator
275 template <
int dim,
typename Number>
277 const ::FEValues<dim> &
285 template <
int dim,
typename Number>
295 template <
int dim,
typename Number>
305 template <
int dim,
typename Number>
315 template <
int dim,
typename Number>
325 template <
int dim,
typename Number>
338 DEAL_II_NAMESPACE_CLOSE
Transformed quadrature weights.
const ::FEValues< dim > & get_fe_values() const
const AlignedVector< Tensor< 2, dim, VectorizedArray< Number > > > & get_inverse_jacobians() const
FE_Nothing< dim > fe_dummy
typename::Triangulation< dim >::cell_iterator present_cell
const Quadrature< 1 > & get_quadrature() const
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
const AlignedVector< VectorizedArray< Number > > & get_JxW_values() const
Transformed quadrature points.
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
::FEValues< dim > fe_values
UpdateFlags get_update_flags() const
const Quadrature< dim > & get_quadrature() const
void reinit(typename::Triangulation< dim >::cell_iterator cell)
const AlignedVector< Tensor< 1, dim, VectorizedArray< Number > > > & get_normal_vectors() const
AlignedVector< Point< dim, VectorizedArray< Number > > > quadrature_points
AlignedVector< Tensor< 2, dim, VectorizedArray< Number > > > inverse_jacobians
#define Assert(cond, exc)
AlignedVector< Tensor< 1, dim, VectorizedArray< Number > > > normal_vectors
const Quadrature< 1 > quadrature_1d
AlignedVector< VectorizedArray< Number > > jxw_values
typename::Triangulation< dim >::cell_iterator get_cell() const
Gradient of volume element.
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
const AlignedVector< Point< dim, VectorizedArray< Number > > > & get_quadrature_points() const
bool is_initialized() const
static::ExceptionBase & ExcNotImplemented()
const Point< spacedim > & quadrature_point(const unsigned int q) const
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell)
double JxW(const unsigned int quadrature_point) const
MappingDataOnTheFly(const Mapping< dim > &mapping, const Quadrature< 1 > &quadrature, const UpdateFlags update_flags)