Reference documentation for deal.II version GIT c1ddcf411d 2023-11-30 16:30:02+00:00
\(\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\}}\)
integration_info.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2022 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_mesh_worker_integration_info_h
18 #define dealii_mesh_worker_integration_info_h
19 
20 #include <deal.II/base/config.h>
21 
23 
25 
26 #include <deal.II/fe/fe_values.h>
27 
31 
32 #include <memory>
33 
35 
36 namespace MeshWorker
37 {
76  template <int dim, int spacedim = dim>
78  {
79  private:
81  std::vector<std::shared_ptr<FEValuesBase<dim, spacedim>>> fevalv;
82 
83  public:
84  static constexpr unsigned int dimension = dim;
85  static constexpr unsigned int space_dimension = spacedim;
86 
91 
96 
115  template <class FEVALUES>
116  void
118  const Mapping<dim, spacedim> &mapping,
120  const UpdateFlags flags,
121  const BlockInfo *local_block_info = nullptr);
122 
126  void
127  initialize_data(const std::shared_ptr<VectorDataBase<dim, spacedim>> &data);
128 
132  void
133  clear();
134 
140  finite_element() const;
141 
143  bool multigrid;
145 
151  fe_values() const;
152 
154 
159  fe_values(const unsigned int i) const;
160 
169  std::vector<std::vector<std::vector<double>>> values;
170 
179  std::vector<std::vector<std::vector<Tensor<1, spacedim>>>> gradients;
180 
189  std::vector<std::vector<std::vector<Tensor<2, spacedim>>>> hessians;
190 
194  template <typename number>
195  void
197 
202  template <typename number>
203  void
205  bool split_fevalues);
206 
211  std::shared_ptr<VectorDataBase<dim, spacedim>> global_data;
212 
216  std::size_t
218 
219  private:
226 
232  template <typename TYPE>
233  void
234  fill_local_data(std::vector<std::vector<std::vector<TYPE>>> &data,
236  bool split_fevalues) const;
240  unsigned int n_components;
241  };
242 
296  template <int dim, int spacedim = dim>
298  {
299  public:
304 
309 
317  void
319  const Mapping<dim, spacedim> &mapping,
320  const BlockInfo *block_info = nullptr);
321 
329  template <typename VectorType>
330  void
332  const Mapping<dim, spacedim> &mapping,
333  const AnyData &data,
334  const VectorType &dummy,
335  const BlockInfo *block_info = nullptr);
343  template <typename VectorType>
344  void
346  const Mapping<dim, spacedim> &mapping,
347  const AnyData &data,
348  const MGLevelObject<VectorType> &dummy,
349  const BlockInfo *block_info = nullptr);
365  void
366  initialize_update_flags(bool neighbor_geometry = false);
367 
372  void
373  add_update_flags_all(const UpdateFlags flags);
374 
378  void
379  add_update_flags_cell(const UpdateFlags flags);
380 
384  void
386 
390  void
391  add_update_flags_face(const UpdateFlags flags);
392 
399  void
401  const bool cell = true,
402  const bool boundary = true,
403  const bool face = true,
404  const bool neighbor = true);
405 
415  void
416  initialize_gauss_quadrature(unsigned int n_cell_points,
417  unsigned int n_boundary_points,
418  unsigned int n_face_points,
419  const bool force = true);
420 
424  std::size_t
426 
439 
446 
454 
459 
464 
490 
496 
502 
503  std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim>> cell_data;
504  std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim>> boundary_data;
505  std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim>> face_data;
529  template <class DOFINFO>
530  void
532 
549  template <class DOFINFO>
550  void
552 
574 
576  };
577 
578 
579  //----------------------------------------------------------------------//
580 
581  template <int dim, int sdim>
583  : fevalv(0)
584  , multigrid(false)
585  , global_data(std::make_shared<VectorDataBase<dim, sdim>>())
586  , n_components(numbers::invalid_unsigned_int)
587  {}
588 
589 #ifndef DOXYGEN
590  template <int dim, int sdim>
592  const IntegrationInfo<dim, sdim> &other)
593  : multigrid(other.multigrid)
594  , values(other.values)
595  , gradients(other.gradients)
596  , hessians(other.hessians)
597  , global_data(other.global_data)
598  , fe_pointer(other.fe_pointer)
599  , n_components(other.n_components)
600  {
601  fevalv.resize(other.fevalv.size());
602  for (unsigned int i = 0; i < other.fevalv.size(); ++i)
603  {
604  const FEValuesBase<dim, sdim> &p = *other.fevalv[i];
605  const FEValues<dim, sdim> *pc =
606  dynamic_cast<const FEValues<dim, sdim> *>(&p);
607  const FEFaceValues<dim, sdim> *pf =
608  dynamic_cast<const FEFaceValues<dim, sdim> *>(&p);
609  const FESubfaceValues<dim, sdim> *ps =
610  dynamic_cast<const FESubfaceValues<dim, sdim> *>(&p);
611 
612  if (pc != nullptr)
613  fevalv[i] =
614  std::make_shared<FEValues<dim, sdim>>(pc->get_mapping(),
615  pc->get_fe(),
616  pc->get_quadrature(),
617  pc->get_update_flags());
618  else if (pf != nullptr)
619  fevalv[i] =
620  std::make_shared<FEFaceValues<dim, sdim>>(pf->get_mapping(),
621  pf->get_fe(),
622  pf->get_quadrature(),
623  pf->get_update_flags());
624  else if (ps != nullptr)
625  fevalv[i] = std::make_shared<FESubfaceValues<dim, sdim>>(
626  ps->get_mapping(),
627  ps->get_fe(),
628  ps->get_quadrature(),
629  ps->get_update_flags());
630  else
631  Assert(false, ExcInternalError());
632  }
633  }
634 #endif
635 
636 
637 
638  template <int dim, int sdim>
639  template <class FEVALUES>
640  inline void
642  const FiniteElement<dim, sdim> &el,
643  const Mapping<dim, sdim> &mapping,
645  const UpdateFlags flags,
646  const BlockInfo *block_info)
647  {
648  fe_pointer = &el;
649  if (block_info == nullptr || block_info->local().size() == 0)
650  {
651  fevalv.resize(1);
652  fevalv[0] = std::make_shared<FEVALUES>(mapping, el, quadrature, flags);
653  }
654  else
655  {
656  fevalv.resize(el.n_base_elements());
657  for (unsigned int i = 0; i < fevalv.size(); ++i)
658  fevalv[i] = std::make_shared<FEVALUES>(mapping,
659  el.base_element(i),
660  quadrature,
661  flags);
662  }
663  n_components = el.n_components();
664  }
665 
666 
667  template <int dim, int spacedim>
668  inline const FiniteElement<dim, spacedim> &
670  {
671  Assert(fe_pointer != nullptr, ExcNotInitialized());
672  return *fe_pointer;
673  }
674 
675  template <int dim, int spacedim>
676  inline const FEValuesBase<dim, spacedim> &
678  {
679  AssertDimension(fevalv.size(), 1);
680  return *fevalv[0];
681  }
682 
683 
684  template <int dim, int spacedim>
685  inline const FEValuesBase<dim, spacedim> &
687  {
688  AssertIndexRange(i, fevalv.size());
689  return *fevalv[i];
690  }
691 
692 
693  template <int dim, int spacedim>
694  template <typename number>
695  inline void
697  const DoFInfo<dim, spacedim, number> &info)
698  {
699  for (unsigned int i = 0; i < fevalv.size(); ++i)
700  {
701  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
703  {
704  // This is a subface
706  dynamic_cast<FESubfaceValues<dim, spacedim> &>(febase);
707  fe.reinit(info.cell, info.face_number, info.sub_number);
708  }
710  {
711  // This is a face
713  dynamic_cast<FEFaceValues<dim, spacedim> &>(febase);
714  fe.reinit(info.cell, info.face_number);
715  }
716  else
717  {
718  // This is a cell
720  dynamic_cast<FEValues<dim, spacedim> &>(febase);
721  fe.reinit(info.cell);
722  }
723  }
724 
725  const bool split_fevalues = info.block_info != nullptr;
726  if (!global_data->empty())
727  fill_local_data(info, split_fevalues);
728  }
729 
730 
731 
732  //----------------------------------------------------------------------//
733 
734  template <int dim, int sdim>
735  inline void
737  unsigned int bp,
738  unsigned int fp,
739  bool force)
740  {
741  if (force || cell_quadrature.empty())
742  cell_quadrature = QGauss<dim>(cp);
743  if (force || boundary_quadrature.empty())
744  boundary_quadrature = QGauss<dim - 1>(bp);
745  if (force || face_quadrature.empty())
746  face_quadrature = QGauss<dim - 1>(fp);
747  }
748 
749 
750  template <int dim, int sdim>
751  inline void
753  {
754  add_update_flags(flags, true, true, true, true);
755  }
756 
757 
758  template <int dim, int sdim>
759  inline void
761  {
762  add_update_flags(flags, true, false, false, false);
763  }
764 
765 
766  template <int dim, int sdim>
767  inline void
769  const UpdateFlags flags)
770  {
771  add_update_flags(flags, false, true, false, false);
772  }
773 
774 
775  template <int dim, int sdim>
776  inline void
778  {
779  add_update_flags(flags, false, false, true, true);
780  }
781 
782 
783 #ifndef DOXYGEN
784  template <int dim, int sdim>
785  inline void
787  const Mapping<dim, sdim> &mapping,
788  const BlockInfo *block_info)
789  {
790  initialize_update_flags();
791  initialize_gauss_quadrature((cell_flags & update_values) ?
792  (el.tensor_degree() + 1) :
793  el.tensor_degree(),
794  (boundary_flags & update_values) ?
795  (el.tensor_degree() + 1) :
796  el.tensor_degree(),
797  (face_flags & update_values) ?
798  (el.tensor_degree() + 1) :
799  el.tensor_degree(),
800  false);
801 
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);
812  }
813 #endif
814 
815 
816  template <int dim, int sdim>
817  template <typename VectorType>
818  void
820  const Mapping<dim, sdim> &mapping,
821  const AnyData &data,
822  const VectorType &,
823  const BlockInfo *block_info)
824  {
825  initialize(el, mapping, block_info);
826  std::shared_ptr<VectorData<VectorType, dim, sdim>> p;
828 
829  p = std::make_shared<VectorData<VectorType, dim, sdim>>(cell_selector);
830  // Public member function of parent class was not found without
831  // explicit cast
832  pp = &*p;
833  pp->initialize(data);
834  cell_data = p;
835  cell.initialize_data(p);
836 
837  p = std::make_shared<VectorData<VectorType, dim, sdim>>(boundary_selector);
838  pp = &*p;
839  pp->initialize(data);
840  boundary_data = p;
841  boundary.initialize_data(p);
842 
843  p = std::make_shared<VectorData<VectorType, dim, sdim>>(face_selector);
844  pp = &*p;
845  pp->initialize(data);
846  face_data = p;
847  face.initialize_data(p);
848  subface.initialize_data(p);
849  neighbor.initialize_data(p);
850  }
851 
852  template <int dim, int sdim>
853  template <typename VectorType>
854  void
856  const Mapping<dim, sdim> &mapping,
857  const AnyData &data,
859  const BlockInfo *block_info)
860  {
861  initialize(el, mapping, block_info);
862  std::shared_ptr<MGVectorData<VectorType, dim, sdim>> p;
864 
865  p = std::make_shared<MGVectorData<VectorType, dim, sdim>>(cell_selector);
866  // Public member function of parent class was not found without
867  // explicit cast
868  pp = &*p;
869  pp->initialize(data);
870  cell_data = p;
871  cell.initialize_data(p);
872 
873  p =
874  std::make_shared<MGVectorData<VectorType, dim, sdim>>(boundary_selector);
875  pp = &*p;
876  pp->initialize(data);
877  boundary_data = p;
878  boundary.initialize_data(p);
879 
880  p = std::make_shared<MGVectorData<VectorType, dim, sdim>>(face_selector);
881  pp = &*p;
882  pp->initialize(data);
883  face_data = p;
884  face.initialize_data(p);
885  subface.initialize_data(p);
886  neighbor.initialize_data(p);
887  }
888 
889  template <int dim, int sdim>
890  template <class DOFINFO>
891  void
893  {}
894 
895 
896  template <int dim, int sdim>
897  template <class DOFINFO>
898  void
900  {}
901 
902 
903 } // namespace MeshWorker
904 
906 
907 #endif
unsigned int size() const
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:96
const BlockIndices & local() const
Definition: block_info.h:229
const Quadrature< dim - 1 > & get_quadrature() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const unsigned int face_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const unsigned int face_no, const unsigned int subface_no)
UpdateFlags get_update_flags() const
const Mapping< dim, spacedim > & get_mapping() const
const FiniteElement< dim, spacedim > & get_fe() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell)
const Quadrature< dim > & get_quadrature() const
unsigned int tensor_degree() const
unsigned int n_components() const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
unsigned int n_base_elements() const
unsigned int sub_number
Definition: dof_info.h:98
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:79
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:172
unsigned int face_number
Definition: dof_info.h:90
std::size_t memory_consumption() const
void initialize_update_flags(bool neighbor_geometry=false)
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 > &)
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const BlockInfo *block_info=nullptr)
void add_update_flags_all(const UpdateFlags flags)
MeshWorker::VectorSelector boundary_selector
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< MeshWorker::VectorDataBase< dim, spacedim > > face_data
std::shared_ptr< MeshWorker::VectorDataBase< dim, spacedim > > boundary_data
std::shared_ptr< MeshWorker::VectorDataBase< dim, spacedim > > cell_data
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const AnyData &data, const MGLevelObject< VectorType > &dummy, const BlockInfo *block_info=nullptr)
Quadrature< dim - 1 > boundary_quadrature
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const AnyData &data, const VectorType &dummy, const BlockInfo *block_info=nullptr)
void add_update_flags_cell(const UpdateFlags flags)
Quadrature< dim - 1 > face_quadrature
MeshWorker::VectorSelector face_selector
void post_faces(const DoFInfoBox< dim, DOFINFO > &)
MeshWorker::VectorSelector cell_selector
void add_update_flags_face(const UpdateFlags flags)
void add_update_flags_boundary(const UpdateFlags flags)
static constexpr unsigned int space_dimension
static constexpr unsigned int dimension
std::vector< std::vector< std::vector< Tensor< 2, spacedim > > > > hessians
std::size_t memory_consumption() const
const FEValuesBase< dim, spacedim > & fe_values(const unsigned int i) const
Access to finite elements.
const FiniteElement< dim, spacedim > & finite_element() const
IntegrationInfo(const IntegrationInfo< dim, spacedim > &other)
const FEValuesBase< dim, spacedim > & fe_values() const
Access to finite element.
bool multigrid
This is true if we are assembling for multigrid.
std::shared_ptr< VectorDataBase< dim, spacedim > > global_data
void reinit(const DoFInfo< dim, spacedim, number > &i)
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)
void initialize_data(const std::shared_ptr< VectorDataBase< dim, spacedim >> &data)
std::vector< std::vector< std::vector< Tensor< 1, spacedim > > > > gradients
void fill_local_data(std::vector< std::vector< std::vector< TYPE >>> &data, VectorSelector &selector, bool split_fevalues) const
std::vector< std::shared_ptr< FEValuesBase< dim, spacedim > > > fevalv
vector of FEValues objects
SmartPointer< const FiniteElement< dim, spacedim >, IntegrationInfo< dim, spacedim > > fe_pointer
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)
std::vector< std::vector< std::vector< double > > > values
void initialize(const AnyData &)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
#define Assert(cond, exc)
Definition: exceptions.h:1631
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1820
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
UpdateFlags
@ update_values
Shape function values.
std::vector< unsigned int > selector(const std::vector< unsigned int > &targets, const std::function< RequestType(const unsigned int)> &create_request, const std::function< AnswerType(const unsigned int, const RequestType &)> &answer_request, const std::function< void(const unsigned int, const AnswerType &)> &process_answer, const MPI_Comm comm)
static const unsigned int invalid_unsigned_int
Definition: types.h:221