Reference documentation for deal.II version SVN Revision 32792
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
integration_info.h
1 // ---------------------------------------------------------------------
2 // @f$Id: integration_info.h 32162 2014-01-03 18:03:48Z heister @f$
3 //
4 // Copyright (C) 2006 - 2013 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_cxx1x/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_cxx1x::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_cxx1x::shared_ptr<VectorDataBase<dim,spacedim> > &data);
129 
133  void clear();
134 
136  bool multigrid;
138 
143  const FEValuesBase<dim, spacedim> &fe_values () const;
144 
146 
152  const FEValuesBase<dim, spacedim> &fe_values (const unsigned int i) const;
153 
162  std::vector<std::vector<std::vector<double> > > values;
163 
172  std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
173 
182  std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
183 
187  template <typename number>
188  void reinit(const DoFInfo<dim, spacedim, number> &i);
189 
194  template<typename number>
195  void fill_local_data(const DoFInfo<dim, spacedim, number> &info, bool split_fevalues);
196 
201  std_cxx1x::shared_ptr<VectorDataBase<dim, spacedim> > global_data;
202 
206  std::size_t memory_consumption () const;
207 
208  private:
214  template <typename TYPE>
215  void fill_local_data(
216  std::vector<std::vector<std::vector<TYPE> > > &data,
217  VectorSelector &selector,
218  bool split_fevalues) const;
222  unsigned int n_components;
223  };
224 
281  template <int dim, int spacedim=dim>
283  {
284  public:
285 
291 
296 
305  const Mapping<dim, spacedim> &mapping,
306  const BlockInfo *block_info = 0);
307 
315  template <typename VECTOR>
317  const Mapping<dim, spacedim> &mapping,
318  const NamedData<VECTOR *> &data,
319  const BlockInfo *block_info = 0);
320 
328  template <typename VECTOR>
330  const Mapping<dim, spacedim> &mapping,
331  const NamedData<MGLevelObject<VECTOR>*> &data,
332  const BlockInfo *block_info = 0);
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 
466 
472 
478 
479  std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
480  std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
481  std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
482  /* @} */
483 
487  /* @{ */
506  template <class DOFINFO>
507  void post_cell(const DoFInfoBox<dim, DOFINFO> &);
508 
526  template <class DOFINFO>
527  void post_faces(const DoFInfoBox<dim, DOFINFO> &);
528 
551 
552  /* @} */
553  };
554 
555 
556 //----------------------------------------------------------------------//
557 
558  template<int dim, int sdim>
559  inline
561  :
562  fevalv(0),
563  multigrid(false),
564  global_data(std_cxx1x::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>))
565  {}
566 
567 
568  template<int dim, int sdim>
569  inline
571  :
572  multigrid(other.multigrid),
573  values(other.values),
574  gradients(other.gradients),
575  hessians(other.hessians),
576  global_data(other.global_data),
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_cxx1x::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_cxx1x::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_cxx1x::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  if (block_info == 0 || block_info->local().size() == 0)
615  {
616  fevalv.resize(1);
617  fevalv[0] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
618  new FEVALUES (mapping, el, quadrature, flags));
619  }
620  else
621  {
622  fevalv.resize(el.n_base_elements());
623  for (unsigned int i=0; i<fevalv.size(); ++i)
624  {
625  fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
626  new FEVALUES (mapping, el.base_element(i), quadrature, flags));
627  }
628  }
629  n_components = el.n_components();
630  }
631 
632 
633  template <int dim, int spacedim>
634  inline const FEValuesBase<dim, spacedim> &
636  {
637  AssertDimension(fevalv.size(), 1);
638  return *fevalv[0];
639  }
640 
641 
642  template <int dim, int spacedim>
643  inline const FEValuesBase<dim, spacedim> &
645  {
646  Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
647  return *fevalv[i];
648  }
649 
650 
651  template <int dim, int spacedim>
652  template <typename number>
653  inline void
655  {
656  for (unsigned int i=0; i<fevalv.size(); ++i)
657  {
658  FEValuesBase<dim, spacedim> &febase = *fevalv[i];
659  if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
660  {
661  // This is a subface
662  FESubfaceValues<dim> &fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
663  fe.reinit(info.cell, info.face_number, info.sub_number);
664  }
665  else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
666  {
667  // This is a face
668  FEFaceValues<dim> &fe = dynamic_cast<FEFaceValues<dim>&> (febase);
669  fe.reinit(info.cell, info.face_number);
670  }
671  else
672  {
673  // This is a cell
674  FEValues<dim,spacedim> &fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
675  fe.reinit(info.cell);
676  }
677  }
678 
679  const bool split_fevalues = info.block_info != 0;
680  if (!global_data->empty())
681  fill_local_data(info, split_fevalues);
682  }
683 
684 
685 
686 
687 //----------------------------------------------------------------------//
688 
689  template <int dim, int sdim>
690  inline
691  void
693  unsigned int cp,
694  unsigned int bp,
695  unsigned int fp,
696  bool force)
697  {
698  if (force || cell_quadrature.size() == 0)
699  cell_quadrature = QGauss<dim>(cp);
700  if (force || boundary_quadrature.size() == 0)
701  boundary_quadrature = QGauss<dim-1>(bp);
702  if (force || face_quadrature.size() == 0)
703  face_quadrature = QGauss<dim-1>(fp);
704  }
705 
706 
707  template <int dim, int sdim>
708  inline
709  void
711  {
712  add_update_flags(flags, true, true, true, true);
713  }
714 
715 
716  template <int dim, int sdim>
717  inline
718  void
720  {
721  add_update_flags(flags, true, false, false, false);
722  }
723 
724 
725  template <int dim, int sdim>
726  inline
727  void
729  {
730  add_update_flags(flags, false, true, false, false);
731  }
732 
733 
734  template <int dim, int sdim>
735  inline
736  void
738  {
739  add_update_flags(flags, false, false, true, true);
740  }
741 
742 
743  template <int dim, int sdim>
744  inline
745  void
747  const FiniteElement<dim,sdim> &el,
748  const Mapping<dim,sdim> &mapping,
749  const BlockInfo *block_info)
750  {
751  initialize_update_flags();
752  initialize_gauss_quadrature(
753  (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
754  (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
755  (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
756 
757  cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
758  cell_flags, block_info);
759  boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
760  boundary_flags, block_info);
761  face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
762  face_flags, block_info);
763  subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
764  face_flags, block_info);
765  neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
766  neighbor_flags, block_info);
767  }
768 
769 
770  template <int dim, int sdim>
771  template <typename VECTOR>
772  void
774  const FiniteElement<dim,sdim> &el,
775  const Mapping<dim,sdim> &mapping,
776  const NamedData<VECTOR *> &data,
777  const BlockInfo *block_info)
778  {
779  initialize(el, mapping, block_info);
780  std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
781 
782  p = std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
783  p->initialize(data);
784  cell_data = p;
785  cell.initialize_data(p);
786 
787  p = std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
788  p->initialize(data);
789  boundary_data = p;
790  boundary.initialize_data(p);
791 
792  p = std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
793  p->initialize(data);
794  face_data = p;
795  face.initialize_data(p);
796  subface.initialize_data(p);
797  neighbor.initialize_data(p);
798  }
799 
800 
801  template <int dim, int sdim>
802  template <typename VECTOR>
803  void
805  const FiniteElement<dim,sdim> &el,
806  const Mapping<dim,sdim> &mapping,
807  const NamedData<MGLevelObject<VECTOR>*> &data,
808  const BlockInfo *block_info)
809  {
810  initialize(el, mapping, block_info);
811  std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
812 
813  p = std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (cell_selector));
814  p->initialize(data);
815  cell_data = p;
816  cell.initialize_data(p);
817 
818  p = std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (boundary_selector));
819  p->initialize(data);
820  boundary_data = p;
821  boundary.initialize_data(p);
822 
823  p = std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (face_selector));
824  p->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 <class DOFINFO>
834  void
836  {}
837 
838 
839  template <int dim, int sdim>
840  template <class DOFINFO>
841  void
843  {}
844 
845 
846 }
847 
848 DEAL_II_NAMESPACE_CLOSE
849 
850 #endif
void reinit(const TriaIterator< DoFCellAccessor< DH, level_dof_access > > cell, const unsigned int face_no)
void add_update_flags(const UpdateFlags flags, const bool cell=true, const bool boundary=true, const bool face=true, const bool neighbor=true)
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
#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
const Mapping< dim, spacedim > & get_mapping() const
void add_update_flags_all(const UpdateFlags flags)
unsigned int sub_number
Definition: dof_info.h:99
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)
void add_update_flags_face(const UpdateFlags flags)
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
const FiniteElement< dim, spacedim > & get_fe() const
MeshWorker::VectorSelector cell_selector
#define Assert(cond, exc)
Definition: exceptions.h:299
UpdateFlags
void post_cell(const DoFInfoBox< dim, DOFINFO > &)
UpdateFlags get_update_flags() const
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:173
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void initialize_update_flags(bool neighbor_geometry=false)
unsigned int n_base_elements() const
unsigned int n_components() const
void initialize_data(const std_cxx1x::shared_ptr< VectorDataBase< dim, spacedim > > &data)
void add_update_flags_boundary(const UpdateFlags flags)
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) 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.
unsigned int face_number
Definition: dof_info.h:88
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
unsigned int size() const
IntegrationInfo< dim, spacedim > CellInfo
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::vector< std_cxx1x::shared_ptr< FEValuesBase< dim, spacedim > > > fevalv
vector of FEValues objects
MeshWorker::VectorSelector face_selector
std_cxx1x::shared_ptr< VectorDataBase< dim, spacedim > > global_data
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)