|
Reference documentation for deal.II version 9.2.0
|
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\)
\(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\)
\(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\)
\(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Go to the documentation of this file.
17 #ifndef dealii_mesh_worker_integration_info_h
18 #define dealii_mesh_worker_integration_info_h
77 template <
int dim,
int spacedim = dim>
82 std::vector<std::shared_ptr<FEValuesBase<dim, spacedim>>>
fevalv;
116 template <
class FEVALUES>
122 const BlockInfo *local_block_info =
nullptr);
170 std::vector<std::vector<std::vector<double>>>
values;
180 std::vector<std::vector<std::vector<Tensor<1, spacedim>>>>
gradients;
190 std::vector<std::vector<std::vector<Tensor<2, spacedim>>>>
hessians;
195 template <
typename number>
203 template <
typename number>
206 bool split_fevalues);
233 template <
typename TYPE>
237 bool split_fevalues)
const;
298 template <
int dim,
int spacedim = dim>
331 template <
typename VectorType>
345 template <
typename VectorType>
403 const bool cell =
true,
405 const bool face =
true,
419 unsigned int n_boundary_points,
420 unsigned int n_face_points,
421 const bool force =
true);
505 std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim>>
cell_data;
507 std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim>>
face_data;
531 template <
class DOFINFO>
551 template <
class DOFINFO>
583 template <
int dim,
int sdim>
592 template <
int dim,
int sdim>
595 : multigrid(other.multigrid)
596 , values(other.values)
597 , gradients(other.gradients)
598 , hessians(other.hessians)
599 , global_data(other.global_data)
600 , fe_pointer(other.fe_pointer)
603 fevalv.resize(other.
fevalv.size());
604 for (
unsigned int i = 0; i < other.
fevalv.size(); ++i)
616 std::make_shared<FEValues<dim, sdim>>(pc->get_mapping(),
619 pc->get_update_flags());
620 else if (pf !=
nullptr)
622 std::make_shared<FEFaceValues<dim, sdim>>(pf->get_mapping(),
625 pf->get_update_flags());
626 else if (ps !=
nullptr)
627 fevalv[i] = std::make_shared<FESubfaceValues<dim, sdim>>(
630 ps->get_quadrature(),
631 ps->get_update_flags());
639 template <
int dim,
int sdim>
640 template <
class FEVALUES>
650 if (block_info ==
nullptr || block_info->
local().
size() == 0)
653 fevalv[0] = std::make_shared<FEVALUES>(mapping, el, quadrature, flags);
658 for (
unsigned int i = 0; i < fevalv.size(); ++i)
659 fevalv[i] = std::make_shared<FEVALUES>(mapping,
668 template <
int dim,
int spacedim>
676 template <
int dim,
int spacedim>
685 template <
int dim,
int spacedim>
694 template <
int dim,
int spacedim>
695 template <
typename number>
700 for (
unsigned int i = 0; i < fevalv.size(); ++i)
726 const bool split_fevalues = info.
block_info !=
nullptr;
727 if (!global_data->empty())
728 fill_local_data(info, split_fevalues);
735 template <
int dim,
int sdim>
742 if (force || cell_quadrature.size() == 0)
744 if (force || boundary_quadrature.size() == 0)
745 boundary_quadrature =
QGauss<dim - 1>(bp);
746 if (force || face_quadrature.size() == 0)
751 template <
int dim,
int sdim>
755 add_update_flags(flags,
true,
true,
true,
true);
759 template <
int dim,
int sdim>
763 add_update_flags(flags,
true,
false,
false,
false);
767 template <
int dim,
int sdim>
772 add_update_flags(flags,
false,
true,
false,
false);
776 template <
int dim,
int sdim>
780 add_update_flags(flags,
false,
false,
true,
true);
784 template <
int dim,
int sdim>
790 initialize_update_flags();
802 cell.template initialize<FEValues<dim, sdim>>(
803 el, mapping, cell_quadrature, cell_flags, block_info);
804 boundary.template initialize<FEFaceValues<dim, sdim>>(
805 el, mapping, boundary_quadrature, boundary_flags, block_info);
806 face.template initialize<FEFaceValues<dim, sdim>>(
807 el, mapping, face_quadrature, face_flags, block_info);
808 subface.template initialize<FESubfaceValues<dim, sdim>>(
809 el, mapping, face_quadrature, face_flags, block_info);
810 neighbor.template initialize<FEFaceValues<dim, sdim>>(
811 el, mapping, face_quadrature, neighbor_flags, block_info);
815 template <
int dim,
int sdim>
816 template <
typename VectorType>
824 initialize(el, mapping, block_info);
825 std::shared_ptr<VectorData<VectorType, dim, sdim>> p;
828 p = std::make_shared<VectorData<VectorType, dim, sdim>>(cell_selector);
832 pp->initialize(data);
834 cell.initialize_data(p);
836 p = std::make_shared<VectorData<VectorType, dim, sdim>>(boundary_selector);
838 pp->initialize(data);
840 boundary.initialize_data(p);
842 p = std::make_shared<VectorData<VectorType, dim, sdim>>(face_selector);
844 pp->initialize(data);
846 face.initialize_data(p);
847 subface.initialize_data(p);
848 neighbor.initialize_data(p);
851 template <
int dim,
int sdim>
852 template <
typename VectorType>
860 initialize(el, mapping, block_info);
861 std::shared_ptr<MGVectorData<VectorType, dim, sdim>> p;
864 p = std::make_shared<MGVectorData<VectorType, dim, sdim>>(cell_selector);
868 pp->initialize(data);
870 cell.initialize_data(p);
873 std::make_shared<MGVectorData<VectorType, dim, sdim>>(boundary_selector);
875 pp->initialize(data);
877 boundary.initialize_data(p);
879 p = std::make_shared<MGVectorData<VectorType, dim, sdim>>(face_selector);
881 pp->initialize(data);
883 face.initialize_data(p);
884 subface.initialize_data(p);
885 neighbor.initialize_data(p);
888 template <
int dim,
int sdim>
889 template <
class DOFINFO>
895 template <
int dim,
int sdim>
896 template <
class DOFINFO>
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const BlockInfo *block_info=nullptr)
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)
unsigned int size() const
void initialize_gauss_quadrature(unsigned int n_cell_points, unsigned int n_boundary_points, unsigned int n_face_points, const bool force=true)
void post_cell(const DoFInfoBox< dim, DOFINFO > &)
unsigned int tensor_degree() const
void add_update_flags_boundary(const UpdateFlags flags)
void post_faces(const DoFInfoBox< dim, DOFINFO > &)
std::vector< std::vector< std::vector< double > > > values
const FEValuesBase< dim, spacedim > & fe_values() const
Access to finite element.
void add_update_flags(const UpdateFlags flags, const bool cell=true, const bool boundary=true, const bool face=true, const bool neighbor=true)
std::shared_ptr< VectorDataBase< dim, spacedim > > global_data
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
#define AssertIndexRange(index, range)
unsigned int n_base_elements() const
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
std::vector< std::shared_ptr< FEValuesBase< dim, spacedim > > > fevalv
vector of FEValues objects
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no, const unsigned int subface_no)
@ update_values
Shape function values.
const FiniteElement< dim, spacedim > & finite_element() const
void initialize_update_flags(bool neighbor_geometry=false)
const BlockIndices & local() const
std::vector< std::vector< std::vector< Tensor< 1, spacedim > > > > gradients
MeshWorker::VectorSelector face_selector
void reinit(const DoFInfo< dim, spacedim, number > &i)
Abstract base class for mapping classes.
void initialize_data(const std::shared_ptr< VectorDataBase< dim, spacedim >> &data)
std::size_t memory_consumption() const
Quadrature< dim > cell_quadrature
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
#define DEAL_II_NAMESPACE_OPEN
void add_update_flags_face(const UpdateFlags flags)
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell)
Quadrature< dim - 1 > face_quadrature
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcNotInitialized()
Quadrature< dim - 1 > boundary_quadrature
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no)
std::shared_ptr< MeshWorker::VectorDataBase< dim, spacedim > > boundary_data
const types::global_dof_index * dummy()
static const unsigned int dimension
static ::ExceptionBase & ExcInternalError()
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=nullptr)
bool multigrid
This is true if we are assembling for multigrid.
#define Assert(cond, exc)
SmartPointer< const FiniteElement< dim, spacedim >, IntegrationInfo< dim, spacedim > > fe_pointer
void add_update_flags_all(const UpdateFlags flags)
static const unsigned int invalid_unsigned_int
unsigned int n_components
std::shared_ptr< MeshWorker::VectorDataBase< dim, spacedim > > cell_data
unsigned int n_components() const
MeshWorker::VectorSelector cell_selector
std::size_t memory_consumption() const
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
MeshWorker::VectorSelector boundary_selector
#define DEAL_II_NAMESPACE_CLOSE
std::vector< std::vector< std::vector< Tensor< 2, spacedim > > > > hessians
UpdateFlags boundary_flags
const Quadrature< dim > & get_quadrature() const
static const unsigned int space_dimension
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
void add_update_flags_cell(const UpdateFlags flags)
UpdateFlags neighbor_flags
std::shared_ptr< MeshWorker::VectorDataBase< dim, spacedim > > face_data