From 6081ccfebeedef8ee769fb45ef46d3cff966b3fd Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Tue, 9 Apr 2019 23:04:47 +0200 Subject: [PATCH] Test ScratchData on simple volume computation --- tests/meshworker/scratch_data_01.cc | 141 ++++++++++++++++++++++++ tests/meshworker/scratch_data_01.output | 13 +++ 2 files changed, 154 insertions(+) create mode 100644 tests/meshworker/scratch_data_01.cc create mode 100644 tests/meshworker/scratch_data_01.output diff --git a/tests/meshworker/scratch_data_01.cc b/tests/meshworker/scratch_data_01.cc new file mode 100644 index 0000000000..52dec18793 --- /dev/null +++ b/tests/meshworker/scratch_data_01.cc @@ -0,0 +1,141 @@ +/* --------------------------------------------------------------------- + * + * 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 ScratchData on simple volume computation + +#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>; + + ScratchData scratch(fe, quad, cell_flags, face_quad, face_flags); + CopyData copy(3); // store measures + + auto cell = dh.begin_active(); + auto endc = dh.end(); + + typedef decltype(cell) Iterator; + + auto measures = std::make_tuple(0.0, 0.0, 0.0); + + + auto cell_worker = [](const Iterator &cell, ScratchData &s, CopyData &c) { + const auto &fev = s.reinit(cell); + const auto &JxW = s.get_JxW_values(); + c = 0; + for (auto w : JxW) + c.vectors[0][0] += w; + }; + + auto boundary_worker = [](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; + }; + + auto face_worker = [](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; + }; + + auto copier = [&measures](const CopyData &c) { + std::get<0>(measures) += c.vectors[0][0]; + std::get<1>(measures) += c.vectors[0][1]; + std::get<2>(measures) += c.vectors[0][2]; + }; + + mesh_loop(cell, + endc, + 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/scratch_data_01.output b/tests/meshworker/scratch_data_01.output new file mode 100644 index 0000000000..b9c1ecd3d3 --- /dev/null +++ b/tests/meshworker/scratch_data_01.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