00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #ifndef __deal2__mesh_worker_integration_info_h
00014 #define __deal2__mesh_worker_integration_info_h
00015
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/quadrature_lib.h>
00018 #include <deal.II/base/std_cxx1x/shared_ptr.h>
00019 #include <deal.II/dofs/block_info.h>
00020 #include <deal.II/fe/fe_values.h>
00021 #include <deal.II/meshworker/local_results.h>
00022 #include <deal.II/meshworker/dof_info.h>
00023 #include <deal.II/meshworker/vector_selector.h>
00024
00025 DEAL_II_NAMESPACE_OPEN
00026
00027 namespace MeshWorker
00028 {
00079 template<int dim, int spacedim = dim>
00080 class IntegrationInfo
00081 {
00082 private:
00084 std::vector<std_cxx1x::shared_ptr<FEValuesBase<dim, spacedim> > > fevalv;
00085 public:
00086 static const unsigned int dimension = dim;
00087 static const unsigned int space_dimension = spacedim;
00088
00092 IntegrationInfo();
00093
00099 IntegrationInfo(const IntegrationInfo<dim, spacedim>& other);
00100
00134 template <class FEVALUES>
00135 void initialize(const FiniteElement<dim,spacedim>& el,
00136 const Mapping<dim,spacedim>& mapping,
00137 const Quadrature<FEVALUES::integral_dimension>& quadrature,
00138 const UpdateFlags flags,
00139 const BlockInfo* local_block_info = 0);
00140
00146 void initialize_data(const std_cxx1x::shared_ptr<VectorDataBase<dim,spacedim> > &data);
00147
00151 void clear();
00152
00154 bool multigrid;
00156
00165 const FEValuesBase<dim, spacedim>& fe_values () const;
00166
00168
00176 const FEValuesBase<dim, spacedim>& fe_values (const unsigned int i) const;
00177
00192 std::vector<std::vector<std::vector<double> > > values;
00193
00208 std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
00209
00224 std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
00225
00230 void reinit(const DoFInfo<dim, spacedim>& i);
00231
00239 template<typename number>
00240 void fill_local_data(const DoFInfo<dim, spacedim, number>& info, bool split_fevalues);
00241
00248 std_cxx1x::shared_ptr<VectorDataBase<dim, spacedim> > global_data;
00249
00253 std::size_t memory_consumption () const;
00254
00255 private:
00265 template <typename TYPE>
00266 void fill_local_data(
00267 std::vector<std::vector<std::vector<TYPE> > >& data,
00268 VectorSelector& selector,
00269 bool split_fevalues) const;
00274 unsigned int n_components;
00275 };
00276
00333 template <int dim, int spacedim=dim>
00334 class IntegrationInfoBox
00335 {
00336 public:
00337
00342 typedef IntegrationInfo<dim, spacedim> CellInfo;
00343
00347 IntegrationInfoBox ();
00348
00362 void initialize(const FiniteElement<dim, spacedim>& el,
00363 const Mapping<dim, spacedim>& mapping,
00364 const BlockInfo* block_info = 0);
00365
00379 template <typename VECTOR>
00380 void initialize(const FiniteElement<dim, spacedim>& el,
00381 const Mapping<dim, spacedim>& mapping,
00382 const NamedData<VECTOR*>& data,
00383 const BlockInfo* block_info = 0);
00384
00398 template <typename VECTOR>
00399 void initialize(const FiniteElement<dim, spacedim>& el,
00400 const Mapping<dim, spacedim>& mapping,
00401 const NamedData<MGLevelObject<VECTOR>*>& data,
00402 const BlockInfo* block_info = 0);
00406
00407
00428 void initialize_update_flags(bool neighbor_geometry = false);
00429
00436 void add_update_flags_all (const UpdateFlags flags);
00437
00442 void add_update_flags_cell(const UpdateFlags flags);
00443
00448 void add_update_flags_boundary(const UpdateFlags flags);
00449
00454 void add_update_flags_face(const UpdateFlags flags);
00455
00468 void add_update_flags(const UpdateFlags flags,
00469 const bool cell = true,
00470 const bool boundary = true,
00471 const bool face = true,
00472 const bool neighbor = true);
00473
00490 void initialize_gauss_quadrature(unsigned int n_cell_points,
00491 unsigned int n_boundary_points,
00492 unsigned int n_face_points,
00493 const bool force = true);
00494
00498 std::size_t memory_consumption () const;
00499
00507 UpdateFlags cell_flags;
00516 UpdateFlags boundary_flags;
00517
00526 UpdateFlags face_flags;
00527
00537 UpdateFlags neighbor_flags;
00538
00543 Quadrature<dim> cell_quadrature;
00544
00549 Quadrature<dim-1> boundary_quadrature;
00550
00555 Quadrature<dim-1> face_quadrature;
00556
00557
00561
00562
00587 MeshWorker::VectorSelector cell_selector;
00588
00596 MeshWorker::VectorSelector boundary_selector;
00597
00605 MeshWorker::VectorSelector face_selector;
00606
00607 std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
00608 std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
00609 std_cxx1x::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
00610
00611
00615
00648 template <class DOFINFO>
00649 void post_cell(const DoFInfoBox<dim, DOFINFO>&);
00650
00681 template <class DOFINFO>
00682 void post_faces(const DoFInfoBox<dim, DOFINFO>&);
00683
00687 CellInfo cell;
00692 CellInfo boundary;
00698 CellInfo face;
00705 CellInfo subface;
00711 CellInfo neighbor;
00712
00713
00714 };
00715
00716
00717
00718
00719 template<int dim, int sdim>
00720 inline
00721 IntegrationInfo<dim,sdim>::IntegrationInfo()
00722 :
00723 fevalv(0),
00724 multigrid(false),
00725 global_data(std_cxx1x::shared_ptr<VectorDataBase<dim, sdim> >(new VectorDataBase<dim, sdim>))
00726 {}
00727
00728
00729 template<int dim, int sdim>
00730 inline
00731 IntegrationInfo<dim,sdim>::IntegrationInfo(const IntegrationInfo<dim,sdim>& other)
00732 :
00733 multigrid(other.multigrid),
00734 values(other.values),
00735 gradients(other.gradients),
00736 hessians(other.hessians),
00737 global_data(other.global_data),
00738 n_components(other.n_components)
00739 {
00740 fevalv.resize(other.fevalv.size());
00741 for (unsigned int i=0;i<other.fevalv.size();++i)
00742 {
00743 const FEValuesBase<dim,sdim>& p = *other.fevalv[i];
00744 const FEValues<dim,sdim>* pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
00745 const FEFaceValues<dim,sdim>* pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
00746 const FESubfaceValues<dim,sdim>* ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
00747
00748 if (pc != 0)
00749 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00750 reinterpret_cast<FEFaceValuesBase<dim,sdim>*>(
00751 new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
00752 pc->get_quadrature(), pc->get_update_flags())));
00753 else if (pf != 0)
00754 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00755 new FEFaceValues<dim,sdim> (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags()));
00756 else if (ps != 0)
00757 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00758 new FESubfaceValues<dim,sdim> (ps->get_mapping(), ps->get_fe(), ps->get_quadrature(), ps->get_update_flags()));
00759 else
00760 Assert(false, ExcInternalError());
00761 }
00762 }
00763
00764
00765 template <>
00766 inline
00767 IntegrationInfo<1,1>::IntegrationInfo(const IntegrationInfo<1,1>& other)
00768 :
00769 multigrid(other.multigrid),
00770 values(other.values),
00771 gradients(other.gradients),
00772 hessians(other.hessians),
00773 global_data(other.global_data),
00774 n_components(other.n_components)
00775 {
00776 const int dim = 1;
00777 const int sdim = 1;
00778
00779 fevalv.resize(other.fevalv.size());
00780 for (unsigned int i=0;i<other.fevalv.size();++i)
00781 {
00782 const FEValuesBase<dim,sdim>& p = *other.fevalv[i];
00783 const FEValues<dim,sdim>* pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
00784 const FEFaceValues<dim,sdim>* pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
00785 const FESubfaceValues<dim,sdim>* ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
00786
00787 if (pc != 0)
00788 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00789 reinterpret_cast<FEFaceValuesBase<dim,sdim>*>(
00790 new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
00791 pc->get_quadrature(), pc->get_update_flags())));
00792 else if (pf != 0)
00793 {
00794 Assert (false, ExcImpossibleInDim(1));
00795 fevalv[i].reset ();
00796 }
00797 else if (ps != 0)
00798 {
00799 Assert (false, ExcImpossibleInDim(1));
00800 fevalv[i].reset();
00801 }
00802 else
00803 Assert(false, ExcInternalError());
00804 }
00805 }
00806
00807
00808 template <>
00809 inline
00810 IntegrationInfo<1,2>::IntegrationInfo(const IntegrationInfo<1,2>& other)
00811 :
00812 multigrid(other.multigrid),
00813 values(other.values),
00814 gradients(other.gradients),
00815 hessians(other.hessians),
00816 global_data(other.global_data),
00817 n_components(other.n_components)
00818 {
00819 const int dim = 1;
00820 const int sdim = 2;
00821
00822 fevalv.resize(other.fevalv.size());
00823 for (unsigned int i=0;i<other.fevalv.size();++i)
00824 {
00825 const FEValuesBase<dim,sdim>& p = *other.fevalv[i];
00826 const FEValues<dim,sdim>* pc = dynamic_cast<const FEValues<dim,sdim>*>(&p);
00827 const FEFaceValues<dim,sdim>* pf = dynamic_cast<const FEFaceValues<dim,sdim>*>(&p);
00828 const FESubfaceValues<dim,sdim>* ps = dynamic_cast<const FESubfaceValues<dim,sdim>*>(&p);
00829
00830 if (pc != 0)
00831 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00832 reinterpret_cast<FEFaceValuesBase<dim,sdim>*>(
00833 new FEValues<dim,sdim> (pc->get_mapping(), pc->get_fe(),
00834 pc->get_quadrature(), pc->get_update_flags())));
00835 else if (pf != 0)
00836 {
00837 Assert (false, ExcImpossibleInDim(1));
00838 fevalv[i].reset();
00839 }
00840 else if (ps != 0)
00841 {
00842 Assert (false, ExcImpossibleInDim(1));
00843 fevalv[i].reset();
00844 }
00845 else
00846 Assert(false, ExcInternalError());
00847 }
00848 }
00849
00850
00851
00852 template<int dim, int sdim>
00853 template <class FEVALUES>
00854 inline void
00855 IntegrationInfo<dim,sdim>::initialize(
00856 const FiniteElement<dim,sdim>& el,
00857 const Mapping<dim,sdim>& mapping,
00858 const Quadrature<FEVALUES::integral_dimension>& quadrature,
00859 const UpdateFlags flags,
00860 const BlockInfo* block_info)
00861 {
00862 if (block_info == 0 || block_info->local().size() == 0)
00863 {
00864 fevalv.resize(1);
00865 fevalv[0] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00866 new FEVALUES (mapping, el, quadrature, flags));
00867 }
00868 else
00869 {
00870 fevalv.resize(el.n_base_elements());
00871 for (unsigned int i=0;i<fevalv.size();++i)
00872 {
00873 fevalv[i] = std_cxx1x::shared_ptr<FEValuesBase<dim,sdim> > (
00874 new FEVALUES (mapping, el.base_element(i), quadrature, flags));
00875 }
00876 }
00877 n_components = el.n_components();
00878 }
00879
00880
00881 template <int dim, int spacedim>
00882 inline const FEValuesBase<dim, spacedim>&
00883 IntegrationInfo<dim,spacedim>::fe_values() const
00884 {
00885 AssertDimension(fevalv.size(), 1);
00886 return *fevalv[0];
00887 }
00888
00889
00890 template <int dim, int spacedim>
00891 inline const FEValuesBase<dim, spacedim>&
00892 IntegrationInfo<dim,spacedim>::fe_values(unsigned int i) const
00893 {
00894 Assert (i<fevalv.size(), ExcIndexRange(i,0,fevalv.size()));
00895 return *fevalv[i];
00896 }
00897
00898
00899 template <int dim, int spacedim>
00900 inline void
00901 IntegrationInfo<dim,spacedim>::reinit(const DoFInfo<dim, spacedim>& info)
00902 {
00903 for (unsigned int i=0;i<fevalv.size();++i)
00904 {
00905 FEValuesBase<dim, spacedim>& febase = *fevalv[i];
00906 if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
00907 {
00908
00909 FESubfaceValues<dim>& fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
00910 fe.reinit(info.cell, info.face_number, info.sub_number);
00911 }
00912 else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
00913 {
00914
00915 FEFaceValues<dim>& fe = dynamic_cast<FEFaceValues<dim>&> (febase);
00916 fe.reinit(info.cell, info.face_number);
00917 }
00918 else
00919 {
00920
00921 FEValues<dim,spacedim>& fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
00922 fe.reinit(info.cell);
00923 }
00924 }
00925
00926 const bool split_fevalues = info.block_info != 0;
00927 if (!global_data->empty())
00928 fill_local_data(info, split_fevalues);
00929 }
00930
00931
00932 template <>
00933 inline void
00934 IntegrationInfo<1,1>::reinit(const DoFInfo<1,1>& info)
00935 {
00936 const int dim = 1;
00937 const int spacedim = 1;
00938
00939 for (unsigned int i=0;i<fevalv.size();++i)
00940 {
00941 FEValuesBase<dim, spacedim>& febase = *fevalv[i];
00942 if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
00943 {
00944
00945 Assert (false, ExcImpossibleInDim(1));
00946 }
00947 else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
00948 {
00949
00950 Assert (false, ExcImpossibleInDim(1));
00951 }
00952 else
00953 {
00954
00955 FEValues<dim,spacedim>& fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
00956 fe.reinit(info.cell);
00957 }
00958 }
00959
00960 const bool split_fevalues = info.block_info != 0;
00961 if (!global_data->empty())
00962 fill_local_data(info, split_fevalues);
00963 }
00964
00965
00966 template <>
00967 inline void
00968 IntegrationInfo<1,2>::reinit(const DoFInfo<1,2>& info)
00969 {
00970 const int dim = 1;
00971 const int spacedim = 2;
00972
00973 for (unsigned int i=0;i<fevalv.size();++i)
00974 {
00975 FEValuesBase<dim, spacedim>& febase = *fevalv[i];
00976 if (info.sub_number != deal_II_numbers::invalid_unsigned_int)
00977 {
00978
00979 Assert (false, ExcImpossibleInDim(1));
00980 }
00981 else if (info.face_number != deal_II_numbers::invalid_unsigned_int)
00982 {
00983
00984 Assert (false, ExcImpossibleInDim(1));
00985 }
00986 else
00987 {
00988
00989 FEValues<dim,spacedim>& fe = dynamic_cast<FEValues<dim,spacedim>&> (febase);
00990 fe.reinit(info.cell);
00991 }
00992 }
00993
00994 const bool split_fevalues = info.block_info != 0;
00995 if (!global_data->empty())
00996 fill_local_data(info, split_fevalues);
00997 }
00998
00999
01000
01001
01002 template <>
01003 inline
01004 void
01005 IntegrationInfoBox<1,1>::initialize_gauss_quadrature(
01006 const unsigned int cp,
01007 const unsigned int,
01008 const unsigned int,
01009 bool force)
01010 {
01011 if (force || cell_quadrature.size() == 0)
01012 cell_quadrature = QGauss<1>(cp);
01013 }
01014
01015
01016 template <>
01017 inline
01018 void
01019 IntegrationInfoBox<1,2>::initialize_gauss_quadrature(
01020 const unsigned int cp,
01021 const unsigned int,
01022 const unsigned int,
01023 bool force)
01024 {
01025 if (force || cell_quadrature.size() == 0)
01026 cell_quadrature = QGauss<1>(cp);
01027 }
01028
01029
01030 template <int dim, int sdim>
01031 inline
01032 void
01033 IntegrationInfoBox<dim,sdim>::initialize_gauss_quadrature(
01034 unsigned int cp,
01035 unsigned int bp,
01036 unsigned int fp,
01037 bool force)
01038 {
01039 if (force || cell_quadrature.size() == 0)
01040 cell_quadrature = QGauss<dim>(cp);
01041 if (force || boundary_quadrature.size() == 0)
01042 boundary_quadrature = QGauss<dim-1>(bp);
01043 if (force || face_quadrature.size() == 0)
01044 face_quadrature = QGauss<dim-1>(fp);
01045 }
01046
01047
01048 template <int dim, int sdim>
01049 inline
01050 void
01051 IntegrationInfoBox<dim,sdim>::add_update_flags_all (const UpdateFlags flags)
01052 {
01053 add_update_flags(flags, true, true, true, true);
01054 }
01055
01056
01057 template <int dim, int sdim>
01058 inline
01059 void
01060 IntegrationInfoBox<dim,sdim>::add_update_flags_cell (const UpdateFlags flags)
01061 {
01062 add_update_flags(flags, true, false, false, false);
01063 }
01064
01065
01066 template <int dim, int sdim>
01067 inline
01068 void
01069 IntegrationInfoBox<dim,sdim>::add_update_flags_boundary (const UpdateFlags flags)
01070 {
01071 add_update_flags(flags, false, true, false, false);
01072 }
01073
01074
01075 template <int dim, int sdim>
01076 inline
01077 void
01078 IntegrationInfoBox<dim,sdim>::add_update_flags_face (const UpdateFlags flags)
01079 {
01080 add_update_flags(flags, false, false, true, true);
01081 }
01082
01083
01084 template <>
01085 inline
01086 void
01087 IntegrationInfoBox<1,1>::initialize(
01088 const FiniteElement<1,1>& el,
01089 const Mapping<1,1>& mapping,
01090 const BlockInfo* block_info)
01091 {
01092 initialize_update_flags();
01093 initialize_gauss_quadrature(
01094 (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), 1, 1, false);
01095
01096 const int dim = 1;
01097 const int sdim = 1;
01098
01099 cell.initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
01100 cell_flags, block_info);
01101 }
01102
01103
01104
01105 template <>
01106 inline
01107 void
01108 IntegrationInfoBox<1,2>::initialize(
01109 const FiniteElement<1,2>& el,
01110 const Mapping<1,2>& mapping,
01111 const BlockInfo* block_info)
01112 {
01113 initialize_update_flags();
01114 initialize_gauss_quadrature(
01115 (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), 1, 1, false);
01116
01117 const int dim = 1;
01118 const int sdim = 2;
01119
01120 cell.initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
01121 cell_flags, block_info);
01122 }
01123
01124
01125 template <int dim, int sdim>
01126 inline
01127 void
01128 IntegrationInfoBox<dim,sdim>::initialize(
01129 const FiniteElement<dim,sdim>& el,
01130 const Mapping<dim,sdim>& mapping,
01131 const BlockInfo* block_info)
01132 {
01133 initialize_update_flags();
01134 initialize_gauss_quadrature(
01135 (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
01136 (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(),
01137 (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false);
01138
01139 cell.template initialize<FEValues<dim,sdim> >(el, mapping, cell_quadrature,
01140 cell_flags, block_info);
01141 boundary.template initialize<FEFaceValues<dim,sdim> >(el, mapping, boundary_quadrature,
01142 boundary_flags, block_info);
01143 face.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
01144 face_flags, block_info);
01145 subface.template initialize<FESubfaceValues<dim,sdim> >(el, mapping, face_quadrature,
01146 face_flags, block_info);
01147 neighbor.template initialize<FEFaceValues<dim,sdim> >(el, mapping, face_quadrature,
01148 neighbor_flags, block_info);
01149 }
01150
01151
01152 template <int dim, int sdim>
01153 template <typename VECTOR>
01154 void
01155 IntegrationInfoBox<dim,sdim>::initialize(
01156 const FiniteElement<dim,sdim>& el,
01157 const Mapping<dim,sdim>& mapping,
01158 const NamedData<VECTOR*>& data,
01159 const BlockInfo* block_info)
01160 {
01161 initialize(el, mapping, block_info);
01162 std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
01163
01164 p = std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (cell_selector));
01165 p->initialize(data);
01166 cell_data = p;
01167 cell.initialize_data(p);
01168
01169 p = std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (boundary_selector));
01170 p->initialize(data);
01171 boundary_data = p;
01172 boundary.initialize_data(p);
01173
01174 p = std_cxx1x::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (face_selector));
01175 p->initialize(data);
01176 face_data = p;
01177 face.initialize_data(p);
01178 subface.initialize_data(p);
01179 neighbor.initialize_data(p);
01180 }
01181
01182
01183 template <int dim, int sdim>
01184 template <typename VECTOR>
01185 void
01186 IntegrationInfoBox<dim,sdim>::initialize(
01187 const FiniteElement<dim,sdim>& el,
01188 const Mapping<dim,sdim>& mapping,
01189 const NamedData<MGLevelObject<VECTOR>*>& data,
01190 const BlockInfo* block_info)
01191 {
01192 initialize(el, mapping, block_info);
01193 std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> > p;
01194
01195 p = std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (cell_selector));
01196 p->initialize(data);
01197 cell_data = p;
01198 cell.initialize_data(p);
01199
01200 p = std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (boundary_selector));
01201 p->initialize(data);
01202 boundary_data = p;
01203 boundary.initialize_data(p);
01204
01205 p = std_cxx1x::shared_ptr<MGVectorData<VECTOR, dim, sdim> >(new MGVectorData<VECTOR, dim, sdim> (face_selector));
01206 p->initialize(data);
01207 face_data = p;
01208 face.initialize_data(p);
01209 subface.initialize_data(p);
01210 neighbor.initialize_data(p);
01211 }
01212
01213
01214 template <int dim, int sdim>
01215 template <class DOFINFO>
01216 void
01217 IntegrationInfoBox<dim,sdim>::post_cell(const DoFInfoBox<dim, DOFINFO>&)
01218 {}
01219
01220
01221 template <int dim, int sdim>
01222 template <class DOFINFO>
01223 void
01224 IntegrationInfoBox<dim,sdim>::post_faces(const DoFInfoBox<dim, DOFINFO>&)
01225 {}
01226
01227
01228 }
01229
01230 DEAL_II_NAMESPACE_CLOSE
01231
01232 #endif
01233