Reference documentation for deal.II version Git a670ecb 2017-04-28 20:25:58 -0400
integration_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2016 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 at
12 // the top level of the deal.II distribution.
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 #include <deal.II/base/quadrature_lib.h>
22 #include <deal.II/dofs/block_info.h>
23 #include <deal.II/fe/fe_values.h>
24 #include <deal.II/meshworker/local_results.h>
25 #include <deal.II/meshworker/dof_info.h>
26 #include <deal.II/meshworker/vector_selector.h>
27 
28 #include <memory>
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 namespace MeshWorker
33 {
73  template<int dim, int spacedim = dim>
75  {
76  private:
78  std::vector<std::shared_ptr<FEValuesBase<dim, spacedim> > > fevalv;
79  public:
80  static const unsigned int dimension = dim;
81  static const unsigned int space_dimension = spacedim;
82 
87 
92 
111  template <class FEVALUES>
113  const Mapping<dim,spacedim> &mapping,
115  const UpdateFlags flags,
116  const BlockInfo *local_block_info = nullptr);
117 
121  void initialize_data(const std::shared_ptr<VectorDataBase<dim,spacedim> > &data);
122 
126  void clear();
127 
133 
135  bool multigrid;
137 
142  const FEValuesBase<dim, spacedim> &fe_values () const;
143 
145 
149  const FEValuesBase<dim, spacedim> &fe_values (const unsigned int i) const;
150 
159  std::vector<std::vector<std::vector<double> > > values;
160 
169  std::vector<std::vector<std::vector<Tensor<1,spacedim> > > > gradients;
170 
179  std::vector<std::vector<std::vector<Tensor<2,spacedim> > > > hessians;
180 
184  template <typename number>
185  void reinit(const DoFInfo<dim, spacedim, number> &i);
186 
191  template<typename number>
192  void fill_local_data(const DoFInfo<dim, spacedim, number> &info, bool split_fevalues);
193 
198  std::shared_ptr<VectorDataBase<dim, spacedim> > global_data;
199 
203  std::size_t memory_consumption () const;
204 
205  private:
210 
216  template <typename TYPE>
217  void fill_local_data(
218  std::vector<std::vector<std::vector<TYPE> > > &data,
219  VectorSelector &selector,
220  bool split_fevalues) const;
224  unsigned int n_components;
225  };
226 
281  template <int dim, int spacedim=dim>
283  {
284  public:
285 
290 
295 
304  const Mapping<dim, spacedim> &mapping,
305  const BlockInfo *block_info = nullptr);
306 
314  template <typename VectorType>
316  const Mapping<dim, spacedim> &mapping,
317  const AnyData &data,
318  const VectorType &dummy,
319  const BlockInfo *block_info = nullptr);
327  template <typename VectorType>
329  const Mapping<dim, spacedim> &mapping,
330  const AnyData &data,
331  const MGLevelObject<VectorType> &dummy,
332  const BlockInfo *block_info = nullptr);
336  /* @{ */
337 
348  void initialize_update_flags(bool neighbor_geometry = false);
349 
354  void add_update_flags_all (const UpdateFlags flags);
355 
359  void add_update_flags_cell(const UpdateFlags flags);
360 
364  void add_update_flags_boundary(const UpdateFlags flags);
365 
369  void add_update_flags_face(const UpdateFlags flags);
370 
377  void add_update_flags(const UpdateFlags flags,
378  const bool cell = true,
379  const bool boundary = true,
380  const bool face = true,
381  const bool neighbor = true);
382 
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);
396 
400  std::size_t memory_consumption () const;
401 
414 
421 
429 
434 
439 
444  /* @} */
445 
449  /* @{ */
450 
465 
471 
477 
478  std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
479  std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
480  std::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
481  /* @} */
482 
486  /* @{ */
504  template <class DOFINFO>
505  void post_cell(const DoFInfoBox<dim, DOFINFO> &);
506 
523  template <class DOFINFO>
524  void post_faces(const DoFInfoBox<dim, DOFINFO> &);
525 
529  CellInfo cell;
533  CellInfo boundary;
537  CellInfo face;
542  CellInfo subface;
546  CellInfo neighbor;
547 
548  /* @} */
549  };
550 
551 
552 //----------------------------------------------------------------------//
553 
554  template<int dim, int sdim>
555  inline
557  :
558  fevalv(0),
559  multigrid(false),
560  global_data(std::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>)),
561  n_components(numbers::invalid_unsigned_int)
562  {}
563 
564 
565  template<int dim, int sdim>
566  inline
568  :
569  multigrid(other.multigrid),
570  values(other.values),
571  gradients(other.gradients),
572  hessians(other.hessians),
573  global_data(other.global_data),
574  fe_pointer(other.fe_pointer),
575  n_components(other.n_components)
576  {
577  fevalv.resize(other.fevalv.size());
578  for (unsigned int i=0; i<other.fevalv.size(); ++i)
579  {
580  const FEValuesBase<dim,sdim> &p = *other.fevalv[i];
581  const FEValues<dim,sdim> *pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
582  const FEFaceValues<dim,sdim> *pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
583  const FESubfaceValues<dim,sdim> *ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
584 
585  if (pc != nullptr)
586  fevalv[i] = std::shared_ptr<FEValuesBase<dim,sdim> > (
587  new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
588  pc->get_quadrature(), pc->get_update_flags()));
589  else if (pf != nullptr)
590  fevalv[i] = std::shared_ptr<FEValuesBase<dim,sdim> > (
591  new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
592  else if (ps != nullptr)
593  fevalv[i] = std::shared_ptr<FEValuesBase<dim,sdim> > (
595  else
596  Assert(false, ExcInternalError());
597  }
598  }
599 
600 
601 
602  template<int dim, int sdim>
603  template <class FEVALUES>
604  inline void
606  const FiniteElement<dim,sdim> &el,
607  const Mapping<dim,sdim> &mapping,
609  const UpdateFlags flags,
610  const BlockInfo *block_info)
611  {
612  fe_pointer = &el;
613  if (block_info == nullptr || block_info->local().size() == 0)
614  {
615  fevalv.resize(1);
616  fevalv[0] = std::shared_ptr<FEValuesBase<dim,sdim> > (
617  new FEVALUES (mapping, el, quadrature, flags));
618  }
619  else
620  {
621  fevalv.resize(el.n_base_elements());
622  for (unsigned int i=0; i<fevalv.size(); ++i)
623  {
624  fevalv[i] = std::shared_ptr<FEValuesBase<dim,sdim> > (
625  new FEVALUES (mapping, el.base_element(i), quadrature, flags));
626  }
627  }
628  n_components = el.n_components();
629  }
630 
631 
632  template <int dim, int spacedim>
633  inline const FiniteElement<dim, spacedim> &
635  {
636  Assert (fe_pointer != nullptr, ExcNotInitialized());
637  return *fe_pointer;
638  }
639 
640  template <int dim, int spacedim>
641  inline const FEValuesBase<dim, spacedim> &
643  {
644  AssertDimension(fevalv.size(), 1);
645  return *fevalv[0];
646  }
647 
648 
649  template <int dim, int spacedim>
650  inline const FEValuesBase<dim, spacedim> &
652  {
653  Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
654  return *fevalv[i];
655  }
656 
657 
658  template <int dim, int spacedim>
659  template <typename number>
660  inline void
662  {
663  for (unsigned int i=0; i<fevalv.size(); ++i)
664  {
665  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
667  {
668  // This is a subface
669  FESubfaceValues<dim,spacedim> &fe = dynamic_cast<FESubfaceValues<dim,spacedim>&> (febase);
670  fe.reinit(info.cell, info.face_number, info.sub_number);
671  }
673  {
674  // This is a face
675  FEFaceValues<dim,spacedim> &fe = dynamic_cast<FEFaceValues<dim,spacedim>&> (febase);
676  fe.reinit(info.cell, info.face_number);
677  }
678  else
679  {
680  // This is a cell
681  FEValues<dim,spacedim> &fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
682  fe.reinit(info.cell);
683  }
684  }
685 
686  const bool split_fevalues = info.block_info != nullptr;
687  if (!global_data->empty())
688  fill_local_data(info, split_fevalues);
689  }
690 
691 
692 
693 
694 //----------------------------------------------------------------------//
695 
696  template <int dim, int sdim>
697  inline
698  void
700  unsigned int cp,
701  unsigned int bp,
702  unsigned int fp,
703  bool force)
704  {
705  if (force || cell_quadrature.size() == 0)
706  cell_quadrature = QGauss<dim>(cp);
707  if (force || boundary_quadrature.size() == 0)
708  boundary_quadrature = QGauss<dim-1>(bp);
709  if (force || face_quadrature.size() == 0)
710  face_quadrature = QGauss<dim-1>(fp);
711  }
712 
713 
714  template <int dim, int sdim>
715  inline
716  void
718  {
719  add_update_flags(flags, true, true, true, true);
720  }
721 
722 
723  template <int dim, int sdim>
724  inline
725  void
727  {
728  add_update_flags(flags, true, false, false, false);
729  }
730 
731 
732  template <int dim, int sdim>
733  inline
734  void
736  {
737  add_update_flags(flags, false, true, false, false);
738  }
739 
740 
741  template <int dim, int sdim>
742  inline
743  void
745  {
746  add_update_flags(flags, false, false, true, true);
747  }
748 
749 
750  template <int dim, int sdim>
751  inline
752  void
754  const FiniteElement<dim,sdim> &el,
755  const Mapping<dim,sdim> &mapping,
756  const BlockInfo *block_info)
757  {
758  initialize_update_flags();
759  initialize_gauss_quadrature(
760  (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
761  (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
762  (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
763 
764  cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
765  cell_flags, block_info);
766  boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
767  boundary_flags, block_info);
768  face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
769  face_flags, block_info);
770  subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
771  face_flags, block_info);
772  neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
773  neighbor_flags, block_info);
774  }
775 
776 
777  template <int dim, int sdim>
778  template <typename VectorType>
779  void
782  const Mapping<dim,sdim> &mapping,
783  const AnyData &data,
784  const VectorType &,
785  const BlockInfo *block_info)
786  {
787  initialize(el, mapping, block_info);
788  std::shared_ptr<VectorData<VectorType, dim, sdim> > p;
790 
791  p = std::shared_ptr<VectorData<VectorType, dim, sdim> >(new VectorData<VectorType, dim, sdim> (cell_selector));
792  // Public member function of parent class was not found without
793  // explicit cast
794  pp = &*p;
795  pp->initialize(data);
796  cell_data = p;
797  cell.initialize_data(p);
798 
799  p = std::shared_ptr<VectorData<VectorType, dim, sdim> >(new VectorData<VectorType, dim, sdim> (boundary_selector));
800  pp = &*p;
801  pp->initialize(data);
802  boundary_data = p;
803  boundary.initialize_data(p);
804 
805  p = std::shared_ptr<VectorData<VectorType, dim, sdim> >(new VectorData<VectorType, dim, sdim> (face_selector));
806  pp = &*p;
807  pp->initialize(data);
808  face_data = p;
809  face.initialize_data(p);
810  subface.initialize_data(p);
811  neighbor.initialize_data(p);
812  }
813 
814  template <int dim, int sdim>
815  template <typename VectorType>
816  void
819  const Mapping<dim,sdim> &mapping,
820  const AnyData &data,
822  const BlockInfo *block_info)
823  {
824  initialize(el, mapping, block_info);
825  std::shared_ptr<MGVectorData<VectorType, dim, sdim> > p;
827 
828  p = std::shared_ptr<MGVectorData<VectorType, dim, sdim> >(new MGVectorData<VectorType, dim, sdim> (cell_selector));
829  // Public member function of parent class was not found without
830  // explicit cast
831  pp = &*p;
832  pp->initialize(data);
833  cell_data = p;
834  cell.initialize_data(p);
835 
836  p = std::shared_ptr<MGVectorData<VectorType, dim, sdim> >(new MGVectorData<VectorType, dim, sdim> (boundary_selector));
837  pp = &*p;
838  pp->initialize(data);
839  boundary_data = p;
840  boundary.initialize_data(p);
841 
842  p = std::shared_ptr<MGVectorData<VectorType, dim, sdim> >(new MGVectorData<VectorType, dim, sdim> (face_selector));
843  pp = &*p;
844  pp->initialize(data);
845  face_data = p;
846  face.initialize_data(p);
847  subface.initialize_data(p);
848  neighbor.initialize_data(p);
849  }
850 
851  template <int dim, int sdim>
852  template <class DOFINFO>
853  void
855  {}
856 
857 
858  template <int dim, int sdim>
859  template <class DOFINFO>
860  void
862  {}
863 
864 
865 }
866 
867 DEAL_II_NAMESPACE_CLOSE
868 
869 #endif
Shape function values.
std::vector< std::vector< std::vector< double > > > values
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:73
static const unsigned int invalid_unsigned_int
Definition: types.h:170
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)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1146
void post_faces(const DoFInfoBox< dim, DOFINFO > &)
const Quadrature< dim-1 > & get_quadrature() const
Quadrature< dim-1 > boundary_quadrature
void add_update_flags_cell(const UpdateFlags flags)
std::shared_ptr< VectorDataBase< dim, spacedim > > global_data
std::vector< std::shared_ptr< FEValuesBase< dim, spacedim > > > fevalv
vector of FEValues objects
void add_update_flags_all(const UpdateFlags flags)
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)
STL namespace.
unsigned int sub_number
Definition: dof_info.h:92
static::ExceptionBase & ExcNotInitialized()
std::vector< std::vector< std::vector< Tensor< 2, spacedim > > > > hessians
void initialize_update_flags(bool neighbor_geometry=false)
const FiniteElement< dim, spacedim > & finite_element() const
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void add_update_flags_face(const UpdateFlags flags)
UpdateFlags get_update_flags() const
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell, const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:4089
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
std::vector< std::vector< std::vector< Tensor< 1, spacedim > > > > gradients
MeshWorker::VectorSelector cell_selector
void initialize(const FiniteElement< dim, spacedim > &el, const Mapping< dim, spacedim > &mapping, const BlockInfo *block_info=nullptr)
#define Assert(cond, exc)
Definition: exceptions.h:313
UpdateFlags
void post_cell(const DoFInfoBox< dim, DOFINFO > &)
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell, const unsigned int face_no)
Definition: fe_values.cc:3923
Abstract base class for mapping classes.
Definition: dof_tools.h:46
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:161
unsigned int n_base_elements() const
Definition: fe.h:2808
unsigned int n_components() const
const Mapping< dim, spacedim > & get_mapping() const
void add_update_flags_boundary(const UpdateFlags flags)
std::size_t memory_consumption() const
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.
void add_update_flags(const UpdateFlags flags, const bool cell=true, const bool boundary=true, const bool face=true, const bool neighbor=true)
Definition: fe.h:30
unsigned int face_number
Definition: dof_info.h:84
const BlockIndices & local() const
Definition: block_info.h:217
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:93
unsigned int size() const
const FiniteElement< dim, spacedim > & get_fe() const
IntegrationInfo< dim, spacedim > CellInfo
SmartPointer< const FiniteElement< dim, spacedim >, IntegrationInfo< dim, spacedim > > fe_pointer
void initialize(const AnyData &)
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
Definition: fe.cc:1182
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell)
Definition: fe_values.cc:3725
std::size_t memory_consumption() const
void initialize_data(const std::shared_ptr< VectorDataBase< dim, spacedim > > &data)
MeshWorker::VectorSelector face_selector
static::ExceptionBase & ExcInternalError()