From a13905aa301099c9d52e409c01167b1d413c76d6 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 6 Aug 2010 02:56:35 +0000 Subject: [PATCH] Make MeshWorker components work in 1d as well. git-svn-id: https://svn.dealii.org/trunk@21613 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/numerics/mesh_worker_info.h | 296 +++++++++++++++++- .../numerics/mesh_worker_info.templates.h | 82 +---- .../source/numerics/mesh_worker_info.cc | 12 +- .../numerics/mesh_worker_vector_selector.cc | 6 +- 4 files changed, 304 insertions(+), 92 deletions(-) diff --git a/deal.II/deal.II/include/numerics/mesh_worker_info.h b/deal.II/deal.II/include/numerics/mesh_worker_info.h index 3979448325..d9d4efb272 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.h @@ -14,6 +14,7 @@ #define __deal2__mesh_worker_info_h #include +#include #include #include #include @@ -943,6 +944,139 @@ namespace MeshWorker //----------------------------------------------------------------------// + template + inline + IntegrationInfo::IntegrationInfo() + : + fevalv(0), + multigrid(false), + global_data(boost::shared_ptr >(new VectorDataBase)) + {} + + + template + inline + IntegrationInfo::IntegrationInfo(const IntegrationInfo& other) + : + multigrid(other.multigrid), + values(other.values), + gradients(other.gradients), + hessians(other.hessians), + global_data(other.global_data), + n_components(other.n_components) + { + fevalv.resize(other.fevalv.size()); + for (unsigned int i=0;i& p = *other.fevalv[i]; + const FEValues* pc = dynamic_cast*>(&p); + const FEFaceValues* pf = dynamic_cast*>(&p); + const FESubfaceValues* ps = dynamic_cast*>(&p); + + if (pc != 0) + fevalv[i] = boost::shared_ptr > ( + reinterpret_cast*>( + new FEValues (pc->get_mapping(), pc->get_fe(), + pc->get_quadrature(), pc->get_update_flags()))); + else if (pf != 0) + fevalv[i] = boost::shared_ptr > ( + new FEFaceValues (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags())); + else if (ps != 0) + fevalv[i] = boost::shared_ptr > ( + new FESubfaceValues (ps->get_mapping(), ps->get_fe(), ps->get_quadrature(), ps->get_update_flags())); + else + Assert(false, ExcInternalError()); + } + } + + + template <> + inline + IntegrationInfo<1,1>::IntegrationInfo(const IntegrationInfo<1,1>& other) + : + multigrid(other.multigrid), + values(other.values), + gradients(other.gradients), + hessians(other.hessians), + global_data(other.global_data), + n_components(other.n_components) + { + const int dim = 1; + const int sdim = 1; + + fevalv.resize(other.fevalv.size()); + for (unsigned int i=0;i& p = *other.fevalv[i]; + const FEValues* pc = dynamic_cast*>(&p); + const FEFaceValues* pf = dynamic_cast*>(&p); + const FESubfaceValues* ps = dynamic_cast*>(&p); + + if (pc != 0) + fevalv[i] = boost::shared_ptr > ( + reinterpret_cast*>( + new FEValues (pc->get_mapping(), pc->get_fe(), + pc->get_quadrature(), pc->get_update_flags()))); + else if (pf != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset (); + } + else if (ps != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset(); + } + else + Assert(false, ExcInternalError()); + } + } + + + template <> + inline + IntegrationInfo<1,2>::IntegrationInfo(const IntegrationInfo<1,2>& other) + : + multigrid(other.multigrid), + values(other.values), + gradients(other.gradients), + hessians(other.hessians), + global_data(other.global_data), + n_components(other.n_components) + { + const int dim = 1; + const int sdim = 2; + + fevalv.resize(other.fevalv.size()); + for (unsigned int i=0;i& p = *other.fevalv[i]; + const FEValues* pc = dynamic_cast*>(&p); + const FEFaceValues* pf = dynamic_cast*>(&p); + const FESubfaceValues* ps = dynamic_cast*>(&p); + + if (pc != 0) + fevalv[i] = boost::shared_ptr > ( + reinterpret_cast*>( + new FEValues (pc->get_mapping(), pc->get_fe(), + pc->get_quadrature(), pc->get_update_flags()))); + else if (pf != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset(); + } + else if (ps != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset(); + } + else + Assert(false, ExcInternalError()); + } + } + + + template template inline void @@ -1012,7 +1146,41 @@ namespace MeshWorker else { // This is a cell - FEValues& fe = dynamic_cast&> (febase); + FEValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell); + } + } + + const bool split_fevalues = info.block_info != 0; + if (!global_data->empty()) + fill_local_data(info, split_fevalues); + } + + + template <> + inline void + IntegrationInfo<1,1>::reinit(const DoFInfo<1,1>& info) + { + const int dim = 1; + const int spacedim = 1; + + for (unsigned int i=0;i& febase = *fevalv[i]; + if (info.sub_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a subface + Assert (false, ExcImpossibleInDim(1)); + } + else if (info.face_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a face + Assert (false, ExcImpossibleInDim(1)); + } + else + { + // This is a cell + FEValues& fe = dynamic_cast&> (febase); fe.reinit(info.cell); } } @@ -1023,6 +1191,94 @@ namespace MeshWorker } + template <> + inline void + IntegrationInfo<1,2>::reinit(const DoFInfo<1,2>& info) + { + const int dim = 1; + const int spacedim = 2; + + for (unsigned int i=0;i& febase = *fevalv[i]; + if (info.sub_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a subface + Assert (false, ExcImpossibleInDim(1)); + } + else if (info.face_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a face + Assert (false, ExcImpossibleInDim(1)); + } + else + { + // This is a cell + FEValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell); + } + } + + const bool split_fevalues = info.block_info != 0; + if (!global_data->empty()) + fill_local_data(info, split_fevalues); + } + + + template + inline + void + IntegrationInfoBox:: + initialize(const FiniteElement& el, + const Mapping& mapping, + const BlockInfo* block_info) + { + cell.template initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + boundary.template initialize >(el, mapping, boundary_quadrature, + boundary_flags, block_info); + face.template initialize >(el, mapping, face_quadrature, + face_flags, block_info); + subface.template initialize >(el, mapping, face_quadrature, + face_flags, block_info); + neighbor.template initialize >(el, mapping, face_quadrature, + neighbor_flags, block_info); + } + + + template <> + inline + void + IntegrationInfoBox<1,1>:: + initialize(const FiniteElement<1,1>& el, + const Mapping<1,1>& mapping, + const BlockInfo* block_info) + { + const int dim = 1; + const int sdim = 1; + + cell.initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + } + + + + template <> + inline + void + IntegrationInfoBox<1,2>:: + initialize(const FiniteElement<1,2>& el, + const Mapping<1,2>& mapping, + const BlockInfo* block_info) + { + const int dim = 1; + const int sdim = 2; + + cell.initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + } + + //----------------------------------------------------------------------// template @@ -1099,6 +1355,44 @@ namespace MeshWorker void IntegrationInfoBox::post_faces(const DoFInfoBox&) {} + + + template + inline + void + IntegrationInfoBox::initialize_gauss_quadrature( + unsigned int cp, + unsigned int bp, + unsigned int fp) + { + cell_quadrature = QGauss(cp); + boundary_quadrature = QGauss(bp); + face_quadrature = QGauss(fp); + } + + + template <> + inline + void + IntegrationInfoBox<1,1>::initialize_gauss_quadrature( + unsigned int cp, + unsigned int bp, + unsigned int fp) + { + cell_quadrature = QGauss<1>(cp); + } + + + template <> + inline + void + IntegrationInfoBox<1,2>::initialize_gauss_quadrature( + unsigned int cp, + unsigned int bp, + unsigned int fp) + { + cell_quadrature = QGauss<1>(cp); + } } DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h b/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h index 40f12d38b3..bf0dfd398d 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h @@ -76,50 +76,6 @@ namespace MeshWorker //----------------------------------------------------------------------// - template - IntegrationInfo::IntegrationInfo() - : - fevalv(0), - multigrid(false), - global_data(boost::shared_ptr >(new VectorDataBase)) - {} - - - template - IntegrationInfo::IntegrationInfo(const IntegrationInfo& other) - : - multigrid(other.multigrid), - values(other.values), - gradients(other.gradients), - hessians(other.hessians), - global_data(other.global_data), - n_components(other.n_components) - { - fevalv.resize(other.fevalv.size()); - for (unsigned int i=0;i& p = *other.fevalv[i]; - const FEValues* pc = dynamic_cast*>(&p); - const FEFaceValues* pf = dynamic_cast*>(&p); - const FESubfaceValues* ps = dynamic_cast*>(&p); - - if (pc != 0) - fevalv[i] = boost::shared_ptr > ( - reinterpret_cast*>( - new FEValues (pc->get_mapping(), pc->get_fe(), - pc->get_quadrature(), pc->get_update_flags()))); - else if (pf != 0) - fevalv[i] = boost::shared_ptr > ( - new FEFaceValues (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags())); - else if (ps != 0) - fevalv[i] = boost::shared_ptr > ( - new FESubfaceValues (ps->get_mapping(), ps->get_fe(), ps->get_quadrature(), ps->get_update_flags())); - else - Assert(false, ExcInternalError()); - } - } - - template void IntegrationInfo::initialize_data( @@ -195,7 +151,7 @@ namespace MeshWorker for (unsigned int b=0;blocal().size();++b) { const unsigned int fe_no = info.block_info->base_element(b); - const FEValuesBase& fe = this->fe_values(fe_no); + const FEValuesBase& fe = this->fe_values(fe_no); const unsigned int n_comp = fe.get_fe().n_components(); const unsigned int block_start = info.block_info->local().block_start(b); const unsigned int block_size = info.block_info->local().block_size(b); @@ -211,7 +167,7 @@ namespace MeshWorker } else { - const FEValuesBase& fe = this->fe_values(0); + const FEValuesBase& fe = this->fe_values(0); const unsigned int n_comp = fe.get_fe().n_components(); if(info.level_cell) this->global_data->mg_fill(values, gradients, hessians, fe, info.cell->level(), info.indices, @@ -267,40 +223,6 @@ namespace MeshWorker if (face) face_flags |= flags; if (neighbor) neighbor_flags |= flags; } - - - template - void - IntegrationInfoBox::initialize_gauss_quadrature( - unsigned int cp, - unsigned int bp, - unsigned int fp) - { - cell_quadrature = QGauss(cp); - boundary_quadrature = QGauss(bp); - face_quadrature = QGauss(fp); - - } - - - template - void - IntegrationInfoBox:: - initialize(const FiniteElement& el, - const Mapping& mapping, - const BlockInfo* block_info) - { - cell.template initialize >(el, mapping, cell_quadrature, - cell_flags, block_info); - boundary.template initialize >(el, mapping, boundary_quadrature, - boundary_flags, block_info); - face.template initialize >(el, mapping, face_quadrature, - face_flags, block_info); - subface.template initialize >(el, mapping, face_quadrature, - face_flags, block_info); - neighbor.template initialize >(el, mapping, face_quadrature, - neighbor_flags, block_info); - } } diff --git a/deal.II/deal.II/source/numerics/mesh_worker_info.cc b/deal.II/deal.II/source/numerics/mesh_worker_info.cc index 3d6d8f6d83..44fe4ad0c8 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker_info.cc +++ b/deal.II/deal.II/source/numerics/mesh_worker_info.cc @@ -22,27 +22,26 @@ DEAL_II_NAMESPACE_OPEN -#if deal_II_dimension > 1 namespace MeshWorker { template class IntegrationInfo; template class IntegrationInfoBox; - + template class DoFInfo; template class DoFInfoBox >; - + template void IntegrationInfo::fill_local_data( const DoFInfo&, bool); - + template class DoFInfo; template class DoFInfoBox >; - + template void IntegrationInfo::fill_local_data( const DoFInfo&, bool); - + // template void IntegrationInfo // ::initialize >( // const FiniteElement&, const Mapping&, @@ -57,7 +56,6 @@ namespace MeshWorker // const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); } -#endif DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.cc b/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.cc index 41d651e01b..90c2b12642 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.cc +++ b/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2009 by the deal.II authors +// Copyright (C) 2009, 2010 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -22,8 +22,6 @@ DEAL_II_NAMESPACE_OPEN -#if deal_II_dimension > 1 - namespace MeshWorker { template class VectorDataBase; @@ -31,6 +29,6 @@ namespace MeshWorker #include "mesh_worker_vector_selector.inst" } -#endif + DEAL_II_NAMESPACE_CLOSE -- 2.39.5