Reference documentation for deal.II version Git 9e5c6cc 2014-09-18 15:25:08 -0400
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
integration_info.h
1 // ---------------------------------------------------------------------
2 // @f$Id@f$
3 //
4 // Copyright (C) 2006 - 2014 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 
18 #ifndef __deal2__mesh_worker_integration_info_h
19 #define __deal2__mesh_worker_integration_info_h
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/quadrature_lib.h>
23 #include <deal.II/base/std_cxx11/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>
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 namespace MeshWorker
33 {
77  template<int dim, int spacedim = dim>
79  {
80  private:
82  std::vector<std_cxx11::shared_ptr<FEValuesBase<dim, spacedim> > > fevalv;
83  public:
84  static const unsigned int dimension = dim;
85  static const unsigned int space_dimension = spacedim;
86 
91 
97 
118  template <class FEVALUES>
120  const Mapping<dim,spacedim> &mapping,
122  const UpdateFlags flags,
123  const BlockInfo *local_block_info = 0);
124 
128  void initialize_data(const std_cxx11::shared_ptr<VectorDataBase<dim,spacedim> > &data);
129 
133  void clear();
134 
140 
142  bool multigrid;
144 
149  const FEValuesBase<dim, spacedim> &fe_values () const;
150 
152 
156  const FEValuesBase<dim, spacedim> &fe_values (const unsigned int i) const;
157 
166  std::vector<std::vector<std::vector<double> > > values;
167 
176  std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
177 
186  std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
187 
191  template <typename number>
192  void reinit(const DoFInfo<dim, spacedim, number> &i);
193 
198  template<typename number>
199  void fill_local_data(const DoFInfo<dim, spacedim, number> &info, bool split_fevalues);
200 
205  std_cxx11::shared_ptr<VectorDataBase<dim, spacedim> > global_data;
206 
210  std::size_t memory_consumption () const;
211 
212  private:
217 
223  template <typename TYPE>
224  void fill_local_data(
225  std::vector<std::vector<std::vector<TYPE> > > &data,
226  VectorSelector &selector,
227  bool split_fevalues) const;
231  unsigned int n_components;
232  };
233 
290  template <int dim, int spacedim=dim>
292  {
293  public:
294 
300 
305 
314  const Mapping<dim, spacedim> &mapping,
315  const BlockInfo *block_info = 0);
316 
324  template <typename VECTOR>
326  const Mapping<dim, spacedim> &mapping,
327  const AnyData &data,
328  const VECTOR &dummy,
329  const BlockInfo *block_info = 0);
333  template <typename VECTOR>
335  const Mapping<dim, spacedim> &mapping,
336  const NamedData<VECTOR *> &data,
337  const BlockInfo *block_info = 0);
338 
342  template <typename VECTOR>
344  const Mapping<dim, spacedim> &mapping,
345  const NamedData<MGLevelObject<VECTOR>*> &data,
346  const BlockInfo *block_info = 0);
350  /* @{ */
351 
362  void initialize_update_flags(bool neighbor_geometry = false);
363 
368  void add_update_flags_all (const UpdateFlags flags);
369 
373  void add_update_flags_cell(const UpdateFlags flags);
374 
378  void add_update_flags_boundary(const UpdateFlags flags);
379 
383  void add_update_flags_face(const UpdateFlags flags);
384 
391  void add_update_flags(const UpdateFlags flags,
392  const bool cell = true,
393  const bool boundary = true,
394  const bool face = true,
395  const bool neighbor = true);
396 
406  void initialize_gauss_quadrature(unsigned int n_cell_points,
407  unsigned int n_boundary_points,
408  unsigned int n_face_points,
409  const bool force = true);
410 
414  std::size_t memory_consumption () const;
415 
428 
435 
443 
448 
453 
458  /* @} */
459 
463  /* @{ */
464 
480 
486 
492 
493  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
494  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
495  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
496  /* @} */
497 
501  /* @{ */
520  template <class DOFINFO>
521  void post_cell(const DoFInfoBox<dim, DOFINFO> &);
522 
540  template <class DOFINFO>
541  void post_faces(const DoFInfoBox<dim, DOFINFO> &);
542 
565 
566  /* @} */
567  };
568 
569 
570 //----------------------------------------------------------------------//
571 
572  template<int dim, int sdim>
573  inline
575  :
576  fevalv(0),
577  multigrid(false),
578  global_data(std_cxx11::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>))
579  {}
580 
581 
582  template<int dim, int sdim>
583  inline
585  :
586  multigrid(other.multigrid),
587  values(other.values),
588  gradients(other.gradients),
589  hessians(other.hessians),
590  global_data(other.global_data),
591  fe_pointer(other.fe_pointer),
592  n_components(other.n_components)
593  {
594  fevalv.resize(other.fevalv.size());
595  for (unsigned int i=0; i<other.fevalv.size(); ++i)
596  {
597  const FEValuesBase<dim,sdim> &p = *other.fevalv[i];
598  const FEValues<dim,sdim> *pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
599  const FEFaceValues<dim,sdim> *pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
600  const FESubfaceValues<dim,sdim> *ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
601 
602  if (pc != 0)
603  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
604  new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
605  pc->get_quadrature(), pc->get_update_flags()));
606  else if (pf != 0)
607  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
608  new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
609  else if (ps != 0)
610  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
612  else
613  Assert(false, ExcInternalError());
614  }
615  }
616 
617 
618 
619  template<int dim, int sdim>
620  template <class FEVALUES>
621  inline void
623  const FiniteElement<dim,sdim> &el,
624  const Mapping<dim,sdim> &mapping,
626  const UpdateFlags flags,
627  const BlockInfo *block_info)
628  {
629  fe_pointer = &el;
630  if (block_info == 0 || block_info->local().size() == 0)
631  {
632  fevalv.resize(1);
633  fevalv[0] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
634  new FEVALUES (mapping, el, quadrature, flags));
635  }
636  else
637  {
638  fevalv.resize(el.n_base_elements());
639  for (unsigned int i=0; i<fevalv.size(); ++i)
640  {
641  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
642  new FEVALUES (mapping, el.base_element(i), quadrature, flags));
643  }
644  }
645  n_components = el.n_components();
646  }
647 
648 
649  template <int dim, int spacedim>
650  inline const FiniteElement<dim, spacedim> &
652  {
653  Assert (fe_pointer !=0, ExcNotInitialized());
654  return *fe_pointer;
655  }
656 
657  template <int dim, int spacedim>
658  inline const FEValuesBase<dim, spacedim> &
660  {
661  AssertDimension(fevalv.size(), 1);
662  return *fevalv[0];
663  }
664 
665 
666  template <int dim, int spacedim>
667  inline const FEValuesBase<dim, spacedim> &
669  {
670  Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
671  return *fevalv[i];
672  }
673 
674 
675  template <int dim, int spacedim>
676  template <typename number>
677  inline void
679  {
680  for (unsigned int i=0; i<fevalv.size(); ++i)
681  {
682  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
683  if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
684  {
685  // This is a subface
686  FESubfaceValues<dim> &fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
687  fe.reinit(info.cell, info.face_number, info.sub_number);
688  }
689  else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
690  {
691  // This is a face
692  FEFaceValues<dim> &fe = dynamic_cast<FEFaceValues<dim>&> (febase);
693  fe.reinit(info.cell, info.face_number);
694  }
695  else
696  {
697  // This is a cell
698  FEValues<dim,spacedim> &fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
699  fe.reinit(info.cell);
700  }
701  }
702 
703  const bool split_fevalues = info.block_info != 0;
704  if (!global_data->empty())
705  fill_local_data(info, split_fevalues);
706  }
707 
708 
709 
710 
711 //----------------------------------------------------------------------//
712 
713  template <int dim, int sdim>
714  inline
715  void
717  unsigned int cp,
718  unsigned int bp,
719  unsigned int fp,
720  bool force)
721  {
722  if (force || cell_quadrature.size() == 0)
723  cell_quadrature = QGauss<dim>(cp);
724  if (force || boundary_quadrature.size() == 0)
725  boundary_quadrature = QGauss<dim-1>(bp);
726  if (force || face_quadrature.size() == 0)
727  face_quadrature = QGauss<dim-1>(fp);
728  }
729 
730 
731  template <int dim, int sdim>
732  inline
733  void
735  {
736  add_update_flags(flags, true, true, true, true);
737  }
738 
739 
740  template <int dim, int sdim>
741  inline
742  void
744  {
745  add_update_flags(flags, true, false, false, false);
746  }
747 
748 
749  template <int dim, int sdim>
750  inline
751  void
753  {
754  add_update_flags(flags, false, true, false, false);
755  }
756 
757 
758  template <int dim, int sdim>
759  inline
760  void
762  {
763  add_update_flags(flags, false, false, true, true);
764  }
765 
766 
767  template <int dim, int sdim>
768  inline
769  void
771  const FiniteElement<dim,sdim> &el,
772  const Mapping<dim,sdim> &mapping,
773  const BlockInfo *block_info)
774  {
775  initialize_update_flags();
776  initialize_gauss_quadrature(
777  (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
778  (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
779  (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
780 
781  cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
782  cell_flags, block_info);
783  boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
784  boundary_flags, block_info);
785  face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
786  face_flags, block_info);
787  subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
788  face_flags, block_info);
789  neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
790  neighbor_flags, block_info);
791  }
792 
793 
794  template <int dim, int sdim>
795  template <typename VECTOR>
796  void
798  const FiniteElement<dim,sdim> &el,
799  const Mapping<dim,sdim> &mapping,
800  const AnyData &data,
801  const VECTOR &dummy,
802  const BlockInfo *block_info)
803  {
804  initialize(el, mapping, block_info);
805  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
807 
808  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
809  // Public member function of parent class was not found without
810  // explicit cast
811  pp = &*p;
812  pp->initialize(data);
813  cell_data = p;
814  cell.initialize_data(p);
815 
816  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
817  pp = &*p;
818  pp->initialize(data);
819  boundary_data = p;
820  boundary.initialize_data(p);
821 
822  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
823  pp = &*p;
824  pp->initialize(data);
825  face_data = p;
826  face.initialize_data(p);
827  subface.initialize_data(p);
828  neighbor.initialize_data(p);
829  }
830 
831 
832  template <int dim, int sdim>
833  template <typename VECTOR>
834  void
836  const FiniteElement<dim,sdim> &el,
837  const Mapping<dim,sdim> &mapping,
838  const NamedData<VECTOR *> &data,
839  const BlockInfo *block_info)
840  {
841  initialize(el, mapping, block_info);
842  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
843 
844  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
845  p->initialize(data);
846  cell_data = p;
847  cell.initialize_data(p);
848 
849  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
850  p->initialize(data);
851  boundary_data = p;
852  boundary.initialize_data(p);
853 
854  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
855  p->initialize(data);
856  face_data = p;
857  face.initialize_data(p);
858  subface.initialize_data(p);
859  neighbor.initialize_data(p);
860  }
861 
862 
863  template <int dim, int sdim>
864  template <typename VECTOR>
865  void
867  const FiniteElement<dim,sdim> &el,
868  const Mapping<dim,sdim> &mapping,
869  const NamedData<MGLevelObject<VECTOR>*> &data,
870  const BlockInfo *block_info)
871  {
872  initialize(el, mapping, block_info);
873  std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
874 
875  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (cell_selector));
876  p->initialize(data);
877  cell_data = p;
878  cell.initialize_data(p);
879 
880  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (boundary_selector));
881  p->initialize(data);
882  boundary_data = p;
883  boundary.initialize_data(p);
884 
885  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (face_selector));
886  p->initialize(data);
887  face_data = p;
888  face.initialize_data(p);
889  subface.initialize_data(p);
890  neighbor.initialize_data(p);
891  }
892 
893 
894  template <int dim, int sdim>
895  template <class DOFINFO>
896  void
898  {}
899 
900 
901  template <int dim, int sdim>
902  template <class DOFINFO>
903  void
905  {}
906 
907 
908 }
909 
910 DEAL_II_NAMESPACE_CLOSE
911 
912 #endif
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell, const unsigned int face_no)
Shape function values.
std::vector< std::vector< std::vector< double > > > values
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:75
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
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:101
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
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:299
UpdateFlags
void post_cell(const DoFInfoBox< dim, DOFINFO > &)
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell)
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:175
void initialize_data(const std_cxx11::shared_ptr< VectorDataBase< dim, spacedim > > &data)
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
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 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:32
unsigned int face_number
Definition: dof_info.h:89
const BlockIndices & local() const
Definition: block_info.h:251
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:93
::ExceptionBase & ExcNotInitialized()
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
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell, const unsigned int face_no, const unsigned int subface_no)
::ExceptionBase & ExcInternalError()
std::size_t memory_consumption() const
std_cxx11::shared_ptr< VectorDataBase< dim, spacedim > > global_data
MeshWorker::VectorSelector face_selector