Reference documentation for deal.II version Git d86a2ae 2015-04-20 23:30:07 +0200
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
integration_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2015 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 __deal2__mesh_worker_integration_info_h
18 #define __deal2__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/base/std_cxx11/shared_ptr.h>
23 #include <deal.II/dofs/block_info.h>
24 #include <deal.II/fe/fe_values.h>
25 #include <deal.II/meshworker/local_results.h>
26 #include <deal.II/meshworker/dof_info.h>
27 #include <deal.II/meshworker/vector_selector.h>
28 
29 DEAL_II_NAMESPACE_OPEN
30 
31 namespace MeshWorker
32 {
72  template<int dim, int spacedim = dim>
74  {
75  private:
77  std::vector<std_cxx11::shared_ptr<FEValuesBase<dim, spacedim> > > fevalv;
78  public:
79  static const unsigned int dimension = dim;
80  static const unsigned int space_dimension = spacedim;
81 
86 
91 
110  template <class FEVALUES>
112  const Mapping<dim,spacedim> &mapping,
114  const UpdateFlags flags,
115  const BlockInfo *local_block_info = 0);
116 
120  void initialize_data(const std_cxx11::shared_ptr<VectorDataBase<dim,spacedim> > &data);
121 
125  void clear();
126 
132 
134  bool multigrid;
136 
141  const FEValuesBase<dim, spacedim> &fe_values () const;
142 
144 
148  const FEValuesBase<dim, spacedim> &fe_values (const unsigned int i) const;
149 
158  std::vector<std::vector<std::vector<double> > > values;
159 
168  std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
169 
178  std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
179 
183  template <typename number>
184  void reinit(const DoFInfo<dim, spacedim, number> &i);
185 
190  template<typename number>
191  void fill_local_data(const DoFInfo<dim, spacedim, number> &info, bool split_fevalues);
192 
197  std_cxx11::shared_ptr<VectorDataBase<dim, spacedim> > global_data;
198 
202  std::size_t memory_consumption () const;
203 
204  private:
209 
215  template <typename TYPE>
216  void fill_local_data(
217  std::vector<std::vector<std::vector<TYPE> > > &data,
218  VectorSelector &selector,
219  bool split_fevalues) const;
223  unsigned int n_components;
224  };
225 
280  template <int dim, int spacedim=dim>
282  {
283  public:
284 
289 
294 
303  const Mapping<dim, spacedim> &mapping,
304  const BlockInfo *block_info = 0);
305 
313  template <typename VECTOR>
315  const Mapping<dim, spacedim> &mapping,
316  const AnyData &data,
317  const VECTOR &dummy,
318  const BlockInfo *block_info = 0);
326  template <typename VECTOR>
328  const Mapping<dim, spacedim> &mapping,
329  const AnyData &data,
330  const MGLevelObject<VECTOR> &dummy,
331  const BlockInfo *block_info = 0);
335  /* @{ */
336 
347  void initialize_update_flags(bool neighbor_geometry = false);
348 
353  void add_update_flags_all (const UpdateFlags flags);
354 
358  void add_update_flags_cell(const UpdateFlags flags);
359 
363  void add_update_flags_boundary(const UpdateFlags flags);
364 
368  void add_update_flags_face(const UpdateFlags flags);
369 
376  void add_update_flags(const UpdateFlags flags,
377  const bool cell = true,
378  const bool boundary = true,
379  const bool face = true,
380  const bool neighbor = true);
381 
391  void initialize_gauss_quadrature(unsigned int n_cell_points,
392  unsigned int n_boundary_points,
393  unsigned int n_face_points,
394  const bool force = true);
395 
399  std::size_t memory_consumption () const;
400 
413 
420 
428 
433 
438 
443  /* @} */
444 
448  /* @{ */
449 
464 
470 
476 
477  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
478  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
479  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
480  /* @} */
481 
485  /* @{ */
503  template <class DOFINFO>
504  void post_cell(const DoFInfoBox<dim, DOFINFO> &);
505 
522  template <class DOFINFO>
523  void post_faces(const DoFInfoBox<dim, DOFINFO> &);
524 
546 
547  /* @} */
548  };
549 
550 
551 //----------------------------------------------------------------------//
552 
553  template<int dim, int sdim>
554  inline
556  :
557  fevalv(0),
558  multigrid(false),
559  global_data(std_cxx11::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>))
560  {}
561 
562 
563  template<int dim, int sdim>
564  inline
566  :
567  multigrid(other.multigrid),
568  values(other.values),
569  gradients(other.gradients),
570  hessians(other.hessians),
571  global_data(other.global_data),
572  fe_pointer(other.fe_pointer),
573  n_components(other.n_components)
574  {
575  fevalv.resize(other.fevalv.size());
576  for (unsigned int i=0; i<other.fevalv.size(); ++i)
577  {
578  const FEValuesBase<dim,sdim> &p = *other.fevalv[i];
579  const FEValues<dim,sdim> *pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
580  const FEFaceValues<dim,sdim> *pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
581  const FESubfaceValues<dim,sdim> *ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
582 
583  if (pc != 0)
584  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
585  new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
586  pc->get_quadrature(), pc->get_update_flags()));
587  else if (pf != 0)
588  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
589  new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
590  else if (ps != 0)
591  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
593  else
594  Assert(false, ExcInternalError());
595  }
596  }
597 
598 
599 
600  template<int dim, int sdim>
601  template <class FEVALUES>
602  inline void
604  const FiniteElement<dim,sdim> &el,
605  const Mapping<dim,sdim> &mapping,
607  const UpdateFlags flags,
608  const BlockInfo *block_info)
609  {
610  fe_pointer = &el;
611  if (block_info == 0 || block_info->local().size() == 0)
612  {
613  fevalv.resize(1);
614  fevalv[0] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
615  new FEVALUES (mapping, el, quadrature, flags));
616  }
617  else
618  {
619  fevalv.resize(el.n_base_elements());
620  for (unsigned int i=0; i<fevalv.size(); ++i)
621  {
622  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
623  new FEVALUES (mapping, el.base_element(i), quadrature, flags));
624  }
625  }
626  n_components = el.n_components();
627  }
628 
629 
630  template <int dim, int spacedim>
631  inline const FiniteElement<dim, spacedim> &
633  {
634  Assert (fe_pointer !=0, ExcNotInitialized());
635  return *fe_pointer;
636  }
637 
638  template <int dim, int spacedim>
639  inline const FEValuesBase<dim, spacedim> &
641  {
642  AssertDimension(fevalv.size(), 1);
643  return *fevalv[0];
644  }
645 
646 
647  template <int dim, int spacedim>
648  inline const FEValuesBase<dim, spacedim> &
650  {
651  Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
652  return *fevalv[i];
653  }
654 
655 
656  template <int dim, int spacedim>
657  template <typename number>
658  inline void
660  {
661  for (unsigned int i=0; i<fevalv.size(); ++i)
662  {
663  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
665  {
666  // This is a subface
667  FESubfaceValues<dim> &fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
668  fe.reinit(info.cell, info.face_number, info.sub_number);
669  }
671  {
672  // This is a face
673  FEFaceValues<dim> &fe = dynamic_cast<FEFaceValues<dim>&> (febase);
674  fe.reinit(info.cell, info.face_number);
675  }
676  else
677  {
678  // This is a cell
679  FEValues<dim,spacedim> &fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
680  fe.reinit(info.cell);
681  }
682  }
683 
684  const bool split_fevalues = info.block_info != 0;
685  if (!global_data->empty())
686  fill_local_data(info, split_fevalues);
687  }
688 
689 
690 
691 
692 //----------------------------------------------------------------------//
693 
694  template <int dim, int sdim>
695  inline
696  void
698  unsigned int cp,
699  unsigned int bp,
700  unsigned int fp,
701  bool force)
702  {
703  if (force || cell_quadrature.size() == 0)
704  cell_quadrature = QGauss<dim>(cp);
705  if (force || boundary_quadrature.size() == 0)
706  boundary_quadrature = QGauss<dim-1>(bp);
707  if (force || face_quadrature.size() == 0)
708  face_quadrature = QGauss<dim-1>(fp);
709  }
710 
711 
712  template <int dim, int sdim>
713  inline
714  void
716  {
717  add_update_flags(flags, true, true, true, true);
718  }
719 
720 
721  template <int dim, int sdim>
722  inline
723  void
725  {
726  add_update_flags(flags, true, false, false, false);
727  }
728 
729 
730  template <int dim, int sdim>
731  inline
732  void
734  {
735  add_update_flags(flags, false, true, false, false);
736  }
737 
738 
739  template <int dim, int sdim>
740  inline
741  void
743  {
744  add_update_flags(flags, false, false, true, true);
745  }
746 
747 
748  template <int dim, int sdim>
749  inline
750  void
752  const FiniteElement<dim,sdim> &el,
753  const Mapping<dim,sdim> &mapping,
754  const BlockInfo *block_info)
755  {
756  initialize_update_flags();
757  initialize_gauss_quadrature(
758  (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
759  (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
760  (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
761 
762  cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
763  cell_flags, block_info);
764  boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
765  boundary_flags, block_info);
766  face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
767  face_flags, block_info);
768  subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
769  face_flags, block_info);
770  neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
771  neighbor_flags, block_info);
772  }
773 
774 
775  template <int dim, int sdim>
776  template <typename VECTOR>
777  void
779  const FiniteElement<dim,sdim> &el,
780  const Mapping<dim,sdim> &mapping,
781  const AnyData &data,
782  const VECTOR &,
783  const BlockInfo *block_info)
784  {
785  initialize(el, mapping, block_info);
786  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
788 
789  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
790  // Public member function of parent class was not found without
791  // explicit cast
792  pp = &*p;
793  pp->initialize(data);
794  cell_data = p;
795  cell.initialize_data(p);
796 
797  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
798  pp = &*p;
799  pp->initialize(data);
800  boundary_data = p;
801  boundary.initialize_data(p);
802 
803  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
804  pp = &*p;
805  pp->initialize(data);
806  face_data = p;
807  face.initialize_data(p);
808  subface.initialize_data(p);
809  neighbor.initialize_data(p);
810  }
811 
812  template <int dim, int sdim>
813  template <typename VECTOR>
814  void
816  const FiniteElement<dim,sdim> &el,
817  const Mapping<dim,sdim> &mapping,
818  const AnyData &data,
819  const MGLevelObject<VECTOR> &,
820  const BlockInfo *block_info)
821  {
822  initialize(el, mapping, block_info);
823  std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
825 
826  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (cell_selector));
827  // Public member function of parent class was not found without
828  // explicit cast
829  pp = &*p;
830  pp->initialize(data);
831  cell_data = p;
832  cell.initialize_data(p);
833 
834  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (boundary_selector));
835  pp = &*p;
836  pp->initialize(data);
837  boundary_data = p;
838  boundary.initialize_data(p);
839 
840  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (face_selector));
841  pp = &*p;
842  pp->initialize(data);
843  face_data = p;
844  face.initialize_data(p);
845  subface.initialize_data(p);
846  neighbor.initialize_data(p);
847  }
848 
849  template <int dim, int sdim>
850  template <class DOFINFO>
851  void
853  {}
854 
855 
856  template <int dim, int sdim>
857  template <class DOFINFO>
858  void
860  {}
861 
862 
863 }
864 
865 DEAL_II_NAMESPACE_CLOSE
866 
867 #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:72
static const unsigned int invalid_unsigned_int
Definition: types.h:153
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > &cell, const unsigned int face_no, const unsigned int subface_no)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1021
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
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)
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)
unsigned int sub_number
Definition: dof_info.h:91
void initialize_update_flags(bool neighbor_geometry=false)
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)
const FiniteElement< dim, spacedim > & finite_element() const
void add_update_flags_face(const UpdateFlags flags)
UpdateFlags get_update_flags() const
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > &cell, const unsigned int face_no)
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)
Definition: exceptions.h:293
UpdateFlags
void post_cell(const DoFInfoBox< dim, DOFINFO > &)
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:160
void initialize_data(const std_cxx11::shared_ptr< VectorDataBase< dim, spacedim > > &data)
unsigned int n_base_elements() const
unsigned int n_components() const
void initialize(const AnyData &)
std::vector< std_cxx11::shared_ptr< FEValuesBase< dim, spacedim > > > fevalv
vector of FEValues objects
const Mapping< dim, spacedim > & get_mapping() const
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > &cell)
void add_update_flags_boundary(const UpdateFlags flags)
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
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:31
unsigned int face_number
Definition: dof_info.h:83
const BlockIndices & local() const
Definition: block_info.h:211
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:87
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
std::vector< std::vector< std::vector< Tensor< 2, dim > > > > hessians
std::size_t memory_consumption() const
std_cxx11::shared_ptr< VectorDataBase< dim, spacedim > > global_data
MeshWorker::VectorSelector face_selector