From: kanschat Date: Tue, 16 Feb 2010 04:00:16 +0000 (+0000) Subject: introduce DoFInfoBox X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8290b7f36d829b2064bc6e31b92489ca9cc05f5f;p=dealii-svn.git introduce DoFInfoBox git-svn-id: https://svn.dealii.org/trunk@20645 0785d39b-7218-0410-832d-ea1e28bc413d --- 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 40d9445f59..8be76611c1 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.h @@ -193,6 +193,7 @@ namespace MeshWorker std::vector > > M2; }; + template class DoFInfoBox; /** * A class containing information on geometry and degrees of freedom @@ -277,8 +278,8 @@ namespace MeshWorker * the #aux_local_indices. */ template - DoFInfo(const DH& dof_handler); - + DoFInfo (const DH& dof_handler); + /** * Set the current cell and * fill #indices. @@ -314,7 +315,18 @@ namespace MeshWorker SmartPointer > block_info; private: - /// Fill index vector with active indices + /** + * Standard constructor, not + * setting any block + * indices. Use of this + * constructor is not + * recommended, but it is + * needed for the arrays in + * DoFInfoBox. + */ + DoFInfo (); + + /// Fill index vector with active indices void get_indices(const typename DoFHandler::cell_iterator& c); /// Fill index vector with level indices @@ -332,8 +344,66 @@ namespace MeshWorker * degrees of freedom per cell. */ BlockIndices aux_local_indices; + + friend class DoFInfoBox; + }; + + + /** + * A class bundling the MeshWorker::DoFInfo objects used on a cell. + * + * @todo Currently, we are storing an object for the cells and two for + * each face. We could gather all face data pertaining to the cell + * itself in one object, saving a bit of memory and a few operations, + * but sacrificing some cleanliness. + * + * @author Guido Kanschat, 2010 + */ + template < int dim, int spacedim=dim> + class DoFInfoBox + { + public: + /** + * Constructor copying the seed + * into all other objects. + */ + DoFInfoBox(const MeshWorker::DoFInfo& seed); + + /** + * Reset all the availability flags. + */ + void reset(); + + /** + * The data for the cell. + */ + MeshWorker::DoFInfo cell; + /** + * The data for the faces from inside. + */ + MeshWorker::DoFInfo interior[GeometryInfo::faces_per_cell]; + /** + * The data for the faces from outside. + */ + MeshWorker::DoFInfo exterior[GeometryInfo::faces_per_cell]; + + /** + * A set of flags, indicating + * whether data on an interior + * face is available. + */ + bool interior_face_available[GeometryInfo::faces_per_cell]; + /** + * A set of flags, indicating + * whether data on an exterior + * face is available. + */ + bool exterior_face_available[GeometryInfo::faces_per_cell]; }; + + + /** * Class for objects handed to local integration functions. * @@ -876,7 +946,35 @@ namespace MeshWorker return aux_local_indices; } +//----------------------------------------------------------------------// + + template < int dim, int spacedim> + inline + DoFInfoBox::DoFInfoBox(const DoFInfo& seed) + : + cell(seed) + { + for (unsigned int i=0;i::faces_per_cell;++i) + { + exterior[i] = seed; + interior[i] = seed; + interior_face_available[i] = false; + exterior_face_available[i] = false; + } + } + + template < int dim, int spacedim> + inline void + DoFInfoBox::reset () + { + for (unsigned int i=0;i::faces_per_cell;++i) + { + interior_face_available[i] = false; + exterior_face_available[i] = false; + } + } + //----------------------------------------------------------------------// template 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 82baecbaa8..3aa4d16b03 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 @@ -25,6 +25,11 @@ namespace MeshWorker {} + template + DoFInfo::DoFInfo() + {} + + template void DoFInfo::get_indices(const typename DoFHandler::cell_iterator& c) diff --git a/deal.II/deal.II/include/numerics/mesh_worker_loop.h b/deal.II/deal.II/include/numerics/mesh_worker_loop.h index 02978105f0..f489e07334 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_loop.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_loop.h @@ -68,14 +68,14 @@ namespace MeshWorker * @ingroup MeshWorker * @author Guido Kanschat, 2009 */ - template + template void loop(ITERATOR begin, typename identity::type end, - DOFINFO cell_dof_info, + DoFInfo dinfo, INFOBOX& info, - const std_cxx1x::function &cell_worker, - const std_cxx1x::function &boundary_worker, - const std_cxx1x::function&, typename INFOBOX::CellInfo &)> &cell_worker, + const std_cxx1x::function&, typename INFOBOX::FaceInfo &)> &boundary_worker, + const std_cxx1x::function&, DoFInfo&, typename INFOBOX::FaceInfo &, typename INFOBOX::FaceInfo &)> &face_worker, ASSEMBLER &assembler, @@ -85,12 +85,14 @@ namespace MeshWorker const bool integrate_boundary = (boundary_worker != 0); const bool integrate_interior_face = (face_worker != 0); - ; - DOFINFO face_dof_info = cell_dof_info; - DOFINFO neighbor_dof_info = cell_dof_info; - assembler.initialize_info(cell_dof_info, false); - assembler.initialize_info(face_dof_info, true); - assembler.initialize_info(neighbor_dof_info, true); + DoFInfoBox dof_info(dinfo); + + assembler.initialize_info(dof_info.cell, false); + for (unsigned int i=0;i::faces_per_cell;++i) + { + assembler.initialize_info(dof_info.interior[i], true); + assembler.initialize_info(dof_info.exterior[i], true); + } // Loop over all cells for (ITERATOR cell = begin; cell != end; ++cell) @@ -100,10 +102,10 @@ namespace MeshWorker // before faces if (integrate_cell && cells_first) { - cell_dof_info.reinit(cell); - info.cell.reinit(cell_dof_info); - cell_worker(cell_dof_info, info.cell); - assembler.assemble(cell_dof_info); + dof_info.cell.reinit(cell); + info.cell.reinit(dof_info.cell); + cell_worker(dof_info.cell, info.cell); + assembler.assemble(dof_info.cell); } if (integrate_interior_face || integrate_boundary) @@ -114,10 +116,11 @@ namespace MeshWorker { if (integrate_boundary) { - cell_dof_info.reinit(cell, face, face_no); - info.boundary.reinit(cell_dof_info); - boundary_worker(cell_dof_info, info.boundary); - assembler.assemble(cell_dof_info); + dof_info.interior_face_available[face_no] = true; + dof_info.interior[face_no].reinit(cell, face, face_no); + info.boundary.reinit(dof_info.interior[face_no]); + boundary_worker(dof_info.interior[face_no], info.boundary); + assembler.assemble(dof_info.interior[face_no]); } } else if (integrate_interior_face) @@ -146,17 +149,21 @@ namespace MeshWorker typename ITERATOR::AccessorType::Container::face_iterator nface = neighbor->face(neighbor_face_no.first); - face_dof_info.reinit(cell, face, face_no); - info.face.reinit(face_dof_info); - neighbor_dof_info.reinit(neighbor, nface, - neighbor_face_no.first, neighbor_face_no.second); - info.subface.reinit(neighbor_dof_info); + dof_info.interior_face_available[face_no] = true; + dof_info.exterior_face_available[face_no] = true; + dof_info.interior[face_no].reinit(cell, face, face_no); + info.face.reinit(dof_info.interior[face_no]); + dof_info.exterior[face_no].reinit( + neighbor, nface, neighbor_face_no.first, neighbor_face_no.second); + info.subface.reinit(dof_info.exterior[face_no]); + // Neighbor // first to // conform to // old version - face_worker(face_dof_info, neighbor_dof_info, info.face, info.subface); - assembler.assemble (face_dof_info, neighbor_dof_info); + face_worker(dof_info.interior[face_no], dof_info.exterior[face_no], + info.face, info.subface); + assembler.assemble (dof_info.interior[face_no], dof_info.exterior[face_no]); } else { @@ -179,13 +186,17 @@ namespace MeshWorker unsigned int neighbor_face_no = cell->neighbor_of_neighbor(face_no); Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError()); // Regular interior face - face_dof_info.reinit(cell, face, face_no); - info.face.reinit(face_dof_info); - neighbor_dof_info.reinit(neighbor, neighbor->face(neighbor_face_no), - neighbor_face_no); - info.neighbor.reinit(neighbor_dof_info); - face_worker(face_dof_info, neighbor_dof_info, info.face, info.neighbor); - assembler.assemble(face_dof_info, neighbor_dof_info); + dof_info.interior_face_available[face_no] = true; + dof_info.exterior_face_available[face_no] = true; + dof_info.interior[face_no].reinit(cell, face, face_no); + info.face.reinit(dof_info.interior[face_no]); + dof_info.exterior[face_no].reinit( + neighbor, neighbor->face(neighbor_face_no), neighbor_face_no); + info.neighbor.reinit(dof_info.exterior[face_no]); + + face_worker(dof_info.interior[face_no], dof_info.exterior[face_no], + info.face, info.neighbor); + assembler.assemble(dof_info.interior[face_no], dof_info.exterior[face_no]); } } } // faces @@ -193,10 +204,10 @@ namespace MeshWorker // have to be handled first if (integrate_cell && !cells_first) { - cell_dof_info.reinit(cell); - info.cell.reinit(cell_dof_info); - cell_worker(cell_dof_info, info.cell); - assembler.assemble (cell_dof_info); + dof_info.cell.reinit(cell); + info.cell.reinit(dof_info.cell); + cell_worker(dof_info.cell, info.cell); + assembler.assemble(dof_info.cell); } } } diff --git a/deal.II/examples/step-38/step-38.cc b/deal.II/examples/step-38/step-38.cc index 43411b0cc6..8eff89bddf 100644 --- a/deal.II/examples/step-38/step-38.cc +++ b/deal.II/examples/step-38/step-38.cc @@ -387,7 +387,7 @@ void Step12::assemble_system () // result of std::bind if the local // integrators were, for example, // non-static member functions. - MeshWorker::loop, MeshWorker::IntegrationInfoBox > + MeshWorker::loop, dim, dim> (dof_handler.begin_active(), dof_handler.end(), dof_info, info_box, &Step12::integrate_cell_term, diff --git a/deal.II/examples/step-39/step-39.cc b/deal.II/examples/step-39/step-39.cc index 12f257628d..c2edd88720 100644 --- a/deal.II/examples/step-39/step-39.cc +++ b/deal.II/examples/step-39/step-39.cc @@ -446,7 +446,7 @@ Step39::assemble_matrix() MeshWorker::IntegrationInfoBox info_box; MeshWorker::DoFInfo dof_info(dof_handler); info_box.initialize(integration_worker, fe, mapping); - MeshWorker::loop, MeshWorker::IntegrationInfoBox >( + MeshWorker::loop, dim, dim>( dof_handler.begin_active(), dof_handler.end(), dof_info, info_box, &MatrixIntegrator::cell, @@ -473,7 +473,7 @@ Step39::assemble_mg_matrix() MeshWorker::IntegrationInfoBox info_box; MeshWorker::DoFInfo dof_info(mg_dof_handler); info_box.initialize(integration_worker, fe, mapping); - MeshWorker::loop, MeshWorker::IntegrationInfoBox > ( + MeshWorker::loop, dim, dim> ( mg_dof_handler.begin(), mg_dof_handler.end(), dof_info, info_box, &MatrixIntegrator::cell, @@ -503,7 +503,7 @@ Step39::assemble_right_hand_side() MeshWorker::DoFInfo dof_info(dof_handler); info_box.initialize(integration_worker, fe, mapping); - MeshWorker::loop, MeshWorker::IntegrationInfoBox >( + MeshWorker::loop, dim, dim>( dof_handler.begin_active(), dof_handler.end(), dof_info, info_box, &RHSIntegrator::cell, @@ -529,7 +529,7 @@ Step39::solve() coarse_matrix.copy_from (mg_matrix[0]); MGCoarseGridHouseholder > mg_coarse; mg_coarse.initialize(coarse_matrix); - typedef PreconditionSOR > RELAXATION; + typedef PreconditionSSOR > RELAXATION; MGSmootherRelaxation, RELAXATION, Vector > mg_smoother(mem); RELAXATION::AdditionalData smoother_data(1.); @@ -624,7 +624,7 @@ Step39::estimate() MeshWorker::IntegrationInfoBox info_box; MeshWorker::DoFInfo dof_info(dof_handler); info_box.initialize(integration_worker, fe, mapping, solution_data); - MeshWorker::loop, MeshWorker::IntegrationInfoBox > ( + MeshWorker::loop, dim, dim> ( dof_handler.begin_active(), dof_handler.end(), dof_info, info_box, &Estimator::cell,