From d4fc1d44af1397d08f7c4564d2b6d8ab32e48a3e Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Mon, 14 Mar 2022 20:10:46 +0100 Subject: [PATCH] Add ScalarType as template parameter to MeshWorker::CopyData The `matrices` and `vectors` can be made to store entries of type double, std::complex, float, etc. --- include/deal.II/meshworker/copy_data.h | 52 ++++++--- tests/meshworker/mesh_loop_09a.cc | 143 +++++++++++++++++++++++++ tests/meshworker/mesh_loop_09a.output | 13 +++ tests/meshworker/mesh_loop_09b.cc | 143 +++++++++++++++++++++++++ tests/meshworker/mesh_loop_09b.output | 13 +++ 5 files changed, 347 insertions(+), 17 deletions(-) create mode 100644 tests/meshworker/mesh_loop_09a.cc create mode 100644 tests/meshworker/mesh_loop_09a.output create mode 100644 tests/meshworker/mesh_loop_09b.cc create mode 100644 tests/meshworker/mesh_loop_09b.output diff --git a/include/deal.II/meshworker/copy_data.h b/include/deal.II/meshworker/copy_data.h index 1945437804..69f857b0bf 100644 --- a/include/deal.II/meshworker/copy_data.h +++ b/include/deal.II/meshworker/copy_data.h @@ -44,10 +44,13 @@ namespace MeshWorker * - @tparam n_matrices: Size of the array of matrices * - @tparam n_vectors: size of the array of vectors * - @tparam n_dof_indices: size of the array of local dof indices + * - @tparam ScalarType: The data type stored by the vectors and matrices. + * This must be a real or complex floating point number. */ - template + template struct CopyData { /** @@ -74,8 +77,8 @@ namespace MeshWorker /** * Deep copy constructor. */ - CopyData(const CopyData &other) = - default; + CopyData(const CopyData + &other) = default; /** * Reinitialize everything the same @p size. This is usually the number of @@ -119,12 +122,12 @@ namespace MeshWorker /** * An array of local matrices. */ - std::array, n_matrices> matrices; + std::array, n_matrices> matrices; /** * An array of local vectors. */ - std::array, n_vectors> vectors; + std::array, n_vectors> vectors; /** * An array of local degrees of freedom indices. @@ -138,8 +141,11 @@ namespace MeshWorker // // Template definitions // - template - CopyData::CopyData( + template + CopyData::CopyData( const unsigned int size) { reinit(size); @@ -147,8 +153,11 @@ namespace MeshWorker - template - CopyData::CopyData( + template + CopyData::CopyData( const ndarray & matrix_sizes, const std::array & vector_sizes, const std::array &dof_indices_sizes) @@ -165,9 +174,12 @@ namespace MeshWorker - template + template void - CopyData::reinit( + CopyData::reinit( const unsigned int size) { for (auto &m : matrices) @@ -180,9 +192,12 @@ namespace MeshWorker - template + template void - CopyData::reinit( + CopyData::reinit( const unsigned int index, const unsigned int size) { @@ -191,9 +206,12 @@ namespace MeshWorker - template + template void - CopyData::reinit( + CopyData::reinit( const unsigned int index, const unsigned int size_rows, const unsigned int size_columns) diff --git a/tests/meshworker/mesh_loop_09a.cc b/tests/meshworker/mesh_loop_09a.cc new file mode 100644 index 0000000000..950b788a57 --- /dev/null +++ b/tests/meshworker/mesh_loop_09a.cc @@ -0,0 +1,143 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2019 - 2020 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 CopyData for complex values: vector +// This test is based off meshworker/scratch_data_04.cc. + +#include + +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include "../tests.h" + +using namespace MeshWorker; + +template +void +test() +{ + Triangulation tria; + FE_Q fe(1); + DoFHandler dh(tria); + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + dh.distribute_dofs(fe); + + QGauss quad(3); + QGauss face_quad(3); + + UpdateFlags cell_flags = update_JxW_values; + UpdateFlags face_flags = update_JxW_values; + + using ScratchData = MeshWorker::ScratchData; + using CopyData = MeshWorker::CopyData<0, 1, 0, std::complex>; + + ScratchData scratch(fe, quad, cell_flags, face_quad, face_flags); + CopyData copy(3); // store measures + + using Iterator = decltype(dh.begin_active()); + + auto measures = std::make_tuple(0.0, 0.0, 0.0); + const std::complex i{0, 1}; + + + auto cell_worker = [&i](const Iterator &cell, ScratchData &s, CopyData &c) { + const auto &fev = s.reinit(cell); + const auto &JxW = s.get_JxW_values(); + for (auto w : JxW) + c.vectors[0][0] += w * i; + }; + + auto boundary_worker = [&i](const Iterator & cell, + const unsigned int &f, + ScratchData & s, + CopyData & c) { + const auto &fev = s.reinit(cell, f); + const auto &JxW = s.get_JxW_values(); + for (auto w : JxW) + c.vectors[0][1] += w * i; + }; + + auto face_worker = [&i](const Iterator & cell, + const unsigned int &f, + const unsigned int &sf, + const Iterator & ncell, + const unsigned int &nf, + const unsigned int &nsf, + ScratchData & s, + CopyData & c) { + const auto &fev = s.reinit(cell, f, sf); + const auto &nfev = s.reinit_neighbor(ncell, nf, nsf); + + const auto &JxW = s.get_JxW_values(); + const auto &nJxW = s.get_neighbor_JxW_values(); + for (auto w : JxW) + c.vectors[0][2] += w * i; + }; + + auto copier = [&measures](const CopyData &c) { + std::get<0>(measures) += std::imag(c.vectors[0][0]); + std::get<1>(measures) += std::imag(c.vectors[0][1]); + std::get<2>(measures) += std::imag(c.vectors[0][2]); + }; + + mesh_loop(dh.active_cell_iterators(), + cell_worker, + copier, + scratch, + copy, + assemble_own_cells | assemble_boundary_faces | + assemble_own_interior_faces_once, + boundary_worker, + face_worker); + + deallog << "Testing <" << dim << ',' << spacedim << '>' << std::endl; + + deallog << "Volume: " << std::get<0>(measures) << std::endl; + + deallog << "Boundary surface: " << std::get<1>(measures) << std::endl; + + deallog << "Interior surface: " << std::get<2>(measures) << std::endl; +} + + + +int +main() +{ + initlog(); + + test<2, 2>(); + test<2, 3>(); + test<3, 3>(); +} diff --git a/tests/meshworker/mesh_loop_09a.output b/tests/meshworker/mesh_loop_09a.output new file mode 100644 index 0000000000..b9c1ecd3d3 --- /dev/null +++ b/tests/meshworker/mesh_loop_09a.output @@ -0,0 +1,13 @@ + +DEAL::Testing <2,2> +DEAL::Volume: 1.00000 +DEAL::Boundary surface: 4.00000 +DEAL::Interior surface: 3.00000 +DEAL::Testing <2,3> +DEAL::Volume: 1.00000 +DEAL::Boundary surface: 4.00000 +DEAL::Interior surface: 3.00000 +DEAL::Testing <3,3> +DEAL::Volume: 1.00000 +DEAL::Boundary surface: 6.00000 +DEAL::Interior surface: 3.75000 diff --git a/tests/meshworker/mesh_loop_09b.cc b/tests/meshworker/mesh_loop_09b.cc new file mode 100644 index 0000000000..8ab2d471a1 --- /dev/null +++ b/tests/meshworker/mesh_loop_09b.cc @@ -0,0 +1,143 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2019 - 2020 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 CopyData for complex values: matrix +// This test is based off meshworker/scratch_data_04.cc. + +#include + +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include "../tests.h" + +using namespace MeshWorker; + +template +void +test() +{ + Triangulation tria; + FE_Q fe(1); + DoFHandler dh(tria); + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + dh.distribute_dofs(fe); + + QGauss quad(3); + QGauss face_quad(3); + + UpdateFlags cell_flags = update_JxW_values; + UpdateFlags face_flags = update_JxW_values; + + using ScratchData = MeshWorker::ScratchData; + using CopyData = MeshWorker::CopyData<1, 0, 0, std::complex>; + + ScratchData scratch(fe, quad, cell_flags, face_quad, face_flags); + CopyData copy(3); // store measures + + using Iterator = decltype(dh.begin_active()); + + auto measures = std::make_tuple(0.0, 0.0, 0.0); + const std::complex i{0, 1}; + + + auto cell_worker = [&i](const Iterator &cell, ScratchData &s, CopyData &c) { + const auto &fev = s.reinit(cell); + const auto &JxW = s.get_JxW_values(); + for (auto w : JxW) + c.matrices[0][0][0] += w * i; + }; + + auto boundary_worker = [&i](const Iterator & cell, + const unsigned int &f, + ScratchData & s, + CopyData & c) { + const auto &fev = s.reinit(cell, f); + const auto &JxW = s.get_JxW_values(); + for (auto w : JxW) + c.matrices[0][0][1] += w * i; + }; + + auto face_worker = [&i](const Iterator & cell, + const unsigned int &f, + const unsigned int &sf, + const Iterator & ncell, + const unsigned int &nf, + const unsigned int &nsf, + ScratchData & s, + CopyData & c) { + const auto &fev = s.reinit(cell, f, sf); + const auto &nfev = s.reinit_neighbor(ncell, nf, nsf); + + const auto &JxW = s.get_JxW_values(); + const auto &nJxW = s.get_neighbor_JxW_values(); + for (auto w : JxW) + c.matrices[0][0][2] += w * i; + }; + + auto copier = [&measures](const CopyData &c) { + std::get<0>(measures) += std::imag(c.matrices[0][0][0]); + std::get<1>(measures) += std::imag(c.matrices[0][0][1]); + std::get<2>(measures) += std::imag(c.matrices[0][0][2]); + }; + + mesh_loop(dh.active_cell_iterators(), + cell_worker, + copier, + scratch, + copy, + assemble_own_cells | assemble_boundary_faces | + assemble_own_interior_faces_once, + boundary_worker, + face_worker); + + deallog << "Testing <" << dim << ',' << spacedim << '>' << std::endl; + + deallog << "Volume: " << std::get<0>(measures) << std::endl; + + deallog << "Boundary surface: " << std::get<1>(measures) << std::endl; + + deallog << "Interior surface: " << std::get<2>(measures) << std::endl; +} + + + +int +main() +{ + initlog(); + + test<2, 2>(); + test<2, 3>(); + test<3, 3>(); +} diff --git a/tests/meshworker/mesh_loop_09b.output b/tests/meshworker/mesh_loop_09b.output new file mode 100644 index 0000000000..b9c1ecd3d3 --- /dev/null +++ b/tests/meshworker/mesh_loop_09b.output @@ -0,0 +1,13 @@ + +DEAL::Testing <2,2> +DEAL::Volume: 1.00000 +DEAL::Boundary surface: 4.00000 +DEAL::Interior surface: 3.00000 +DEAL::Testing <2,3> +DEAL::Volume: 1.00000 +DEAL::Boundary surface: 4.00000 +DEAL::Interior surface: 3.00000 +DEAL::Testing <3,3> +DEAL::Volume: 1.00000 +DEAL::Boundary surface: 6.00000 +DEAL::Interior surface: 3.75000 -- 2.39.5