From: Timo Heister Date: Sun, 13 Aug 2017 18:17:43 +0000 (-0600) Subject: cleanup assemble flags, fix logic, add asserts, adjust tests X-Git-Tag: v9.0.0-rc1~1271^2~4 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=26b3190f96b0aab145bd0deb11fdd8d626921695;p=dealii.git cleanup assemble flags, fix logic, add asserts, adjust tests --- diff --git a/include/deal.II/meshworker/assemble_flags.h b/include/deal.II/meshworker/assemble_flags.h index 4fcd6e86e0..cea5a1902c 100644 --- a/include/deal.II/meshworker/assemble_flags.h +++ b/include/deal.II/meshworker/assemble_flags.h @@ -40,67 +40,68 @@ namespace MeshWorker */ enum AssembleFlags { - //! No update - assemble_default = 0, - //! Own cells /** - * Assemble on cells. + * Do Nothing. + */ + assemble_nothing = 0, + /** + * Assemble on locally owned cells. */ assemble_own_cells = 0x0001, - //! Ghost cells /** * Assemble on ghost cells. */ assemble_ghost_cells = 0x0002, - //! Own faces once /** - * Assemble on own interior faces, visiting each face only once + * Assemble on interior faces between two locally owned cells, + * visiting each face only once */ assemble_own_interior_faces_once = 0x0004, - //! Own faces both /** - * Assemble on own interior faces, visiting each interior face twice + * Assemble on interior faces between two locally owned cells, + * visiting each interior face twice */ assemble_own_interior_faces_both = 0x0008, - //! Ghost faces once /** * Assemble on faces between a locally owned cell and a ghost cell, making * sure that only one of the processes will assemble these faces (from the * finer side or the process with the lower mpi rank) */ assemble_ghost_faces_once = 0x0010, - //! Ghost faces both /** * Assemble on faces between a locally owned cell and a ghost cell, both * processes will assemble these faces. Note that these faces are never * assembled from both sides on a single process. */ assemble_ghost_faces_both = 0x0020, - //! Assemble cells before faces /** - * Assemble cell integrals before face integrals. + * Assemble on boundary faces of the locally owned cells */ - assemble_cells_first = 0x0040, - //! Assemble boundary faces + assemble_boundary_faces = 0x0040, + /** - * Assemble on boundary faces + * Assemble cell integrals before face integrals. */ - assemble_boundary_faces = 0x0080, - //! Assemble own interior aces + cells_first = 0x0080, + /** - * Assemble own interior faces, either interior ones or on the boundary. + * Combination of flags to determine if any work on cells is done */ - assemble_own_interior_faces = assemble_own_interior_faces_both | assemble_own_interior_faces_once, - //! Assemble own faces + work_on_cells = assemble_own_cells | assemble_ghost_cells, + /** - * Assemble own faces, either interior ones or on the boundary. + * Combination of flags to determine if any work is done on faces */ - assemble_own_faces = assemble_own_interior_faces | assemble_boundary_faces, - //! Assemble ghost faces + work_on_faces = assemble_own_interior_faces_once + | assemble_own_interior_faces_both + | assemble_ghost_faces_once + | assemble_ghost_faces_both, + /** - * Assemble ghost faces + * Combination of flags to determine if any work is done on the boundary + * faces */ - assemble_ghost_faces = assemble_ghost_faces_both | assemble_ghost_faces_once, + work_on_boundary = assemble_boundary_faces, }; diff --git a/include/deal.II/meshworker/mesh_loop.h b/include/deal.II/meshworker/mesh_loop.h index f52400de9a..86fcab3542 100644 --- a/include/deal.II/meshworker/mesh_loop.h +++ b/include/deal.II/meshworker/mesh_loop.h @@ -129,6 +129,31 @@ namespace MeshWorker const unsigned int queue_length = 2*MultithreadInfo::n_threads(), const unsigned int chunk_size = 8) { + Assert((!cell_worker) + || (flags & work_on_cells), + ExcMessage("If you specify a cell_worker, you need to set assemble_own_cells or assemble_ghost_cells.")); + + Assert((flags & (assemble_own_interior_faces_once|assemble_own_interior_faces_both)) + != (assemble_own_interior_faces_once|assemble_own_interior_faces_both), + ExcMessage("You can only specify assemble_own_interior_faces_once OR assemble_own_interior_faces_both.")); + + Assert((flags & (assemble_ghost_faces_once|assemble_ghost_faces_both)) + != (assemble_ghost_faces_once|assemble_ghost_faces_both), + ExcMessage("You can only specify assemble_ghost_faces_once OR assemble_ghost_faces_both.")); + + Assert(!(flags & cells_first) || + (flags & (assemble_own_cells | assemble_ghost_cells)), + ExcMessage("The option cells_first only makes sense if you assemble on cells.")); + + Assert((!face_worker) + || (flags & work_on_faces), + ExcMessage("If you specify a boundary_worker, assemble_boundary_faces needs to be set.")); + + Assert((!boundary_worker) + || (flags & assemble_boundary_faces), + ExcMessage("If you specify a boundary_worker, assemble_boundary_faces needs to be set.")); + + auto cell_action = [&] (const CellIteratorType &cell, ScratchData &scratch, CopyData ©) { const bool ignore_subdomain = (cell->get_triangulation().locally_owned_subdomain() @@ -143,12 +168,12 @@ namespace MeshWorker if ((!ignore_subdomain) && (current_subdomain_id == numbers::artificial_subdomain_id)) return; - if ( (flags & (assemble_cells_first)) && + if ( (flags & (cells_first)) && ( ((flags & (assemble_own_cells)) && own_cell) || ( (flags & assemble_ghost_cells) && !own_cell) ) ) cell_worker(cell, scratch, copy); - if (flags & assemble_own_faces) + if (flags & (work_on_faces | work_on_boundary)) for (unsigned int face_no=0; face_no < GeometryInfo::faces_per_cell; ++face_no) { typename CellIteratorType::AccessorType::Container::face_iterator face = cell->face(face_no); @@ -156,13 +181,11 @@ namespace MeshWorker { // only integrate boundary faces of own cells if ( (flags & assemble_boundary_faces) && own_cell) - { - boundary_worker(cell, face_no, scratch, copy); - } + boundary_worker(cell, face_no, scratch, copy); } - else if (flags & assemble_own_interior_faces) + else { - // Interior face + // interior face, potentially assemble TriaIterator neighbor = cell->neighbor_or_periodic_neighbor(face_no); types::subdomain_id neighbor_subdomain_id = numbers::artificial_subdomain_id; @@ -180,11 +203,11 @@ namespace MeshWorker continue; // skip if the user doesn't want faces between own cells - if (own_cell && own_neighbor && !(flags & assemble_own_interior_faces)) + if (own_cell && own_neighbor && !(flags & (assemble_own_interior_faces_both | assemble_own_interior_faces_once))) continue; // skip face to ghost - if (own_cell != own_neighbor && !(flags & assemble_ghost_faces)) + if (own_cell != own_neighbor && !(flags & (assemble_ghost_faces_both | assemble_ghost_faces_once))) continue; // Deal with refinement edges from the refined side. Assuming one-irregular @@ -215,6 +238,10 @@ namespace MeshWorker if (flags & assemble_own_interior_faces_both) { + // If own faces are to be assembled from both sides, call the + // faceworker again with swapped arguments. This is because + // we won't be looking at an adaptively refined edge + // coming from the other side. face_worker(neighbor, neighbor_face_no.first, neighbor_face_no.second, cell, face_no, numbers::invalid_unsigned_int, scratch, copy); @@ -238,8 +265,7 @@ namespace MeshWorker && (neighbor < cell)) continue; - // independent of loop_control.faces_to_ghost, - // we only look at faces to ghost on the same level once + // We only look at faces to ghost on the same level once // (only where own_cell=true and own_neighbor=false) if (!own_cell) continue; @@ -260,20 +286,17 @@ namespace MeshWorker face_worker(cell, face_no, numbers::invalid_unsigned_int, neighbor, neighbor_face_no, numbers::invalid_unsigned_int, scratch, copy); - } } } // faces - // Execute this, if faces - // have to be handled first - if ((flags & assemble_own_cells) && !(flags & assemble_cells_first) && + // Execute the cell_worker if faces are handled before cells + if (!(flags & cells_first) && ( ((flags & assemble_own_cells) && own_cell) || ((flags & assemble_ghost_cells) && !own_cell))) cell_worker(cell, scratch, copy); }; - - // Loop over all cells + // Submit to workstream WorkStream::run(begin, end, cell_action, copier, scratch_data, copy_data, diff --git a/tests/meshworker/mesh_loop_02.cc b/tests/meshworker/mesh_loop_02.cc index bdd5feffe7..9846217ed7 100644 --- a/tests/meshworker/mesh_loop_02.cc +++ b/tests/meshworker/mesh_loop_02.cc @@ -76,44 +76,48 @@ void test() deallog << "copier" << std::endl; }; + std::function + empty_cell_worker; + std::function + empty_boundary_worker; + deallog << "CELLS ONLY" << std::endl << std::endl; mesh_loop(cell, endc, cell_worker, copier, scratch, copy, - assemble_own_cells, - boundary_worker, face_worker); + assemble_own_cells); deallog << "BOUNDARY ONLY" << std::endl << std::endl; - mesh_loop(cell, endc, cell_worker, copier, scratch, copy, + mesh_loop(cell, endc, empty_cell_worker, copier, scratch, copy, assemble_boundary_faces, - boundary_worker, face_worker); + boundary_worker); deallog << "CELLS AND BOUNDARY" << std::endl << std::endl; mesh_loop(cell, endc, cell_worker, copier, scratch, copy, assemble_own_cells | assemble_boundary_faces, - boundary_worker, face_worker); + boundary_worker); deallog << "CELLS FIRST AND BOUNDARY" << std::endl << std::endl; mesh_loop(cell, endc, cell_worker, copier, scratch, copy, - assemble_own_cells | assemble_cells_first | assemble_boundary_faces, - boundary_worker, face_worker); + assemble_own_cells | cells_first | assemble_boundary_faces, + boundary_worker); deallog << "ONLY FACES ONCE" << std::endl << std::endl; - mesh_loop(cell, endc, cell_worker, copier, scratch, copy, + mesh_loop(cell, endc, empty_cell_worker, copier, scratch, copy, assemble_own_interior_faces_once, - boundary_worker, face_worker); + empty_boundary_worker, face_worker); deallog << "ONLY FACES BOTH" << std::endl << std::endl; - mesh_loop(cell, endc, cell_worker, copier, scratch, copy, + mesh_loop(cell, endc, empty_cell_worker, copier, scratch, copy, assemble_own_interior_faces_both, - boundary_worker, face_worker); + empty_boundary_worker, face_worker); }