Reference documentation for deal.II version Git 9deac93 2014-10-24 21:41:49 -0500
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
integration_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2014 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 {
76  template<int dim, int spacedim = dim>
78  {
79  private:
81  std::vector<std_cxx11::shared_ptr<FEValuesBase<dim, spacedim> > > fevalv;
82  public:
83  static const unsigned int dimension = dim;
84  static const unsigned int space_dimension = spacedim;
85 
90 
96 
117  template <class FEVALUES>
119  const Mapping<dim,spacedim> &mapping,
121  const UpdateFlags flags,
122  const BlockInfo *local_block_info = 0);
123 
127  void initialize_data(const std_cxx11::shared_ptr<VectorDataBase<dim,spacedim> > &data);
128 
132  void clear();
133 
139 
141  bool multigrid;
143 
148  const FEValuesBase<dim, spacedim> &fe_values () const;
149 
151 
155  const FEValuesBase<dim, spacedim> &fe_values (const unsigned int i) const;
156 
165  std::vector<std::vector<std::vector<double> > > values;
166 
175  std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
176 
185  std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
186 
190  template <typename number>
191  void reinit(const DoFInfo<dim, spacedim, number> &i);
192 
197  template<typename number>
198  void fill_local_data(const DoFInfo<dim, spacedim, number> &info, bool split_fevalues);
199 
204  std_cxx11::shared_ptr<VectorDataBase<dim, spacedim> > global_data;
205 
209  std::size_t memory_consumption () const;
210 
211  private:
216 
222  template <typename TYPE>
223  void fill_local_data(
224  std::vector<std::vector<std::vector<TYPE> > > &data,
225  VectorSelector &selector,
226  bool split_fevalues) const;
230  unsigned int n_components;
231  };
232 
289  template <int dim, int spacedim=dim>
291  {
292  public:
293 
299 
304 
313  const Mapping<dim, spacedim> &mapping,
314  const BlockInfo *block_info = 0);
315 
323  template <typename VECTOR>
325  const Mapping<dim, spacedim> &mapping,
326  const AnyData &data,
327  const VECTOR &dummy,
328  const BlockInfo *block_info = 0);
332  template <typename VECTOR>
334  const Mapping<dim, spacedim> &mapping,
335  const NamedData<VECTOR *> &data,
336  const BlockInfo *block_info = 0);
337 
341  template <typename VECTOR>
343  const Mapping<dim, spacedim> &mapping,
344  const NamedData<MGLevelObject<VECTOR>*> &data,
345  const BlockInfo *block_info = 0);
349  /* @{ */
350 
361  void initialize_update_flags(bool neighbor_geometry = false);
362 
367  void add_update_flags_all (const UpdateFlags flags);
368 
372  void add_update_flags_cell(const UpdateFlags flags);
373 
377  void add_update_flags_boundary(const UpdateFlags flags);
378 
382  void add_update_flags_face(const UpdateFlags flags);
383 
390  void add_update_flags(const UpdateFlags flags,
391  const bool cell = true,
392  const bool boundary = true,
393  const bool face = true,
394  const bool neighbor = true);
395 
405  void initialize_gauss_quadrature(unsigned int n_cell_points,
406  unsigned int n_boundary_points,
407  unsigned int n_face_points,
408  const bool force = true);
409 
413  std::size_t memory_consumption () const;
414 
427 
434 
442 
447 
452 
457  /* @} */
458 
462  /* @{ */
463 
479 
485 
491 
492  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
493  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
494  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
495  /* @} */
496 
500  /* @{ */
519  template <class DOFINFO>
520  void post_cell(const DoFInfoBox<dim, DOFINFO> &);
521 
539  template <class DOFINFO>
540  void post_faces(const DoFInfoBox<dim, DOFINFO> &);
541 
564 
565  /* @} */
566  };
567 
568 
569 //----------------------------------------------------------------------//
570 
571  template<int dim, int sdim>
572  inline
574  :
575  fevalv(0),
576  multigrid(false),
577  global_data(std_cxx11::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>))
578  {}
579 
580 
581  template<int dim, int sdim>
582  inline
584  :
585  multigrid(other.multigrid),
586  values(other.values),
587  gradients(other.gradients),
588  hessians(other.hessians),
589  global_data(other.global_data),
590  fe_pointer(other.fe_pointer),
591  n_components(other.n_components)
592  {
593  fevalv.resize(other.fevalv.size());
594  for (unsigned int i=0; i<other.fevalv.size(); ++i)
595  {
596  const FEValuesBase<dim,sdim> &p = *other.fevalv[i];
597  const FEValues<dim,sdim> *pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
598  const FEFaceValues<dim,sdim> *pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
599  const FESubfaceValues<dim,sdim> *ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
600 
601  if (pc != 0)
602  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
603  new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
604  pc->get_quadrature(), pc->get_update_flags()));
605  else if (pf != 0)
606  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
607  new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
608  else if (ps != 0)
609  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
611  else
612  Assert(false, ExcInternalError());
613  }
614  }
615 
616 
617 
618  template<int dim, int sdim>
619  template <class FEVALUES>
620  inline void
622  const FiniteElement<dim,sdim> &el,
623  const Mapping<dim,sdim> &mapping,
625  const UpdateFlags flags,
626  const BlockInfo *block_info)
627  {
628  fe_pointer = &el;
629  if (block_info == 0 || block_info->local().size() == 0)
630  {
631  fevalv.resize(1);
632  fevalv[0] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
633  new FEVALUES (mapping, el, quadrature, flags));
634  }
635  else
636  {
637  fevalv.resize(el.n_base_elements());
638  for (unsigned int i=0; i<fevalv.size(); ++i)
639  {
640  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
641  new FEVALUES (mapping, el.base_element(i), quadrature, flags));
642  }
643  }
644  n_components = el.n_components();
645  }
646 
647 
648  template <int dim, int spacedim>
649  inline const FiniteElement<dim, spacedim> &
651  {
652  Assert (fe_pointer !=0, ExcNotInitialized());
653  return *fe_pointer;
654  }
655 
656  template <int dim, int spacedim>
657  inline const FEValuesBase<dim, spacedim> &
659  {
660  AssertDimension(fevalv.size(), 1);
661  return *fevalv[0];
662  }
663 
664 
665  template <int dim, int spacedim>
666  inline const FEValuesBase<dim, spacedim> &
668  {
669  Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
670  return *fevalv[i];
671  }
672 
673 
674  template <int dim, int spacedim>
675  template <typename number>
676  inline void
678  {
679  for (unsigned int i=0; i<fevalv.size(); ++i)
680  {
681  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
682  if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
683  {
684  // This is a subface
685  FESubfaceValues<dim> &fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
686  fe.reinit(info.cell, info.face_number, info.sub_number);
687  }
688  else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
689  {
690  // This is a face
691  FEFaceValues<dim> &fe = dynamic_cast<FEFaceValues<dim>&> (febase);
692  fe.reinit(info.cell, info.face_number);
693  }
694  else
695  {
696  // This is a cell
697  FEValues<dim,spacedim> &fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
698  fe.reinit(info.cell);
699  }
700  }
701 
702  const bool split_fevalues = info.block_info != 0;
703  if (!global_data->empty())
704  fill_local_data(info, split_fevalues);
705  }
706 
707 
708 
709 
710 //----------------------------------------------------------------------//
711 
712  template <int dim, int sdim>
713  inline
714  void
716  unsigned int cp,
717  unsigned int bp,
718  unsigned int fp,
719  bool force)
720  {
721  if (force || cell_quadrature.size() == 0)
722  cell_quadrature = QGauss<dim>(cp);
723  if (force || boundary_quadrature.size() == 0)
724  boundary_quadrature = QGauss<dim-1>(bp);
725  if (force || face_quadrature.size() == 0)
726  face_quadrature = QGauss<dim-1>(fp);
727  }
728 
729 
730  template <int dim, int sdim>
731  inline
732  void
734  {
735  add_update_flags(flags, true, true, true, true);
736  }
737 
738 
739  template <int dim, int sdim>
740  inline
741  void
743  {
744  add_update_flags(flags, true, false, false, false);
745  }
746 
747 
748  template <int dim, int sdim>
749  inline
750  void
752  {
753  add_update_flags(flags, false, true, false, false);
754  }
755 
756 
757  template <int dim, int sdim>
758  inline
759  void
761  {
762  add_update_flags(flags, false, false, true, true);
763  }
764 
765 
766  template <int dim, int sdim>
767  inline
768  void
770  const FiniteElement<dim,sdim> &el,
771  const Mapping<dim,sdim> &mapping,
772  const BlockInfo *block_info)
773  {
774  initialize_update_flags();
775  initialize_gauss_quadrature(
776  (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
777  (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
778  (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
779 
780  cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
781  cell_flags, block_info);
782  boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
783  boundary_flags, block_info);
784  face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
785  face_flags, block_info);
786  subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
787  face_flags, block_info);
788  neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
789  neighbor_flags, block_info);
790  }
791 
792 
793  template <int dim, int sdim>
794  template <typename VECTOR>
795  void
797  const FiniteElement<dim,sdim> &el,
798  const Mapping<dim,sdim> &mapping,
799  const AnyData &data,
800  const VECTOR &dummy,
801  const BlockInfo *block_info)
802  {
803  initialize(el, mapping, block_info);
804  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
806 
807  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
808  // Public member function of parent class was not found without
809  // explicit cast
810  pp = &*p;
811  pp->initialize(data);
812  cell_data = p;
813  cell.initialize_data(p);
814 
815  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
816  pp = &*p;
817  pp->initialize(data);
818  boundary_data = p;
819  boundary.initialize_data(p);
820 
821  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
822  pp = &*p;
823  pp->initialize(data);
824  face_data = p;
825  face.initialize_data(p);
826  subface.initialize_data(p);
827  neighbor.initialize_data(p);
828  }
829 
830 
831  template <int dim, int sdim>
832  template <typename VECTOR>
833  void
835  const FiniteElement<dim,sdim> &el,
836  const Mapping<dim,sdim> &mapping,
837  const NamedData<VECTOR *> &data,
838  const BlockInfo *block_info)
839  {
840  initialize(el, mapping, block_info);
841  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
842 
843  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
844  p->initialize(data);
845  cell_data = p;
846  cell.initialize_data(p);
847 
848  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
849  p->initialize(data);
850  boundary_data = p;
851  boundary.initialize_data(p);
852 
853  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
854  p->initialize(data);
855  face_data = p;
856  face.initialize_data(p);
857  subface.initialize_data(p);
858  neighbor.initialize_data(p);
859  }
860 
861 
862  template <int dim, int sdim>
863  template <typename VECTOR>
864  void
866  const FiniteElement<dim,sdim> &el,
867  const Mapping<dim,sdim> &mapping,
868  const NamedData<MGLevelObject<VECTOR>*> &data,
869  const BlockInfo *block_info)
870  {
871  initialize(el, mapping, block_info);
872  std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
873 
874  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (cell_selector));
875  p->initialize(data);
876  cell_data = p;
877  cell.initialize_data(p);
878 
879  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (boundary_selector));
880  p->initialize(data);
881  boundary_data = p;
882  boundary.initialize_data(p);
883 
884  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (face_selector));
885  p->initialize(data);
886  face_data = p;
887  face.initialize_data(p);
888  subface.initialize_data(p);
889  neighbor.initialize_data(p);
890  }
891 
892 
893  template <int dim, int sdim>
894  template <class DOFINFO>
895  void
897  {}
898 
899 
900  template <int dim, int sdim>
901  template <class DOFINFO>
902  void
904  {}
905 
906 
907 }
908 
909 DEAL_II_NAMESPACE_CLOSE
910 
911 #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:74
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:857
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:100
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:298
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:174
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:31
unsigned int face_number
Definition: dof_info.h:88
const BlockIndices & local() const
Definition: block_info.h:250
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:92
::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