18 #ifndef __deal2__mesh_worker_integration_info_h 19 #define __deal2__mesh_worker_integration_info_h 21 #include <deal.II/base/config.h> 22 #include <deal.II/base/quadrature_lib.h> 23 #include <deal.II/base/std_cxx1x/shared_ptr.h> 24 #include <deal.II/dofs/block_info.h> 25 #include <deal.II/fe/fe_values.h> 26 #include <deal.II/meshworker/local_results.h> 27 #include <deal.II/meshworker/dof_info.h> 28 #include <deal.II/meshworker/vector_selector.h> 77 template<
int dim,
int spacedim = dim>
82 std::vector<std_cxx1x::shared_ptr<FEValuesBase<dim, spacedim> > >
fevalv;
84 static const unsigned int dimension = dim;
85 static const unsigned int space_dimension = spacedim;
118 template <
class FEVALUES>
162 std::vector<std::vector<std::vector<double> > >
values;
172 std::vector<std::vector<std::vector<Tensor<1,dim> > > >
gradients;
182 std::vector<std::vector<std::vector<Tensor<2,dim> > > >
hessians;
187 template <
typename number>
194 template<
typename number>
201 std_cxx1x::shared_ptr<VectorDataBase<dim, spacedim> >
global_data;
214 template <
typename TYPE>
216 std::vector<std::vector<std::vector<TYPE> > > &data,
218 bool split_fevalues)
const;
281 template <
int dim,
int spacedim=dim>
315 template <
typename VECTOR>
328 template <
typename VECTOR>
348 void initialize_update_flags(
bool neighbor_geometry =
false);
354 void add_update_flags_all (
const UpdateFlags flags);
359 void add_update_flags_cell(
const UpdateFlags flags);
364 void add_update_flags_boundary(
const UpdateFlags flags);
369 void add_update_flags_face(
const UpdateFlags flags);
378 const bool cell =
true,
379 const bool boundary =
true,
380 const bool face =
true,
381 const bool neighbor =
true);
392 void initialize_gauss_quadrature(
unsigned int n_cell_points,
393 unsigned int n_boundary_points,
394 unsigned int n_face_points,
395 const bool force =
true);
479 std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
480 std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
481 std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
506 template <
class DOFINFO>
526 template <
class DOFINFO>
558 template<
int dim,
int sdim>
568 template<
int dim,
int sdim>
579 fevalv.resize(other.fevalv.size());
580 for (
unsigned int i=0; i<other.fevalv.size(); ++i)
588 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
592 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
595 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
604 template<
int dim,
int sdim>
605 template <
class FEVALUES>
614 if (block_info == 0 || block_info->
local().
size() == 0)
617 fevalv[0] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
618 new FEVALUES (mapping, el, quadrature, flags));
623 for (
unsigned int i=0; i<
fevalv.size(); ++i)
625 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
626 new FEVALUES (mapping, el.
base_element(i), quadrature, flags));
633 template <
int dim,
int spacedim>
642 template <
int dim,
int spacedim>
651 template <
int dim,
int spacedim>
652 template <
typename number>
656 for (
unsigned int i=0; i<
fevalv.size(); ++i)
659 if (info.
sub_number != deal_II_numbers::invalid_unsigned_int)
665 else if (info.
face_number != deal_II_numbers::invalid_unsigned_int)
679 const bool split_fevalues = info.
block_info != 0;
689 template <
int dim,
int sdim>
698 if (force || cell_quadrature.size() == 0)
700 if (force || boundary_quadrature.size() == 0)
702 if (force || face_quadrature.size() == 0)
707 template <
int dim,
int sdim>
712 add_update_flags(flags,
true,
true,
true,
true);
716 template <
int dim,
int sdim>
721 add_update_flags(flags,
true,
false,
false,
false);
725 template <
int dim,
int sdim>
730 add_update_flags(flags,
false,
true,
false,
false);
734 template <
int dim,
int sdim>
739 add_update_flags(flags,
false,
false,
true,
true);
743 template <
int dim,
int sdim>
751 initialize_update_flags();
752 initialize_gauss_quadrature(
757 cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
758 cell_flags, block_info);
759 boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
760 boundary_flags, block_info);
761 face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
762 face_flags, block_info);
763 subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
764 face_flags, block_info);
765 neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
766 neighbor_flags, block_info);
770 template <
int dim,
int sdim>
771 template <
typename VECTOR>
780 std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
785 cell.initialize_data(p);
790 boundary.initialize_data(p);
795 face.initialize_data(p);
796 subface.initialize_data(p);
797 neighbor.initialize_data(p);
801 template <
int dim,
int sdim>
802 template <
typename VECTOR>
811 std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
816 cell.initialize_data(p);
821 boundary.initialize_data(p);
826 face.initialize_data(p);
827 subface.initialize_data(p);
828 neighbor.initialize_data(p);
832 template <
int dim,
int sdim>
833 template <
class DOFINFO>
839 template <
int dim,
int sdim>
840 template <
class DOFINFO>
848 DEAL_II_NAMESPACE_CLOSE
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell, const unsigned int face_no)
std::vector< std::vector< std::vector< double > > > values
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
#define AssertDimension(dim1, dim2)
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const BlockInfo *block_info=0)
void post_faces(const DoFInfoBox< dim, DOFINFO > &)
const Quadrature< dim-1 > & get_quadrature() const
Quadrature< dim-1 > boundary_quadrature
UpdateFlags boundary_flags
void add_update_flags_cell(const UpdateFlags flags)
std::vector< std::vector< std::vector< Tensor< 1, dim > > > > gradients
void add_update_flags_all(const UpdateFlags flags)
std::size_t memory_consumption() const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const Quadrature< FEVALUES::integral_dimension > &quadrature, const UpdateFlags flags, const BlockInfo *local_block_info=0)
void add_update_flags_face(const UpdateFlags flags)
UpdateFlags get_update_flags() const
const Quadrature< dim > & get_quadrature() const
unsigned int tensor_degree() const
void reinit(const DoFInfo< dim, spacedim, number > &i)
void initialize_gauss_quadrature(unsigned int n_cell_points, unsigned int n_boundary_points, unsigned int n_face_points, const bool force=true)
Quadrature< dim-1 > face_quadrature
MeshWorker::VectorSelector cell_selector
#define Assert(cond, exc)
void post_cell(const DoFInfoBox< dim, DOFINFO > &)
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell)
Quadrature< dim > cell_quadrature
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
UpdateFlags neighbor_flags
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
unsigned int n_base_elements() const
unsigned int n_components() const
void initialize_data(const std_cxx1x::shared_ptr< VectorDataBase< dim, spacedim > > &data)
const Mapping< dim, spacedim > & get_mapping() const
void add_update_flags_boundary(const UpdateFlags flags)
const FEValuesBase< dim, spacedim > & fe_values() const
Access to finite element.
MeshWorker::VectorSelector boundary_selector
bool multigrid
This is true if we are assembling for multigrid.
const BlockIndices & local() const
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
unsigned int size() const
const FiniteElement< dim, spacedim > & get_fe() const
IntegrationInfo< dim, spacedim > CellInfo
std::vector< std::vector< std::vector< Tensor< 2, dim > > > > hessians
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell, const unsigned int face_no, const unsigned int subface_no)
::ExceptionBase & ExcInternalError()
unsigned int n_components
std::vector< std_cxx1x::shared_ptr< FEValuesBase< dim, spacedim > > > fevalv
vector of FEValues objects
MeshWorker::VectorSelector face_selector
std_cxx1x::shared_ptr< VectorDataBase< dim, spacedim > > global_data
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)