From 35aed646e7878f1bb5fd0f1a9948e5ec5b2f3b2d Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Sat, 23 Mar 2019 03:12:19 +0100 Subject: [PATCH] New variant of mesh_loop. --- doc/news/changes/minor/20190323LucaHeltai | 3 + include/deal.II/meshworker/mesh_loop.h | 110 ++++++++++++++++ tests/meshworker/mesh_loop_class_01.cc | 144 +++++++++++++++++++++ tests/meshworker/mesh_loop_class_01.output | 55 ++++++++ 4 files changed, 312 insertions(+) create mode 100644 doc/news/changes/minor/20190323LucaHeltai create mode 100644 tests/meshworker/mesh_loop_class_01.cc create mode 100644 tests/meshworker/mesh_loop_class_01.output diff --git a/doc/news/changes/minor/20190323LucaHeltai b/doc/news/changes/minor/20190323LucaHeltai new file mode 100644 index 0000000000..265ceaaf7c --- /dev/null +++ b/doc/news/changes/minor/20190323LucaHeltai @@ -0,0 +1,3 @@ +New: Added a variant of mesh_loop that takes a class and its member functions as workers and copiers. +
+(Luca Heltai, 2019/03/23) diff --git a/include/deal.II/meshworker/mesh_loop.h b/include/deal.II/meshworker/mesh_loop.h index f2a53f0502..8589fbba06 100644 --- a/include/deal.II/meshworker/mesh_loop.h +++ b/include/deal.II/meshworker/mesh_loop.h @@ -392,6 +392,116 @@ namespace MeshWorker queue_length, chunk_size); } + + /** + * This is a variant of the mesh_loop function, that can be used for worker + * and copier functions that are member functions of a class. + * + * The argument passed as @p end must be convertible to the same type as @p + * begin, but doesn't have to be of the same type itself. This allows to + * write code like mesh_loop(dof_handler.begin_active(), + * dof_handler.end(), ... where the first is of type + * DoFHandler::active_cell_iterator whereas the second is of type + * DoFHandler::raw_cell_iterator. + * + * The @p queue_length argument indicates the number of items that can be + * live at any given time. Each item consists of @p chunk_size elements of + * the input stream that will be worked on by the worker and copier + * functions one after the other on the same thread. + * + * @note If your data objects are large, or their constructors are + * expensive, it is helpful to keep in mind that queue_length + * copies of the ScratchData object and + * queue_length*chunk_size copies of the CopyData object + * are generated. + * + * @ingroup MeshWorker + * @author Luca Heltai, 2019 + */ + template + void + mesh_loop(const CellIteratorType & begin, + const typename identity::type &end, + MainClass & main_class, + void (MainClass::*cell_worker)(const CellIteratorType &, + ScratchData &, + CopyData &), + void (MainClass::*copier)(const CopyData &), + const ScratchData & sample_scratch_data, + const CopyData & sample_copy_data, + const AssembleFlags flags = assemble_own_cells, + void (MainClass::*boundary_worker)(const CellIteratorType &, + const unsigned int, + ScratchData &, + CopyData &) = nullptr, + void (MainClass::*face_worker)(const CellIteratorType &, + const unsigned int, + const unsigned int, + const CellIteratorType &, + const unsigned int, + const unsigned int, + ScratchData &, + CopyData &) = nullptr, + const unsigned int queue_length = 2 * MultithreadInfo::n_threads(), + const unsigned int chunk_size = 8) + { + std::function + f_cell_worker; + + std::function + f_boundary_worker; + + std::function + f_face_worker; + + if (cell_worker != nullptr) + f_cell_worker = std::bind(cell_worker, + std::ref(main_class), + std::placeholders::_1, + std::placeholders::_2, + std::placeholders::_3); + + if (boundary_worker != nullptr) + f_boundary_worker = std::bind(boundary_worker, + std::ref(main_class), + std::placeholders::_1, + std::placeholders::_2, + std::placeholders::_3, + std::placeholders::_4); + + if (face_worker != nullptr) + f_face_worker = std::bind(face_worker, + std::ref(main_class), + std::placeholders::_1, + std::placeholders::_2, + std::placeholders::_3, + std::placeholders::_4, + std::placeholders::_5, + std::placeholders::_6, + std::placeholders::_7, + std::placeholders::_8); + + mesh_loop(begin, + end, + f_cell_worker, + std::bind(copier, main_class, std::placeholders::_1), + sample_scratch_data, + sample_copy_data, + flags, + f_boundary_worker, + f_face_worker); + } } // namespace MeshWorker DEAL_II_NAMESPACE_CLOSE diff --git a/tests/meshworker/mesh_loop_class_01.cc b/tests/meshworker/mesh_loop_class_01.cc new file mode 100644 index 0000000000..bda23d029d --- /dev/null +++ b/tests/meshworker/mesh_loop_class_01.cc @@ -0,0 +1,144 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +// Test assembly using a class and mesh_loop + +#include + +#include +#include +#include + +#include +#include + +#include + +#include "../tests.h" + +using namespace MeshWorker; + +struct ScratchData +{}; +struct CopyData +{}; + +template +class TestClass +{ +public: + using CellIteratorType = + typename Triangulation::active_cell_iterator; + + void + cell_worker(const CellIteratorType &cell, ScratchData &, CopyData &) + { + deallog << "Workgin on cell " << cell << std::endl; + } + + void + copier(const CopyData &) + {} + + void + boundary_worker(const CellIteratorType &cell, + const unsigned int f, + ScratchData &, + CopyData &) + { + deallog << "Boundary worker on : " << cell << ", Face : " << f << std::endl; + } + + void + face_worker(const CellIteratorType &cell, + const unsigned int f, + const unsigned int sf, + const CellIteratorType &ncell, + const unsigned int nf, + const unsigned int nsf, + ScratchData &, + CopyData &) + { + deallog << "Face worker on : " << cell << ", Neighbor cell : " << ncell + << ", Face : " << f << ", Neighbor Face : " << nf + << ", Subface: " << sf << ", Neighbor Subface: " << nsf + << std::endl; + } +}; + +template +void +test() +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + + TestClass test_class; + ScratchData scratch; + CopyData copy; + + auto cell = tria.begin_active(); + auto endc = tria.end(); + + deallog << "CELLS ONLY" << std::endl; + mesh_loop(cell, + endc, + test_class, + &TestClass::cell_worker, + &TestClass::copier, + scratch, + copy, + assemble_own_cells); + + deallog << "CELLS+BOUNDARY" << std::endl; + + mesh_loop(cell, + endc, + test_class, + &TestClass::cell_worker, + &TestClass::copier, + scratch, + copy, + assemble_own_cells | assemble_boundary_faces, + &TestClass::boundary_worker); + + + deallog << "CELLS+BOUNDARY+FACES" << std::endl; + + mesh_loop(cell, + endc, + test_class, + &TestClass::cell_worker, + &TestClass::copier, + scratch, + copy, + assemble_own_cells | assemble_boundary_faces | + assemble_own_interior_faces_once, + &TestClass::boundary_worker, + &TestClass::face_worker); +} + + +int +main() +{ + initlog(); + MultithreadInfo::set_thread_limit(1); // to make output deterministic + + test<2, 2>(); +} diff --git a/tests/meshworker/mesh_loop_class_01.output b/tests/meshworker/mesh_loop_class_01.output new file mode 100644 index 0000000000..63a972588d --- /dev/null +++ b/tests/meshworker/mesh_loop_class_01.output @@ -0,0 +1,55 @@ + +DEAL::CELLS ONLY +DEAL::Workgin on cell 1.1 +DEAL::Workgin on cell 1.2 +DEAL::Workgin on cell 1.3 +DEAL::Workgin on cell 2.0 +DEAL::Workgin on cell 2.1 +DEAL::Workgin on cell 2.2 +DEAL::Workgin on cell 2.3 +DEAL::CELLS+BOUNDARY +DEAL::Workgin on cell 1.1 +DEAL::Boundary worker on : 1.1, Face : 1 +DEAL::Boundary worker on : 1.1, Face : 2 +DEAL::Workgin on cell 1.2 +DEAL::Boundary worker on : 1.2, Face : 0 +DEAL::Boundary worker on : 1.2, Face : 3 +DEAL::Workgin on cell 1.3 +DEAL::Boundary worker on : 1.3, Face : 1 +DEAL::Boundary worker on : 1.3, Face : 3 +DEAL::Workgin on cell 2.0 +DEAL::Boundary worker on : 2.0, Face : 0 +DEAL::Boundary worker on : 2.0, Face : 2 +DEAL::Workgin on cell 2.1 +DEAL::Boundary worker on : 2.1, Face : 2 +DEAL::Workgin on cell 2.2 +DEAL::Boundary worker on : 2.2, Face : 0 +DEAL::Workgin on cell 2.3 +DEAL::CELLS+BOUNDARY+FACES +DEAL::Workgin on cell 1.1 +DEAL::Boundary worker on : 1.1, Face : 1 +DEAL::Boundary worker on : 1.1, Face : 2 +DEAL::Face worker on : 1.1, Neighbor cell : 1.3, Face : 3, Neighbor Face : 2, Subface: 4294967295, Neighbor Subface: 4294967295 +DEAL::Workgin on cell 1.2 +DEAL::Boundary worker on : 1.2, Face : 0 +DEAL::Face worker on : 1.2, Neighbor cell : 1.3, Face : 1, Neighbor Face : 0, Subface: 4294967295, Neighbor Subface: 4294967295 +DEAL::Boundary worker on : 1.2, Face : 3 +DEAL::Workgin on cell 1.3 +DEAL::Boundary worker on : 1.3, Face : 1 +DEAL::Boundary worker on : 1.3, Face : 3 +DEAL::Workgin on cell 2.0 +DEAL::Boundary worker on : 2.0, Face : 0 +DEAL::Face worker on : 2.0, Neighbor cell : 2.1, Face : 1, Neighbor Face : 0, Subface: 4294967295, Neighbor Subface: 4294967295 +DEAL::Boundary worker on : 2.0, Face : 2 +DEAL::Face worker on : 2.0, Neighbor cell : 2.2, Face : 3, Neighbor Face : 2, Subface: 4294967295, Neighbor Subface: 4294967295 +DEAL::Workgin on cell 2.1 +DEAL::Face worker on : 2.1, Neighbor cell : 1.1, Face : 1, Neighbor Face : 0, Subface: 4294967295, Neighbor Subface: 0 +DEAL::Boundary worker on : 2.1, Face : 2 +DEAL::Face worker on : 2.1, Neighbor cell : 2.3, Face : 3, Neighbor Face : 2, Subface: 4294967295, Neighbor Subface: 4294967295 +DEAL::Workgin on cell 2.2 +DEAL::Boundary worker on : 2.2, Face : 0 +DEAL::Face worker on : 2.2, Neighbor cell : 2.3, Face : 1, Neighbor Face : 0, Subface: 4294967295, Neighbor Subface: 4294967295 +DEAL::Face worker on : 2.2, Neighbor cell : 1.2, Face : 3, Neighbor Face : 2, Subface: 4294967295, Neighbor Subface: 0 +DEAL::Workgin on cell 2.3 +DEAL::Face worker on : 2.3, Neighbor cell : 1.1, Face : 1, Neighbor Face : 0, Subface: 4294967295, Neighbor Subface: 1 +DEAL::Face worker on : 2.3, Neighbor cell : 1.2, Face : 3, Neighbor Face : 2, Subface: 4294967295, Neighbor Subface: 1 -- 2.39.5