Reference documentation for deal.II version Git c77b683 2014-12-17 17:11:36 +0100
 All Classes Namespaces 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 {
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);
322  template <typename VECTOR>
324  const Mapping<dim, spacedim> &mapping,
325  const NamedData<VECTOR *> &data,
326  const BlockInfo *block_info = 0);
327 
331  template <typename VECTOR>
333  const Mapping<dim, spacedim> &mapping,
334  const NamedData<MGLevelObject<VECTOR>*> &data,
335  const BlockInfo *block_info = 0);
339  /* @{ */
340 
351  void initialize_update_flags(bool neighbor_geometry = false);
352 
357  void add_update_flags_all (const UpdateFlags flags);
358 
362  void add_update_flags_cell(const UpdateFlags flags);
363 
367  void add_update_flags_boundary(const UpdateFlags flags);
368 
372  void add_update_flags_face(const UpdateFlags flags);
373 
380  void add_update_flags(const UpdateFlags flags,
381  const bool cell = true,
382  const bool boundary = true,
383  const bool face = true,
384  const bool neighbor = true);
385 
395  void initialize_gauss_quadrature(unsigned int n_cell_points,
396  unsigned int n_boundary_points,
397  unsigned int n_face_points,
398  const bool force = true);
399 
403  std::size_t memory_consumption () const;
404 
417 
424 
432 
437 
442 
447  /* @} */
448 
452  /* @{ */
453 
468 
474 
480 
481  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
482  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
483  std_cxx11::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
484  /* @} */
485 
489  /* @{ */
507  template <class DOFINFO>
508  void post_cell(const DoFInfoBox<dim, DOFINFO> &);
509 
526  template <class DOFINFO>
527  void post_faces(const DoFInfoBox<dim, DOFINFO> &);
528 
550 
551  /* @} */
552  };
553 
554 
555 //----------------------------------------------------------------------//
556 
557  template<int dim, int sdim>
558  inline
560  :
561  fevalv(0),
562  multigrid(false),
563  global_data(std_cxx11::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>))
564  {}
565 
566 
567  template<int dim, int sdim>
568  inline
570  :
571  multigrid(other.multigrid),
572  values(other.values),
573  gradients(other.gradients),
574  hessians(other.hessians),
575  global_data(other.global_data),
576  fe_pointer(other.fe_pointer),
577  n_components(other.n_components)
578  {
579  fevalv.resize(other.fevalv.size());
580  for (unsigned int i=0; i<other.fevalv.size(); ++i)
581  {
582  const FEValuesBase<dim,sdim> &p = *other.fevalv[i];
583  const FEValues<dim,sdim> *pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
584  const FEFaceValues<dim,sdim> *pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
585  const FESubfaceValues<dim,sdim> *ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
586 
587  if (pc != 0)
588  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
589  new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
590  pc->get_quadrature(), pc->get_update_flags()));
591  else if (pf != 0)
592  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
593  new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
594  else if (ps != 0)
595  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
597  else
598  Assert(false, ExcInternalError());
599  }
600  }
601 
602 
603 
604  template<int dim, int sdim>
605  template <class FEVALUES>
606  inline void
608  const FiniteElement<dim,sdim> &el,
609  const Mapping<dim,sdim> &mapping,
611  const UpdateFlags flags,
612  const BlockInfo *block_info)
613  {
614  fe_pointer = &el;
615  if (block_info == 0 || block_info->local().size() == 0)
616  {
617  fevalv.resize(1);
618  fevalv[0] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
619  new FEVALUES (mapping, el, quadrature, flags));
620  }
621  else
622  {
623  fevalv.resize(el.n_base_elements());
624  for (unsigned int i=0; i<fevalv.size(); ++i)
625  {
626  fevalv[i] = std_cxx11::shared_ptr<FEValuesBase<dim,sdim> > (
627  new FEVALUES (mapping, el.base_element(i), quadrature, flags));
628  }
629  }
630  n_components = el.n_components();
631  }
632 
633 
634  template <int dim, int spacedim>
635  inline const FiniteElement<dim, spacedim> &
637  {
638  Assert (fe_pointer !=0, ExcNotInitialized());
639  return *fe_pointer;
640  }
641 
642  template <int dim, int spacedim>
643  inline const FEValuesBase<dim, spacedim> &
645  {
646  AssertDimension(fevalv.size(), 1);
647  return *fevalv[0];
648  }
649 
650 
651  template <int dim, int spacedim>
652  inline const FEValuesBase<dim, spacedim> &
654  {
655  Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
656  return *fevalv[i];
657  }
658 
659 
660  template <int dim, int spacedim>
661  template <typename number>
662  inline void
664  {
665  for (unsigned int i=0; i<fevalv.size(); ++i)
666  {
667  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
668  if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
669  {
670  // This is a subface
671  FESubfaceValues<dim> &fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
672  fe.reinit(info.cell, info.face_number, info.sub_number);
673  }
674  else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
675  {
676  // This is a face
677  FEFaceValues<dim> &fe = dynamic_cast<FEFaceValues<dim>&> (febase);
678  fe.reinit(info.cell, info.face_number);
679  }
680  else
681  {
682  // This is a cell
683  FEValues<dim,spacedim> &fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
684  fe.reinit(info.cell);
685  }
686  }
687 
688  const bool split_fevalues = info.block_info != 0;
689  if (!global_data->empty())
690  fill_local_data(info, split_fevalues);
691  }
692 
693 
694 
695 
696 //----------------------------------------------------------------------//
697 
698  template <int dim, int sdim>
699  inline
700  void
702  unsigned int cp,
703  unsigned int bp,
704  unsigned int fp,
705  bool force)
706  {
707  if (force || cell_quadrature.size() == 0)
708  cell_quadrature = QGauss<dim>(cp);
709  if (force || boundary_quadrature.size() == 0)
710  boundary_quadrature = QGauss<dim-1>(bp);
711  if (force || face_quadrature.size() == 0)
712  face_quadrature = QGauss<dim-1>(fp);
713  }
714 
715 
716  template <int dim, int sdim>
717  inline
718  void
720  {
721  add_update_flags(flags, true, true, true, true);
722  }
723 
724 
725  template <int dim, int sdim>
726  inline
727  void
729  {
730  add_update_flags(flags, true, false, false, false);
731  }
732 
733 
734  template <int dim, int sdim>
735  inline
736  void
738  {
739  add_update_flags(flags, false, true, false, false);
740  }
741 
742 
743  template <int dim, int sdim>
744  inline
745  void
747  {
748  add_update_flags(flags, false, false, true, true);
749  }
750 
751 
752  template <int dim, int sdim>
753  inline
754  void
756  const FiniteElement<dim,sdim> &el,
757  const Mapping<dim,sdim> &mapping,
758  const BlockInfo *block_info)
759  {
760  initialize_update_flags();
761  initialize_gauss_quadrature(
762  (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
763  (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
764  (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
765 
766  cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
767  cell_flags, block_info);
768  boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
769  boundary_flags, block_info);
770  face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
771  face_flags, block_info);
772  subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
773  face_flags, block_info);
774  neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
775  neighbor_flags, block_info);
776  }
777 
778 
779  template <int dim, int sdim>
780  template <typename VECTOR>
781  void
783  const FiniteElement<dim,sdim> &el,
784  const Mapping<dim,sdim> &mapping,
785  const AnyData &data,
786  const VECTOR &dummy,
787  const BlockInfo *block_info)
788  {
789  initialize(el, mapping, block_info);
790  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
792 
793  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
794  // Public member function of parent class was not found without
795  // explicit cast
796  pp = &*p;
797  pp->initialize(data);
798  cell_data = p;
799  cell.initialize_data(p);
800 
801  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
802  pp = &*p;
803  pp->initialize(data);
804  boundary_data = p;
805  boundary.initialize_data(p);
806 
807  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
808  pp = &*p;
809  pp->initialize(data);
810  face_data = p;
811  face.initialize_data(p);
812  subface.initialize_data(p);
813  neighbor.initialize_data(p);
814  }
815 
816 
817  template <int dim, int sdim>
818  template <typename VECTOR>
819  void
821  const FiniteElement<dim,sdim> &el,
822  const Mapping<dim,sdim> &mapping,
823  const NamedData<VECTOR *> &data,
824  const BlockInfo *block_info)
825  {
826  initialize(el, mapping, block_info);
827  std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
828 
829  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
830  p->initialize(data);
831  cell_data = p;
832  cell.initialize_data(p);
833 
834  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
835  p->initialize(data);
836  boundary_data = p;
837  boundary.initialize_data(p);
838 
839  p = std_cxx11::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
840  p->initialize(data);
841  face_data = p;
842  face.initialize_data(p);
843  subface.initialize_data(p);
844  neighbor.initialize_data(p);
845  }
846 
847 
848  template <int dim, int sdim>
849  template <typename VECTOR>
850  void
852  const FiniteElement<dim,sdim> &el,
853  const Mapping<dim,sdim> &mapping,
854  const NamedData<MGLevelObject<VECTOR>*> &data,
855  const BlockInfo *block_info)
856  {
857  initialize(el, mapping, block_info);
858  std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
859 
860  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (cell_selector));
861  p->initialize(data);
862  cell_data = p;
863  cell.initialize_data(p);
864 
865  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (boundary_selector));
866  p->initialize(data);
867  boundary_data = p;
868  boundary.initialize_data(p);
869 
870  p = std_cxx11::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (face_selector));
871  p->initialize(data);
872  face_data = p;
873  face.initialize_data(p);
874  subface.initialize_data(p);
875  neighbor.initialize_data(p);
876  }
877 
878 
879  template <int dim, int sdim>
880  template <class DOFINFO>
881  void
883  {}
884 
885 
886  template <int dim, int sdim>
887  template <class DOFINFO>
888  void
890  {}
891 
892 
893 }
894 
895 DEAL_II_NAMESPACE_CLOSE
896 
897 #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
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:846
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:291
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)
::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 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
::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
::ExceptionBase & ExcInternalError()
std::size_t memory_consumption() const
std_cxx11::shared_ptr< VectorDataBase< dim, spacedim > > global_data
MeshWorker::VectorSelector face_selector