From 2f926fcdc75661c47856e3da93d0cec4c3c2479a Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Fri, 8 May 2020 11:39:30 +0200 Subject: [PATCH] Test case for new feature --- tests/matrix_free/pre_and_post_loops_06.cc | 212 ++++++++++++++++++ ...h_mpi=true.with_p4est=true.mpirun=1.output | 10 + ...h_mpi=true.with_p4est=true.mpirun=3.output | 10 + 3 files changed, 232 insertions(+) create mode 100644 tests/matrix_free/pre_and_post_loops_06.cc create mode 100644 tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=1.output create mode 100644 tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=3.output diff --git a/tests/matrix_free/pre_and_post_loops_06.cc b/tests/matrix_free/pre_and_post_loops_06.cc new file mode 100644 index 0000000000..97c140dfba --- /dev/null +++ b/tests/matrix_free/pre_and_post_loops_06.cc @@ -0,0 +1,212 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 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 the correctness of the matrix-free cell loop with additional +// operation_before_loop and operation_after_loop lambdas when using a +// threaded matrix-free loop + +#include + +#include +#include + +#include +#include + +#include +#include + +#include +#include +#include +#include + +#include + +#include "../tests.h" + +#include "matrix_vector_mf.h" + + + +template +class Matrix +{ +public: + Matrix(const MatrixFree &data_in) + : data(data_in) + {} + + void + vmult(LinearAlgebra::distributed::Vector & dst, + const LinearAlgebra::distributed::Vector &src) const + { + const std::function &, + LinearAlgebra::distributed::Vector &, + const LinearAlgebra::distributed::Vector &, + const std::pair &)> + wrap = helmholtz_operator, + n_q_points_1d>; + dst = 0; + data.cell_loop(wrap, dst, src); + for (auto i : data.get_constrained_dofs()) + dst.local_element(i) = src.local_element(i); + } + + void + vmult_with_update_basic( + LinearAlgebra::distributed::Vector &dst, + LinearAlgebra::distributed::Vector &src, + LinearAlgebra::distributed::Vector &other_vector) const + { + src.add(1.5, other_vector); + other_vector.add(0.7, dst); + dst = 0; + vmult(dst, src); + dst.sadd(0.63, -1.3, other_vector); + } + + void + vmult_with_update_merged( + LinearAlgebra::distributed::Vector &dst, + LinearAlgebra::distributed::Vector &src, + LinearAlgebra::distributed::Vector &other_vector) const + { + const std::function &, + LinearAlgebra::distributed::Vector &, + const LinearAlgebra::distributed::Vector &, + const std::pair &)> + wrap = helmholtz_operator, + n_q_points_1d>; + data.cell_loop( + wrap, + dst, + src, + // operation before cell operation + [&](const unsigned int start_range, const unsigned int end_range) { + for (unsigned int i = start_range; i < end_range; ++i) + { + src.local_element(i) += 1.5 * other_vector.local_element(i); + other_vector.local_element(i) += 0.7 * dst.local_element(i); + dst.local_element(i) = 0; + } + }, + // operation after cell operation + [&](const unsigned int start_range, const unsigned int end_range) { + for (unsigned int i = start_range; i < end_range; ++i) + { + dst.local_element(i) = + 0.63 * dst.local_element(i) - 1.3 * other_vector.local_element(i); + } + }); + } + +private: + const MatrixFree &data; +}; + + + +template +void +test() +{ + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tria); + tria.refine_global(7 - dim); + + FE_Q fe(fe_degree); + DoFHandler dof(tria); + dof.distribute_dofs(fe); + AffineConstraints constraints; + constraints.close(); + + deallog << "Testing " << dof.get_fe().get_name() << std::endl; + + MatrixFree mf_data; + { + const QGauss<1> quad(fe_degree + 1); + typename MatrixFree::AdditionalData data; + data.tasks_parallel_scheme = + MatrixFree::AdditionalData::partition_partition; + mf_data.reinit(dof, constraints, quad, data); + } + + Matrix mf(mf_data); + LinearAlgebra::distributed::Vector vec1, vec2, vec3; + mf_data.initialize_dof_vector(vec1); + mf_data.initialize_dof_vector(vec2); + mf_data.initialize_dof_vector(vec3); + + for (unsigned int i = 0; i < vec1.local_size(); ++i) + { + // Multiply by 0.01 to make float error with roundoff less than the + // numdiff absolute tolerance + double entry = 0.01 * random_value(); + vec1.local_element(i) = entry; + entry = 0.01 * random_value(); + vec2.local_element(i) = entry; + entry = 0.01 * random_value(); + vec3.local_element(i) = entry; + } + + LinearAlgebra::distributed::Vector ref1 = vec1; + LinearAlgebra::distributed::Vector ref2 = vec2; + LinearAlgebra::distributed::Vector ref3 = vec3; + + mf.vmult_with_update_basic(ref3, ref2, ref1); + mf.vmult_with_update_merged(vec3, vec2, vec1); + + vec3 -= ref3; + deallog << "Error in 1x merged operation: " << vec3.linfty_norm() + << std::endl; + + ref3 = 0; + + mf.vmult_with_update_basic(ref1, ref2, ref3); + mf.vmult_with_update_merged(vec1, vec2, vec3); + + mf.vmult_with_update_basic(ref1, ref2, ref3); + mf.vmult_with_update_merged(vec1, vec2, vec3); + + vec3 -= ref3; + deallog << "Error in 2x merged operation: " << vec3.linfty_norm() + << std::endl; +} + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_init(argc, + argv, + testing_max_num_threads()); + + mpi_initlog(); + deallog << std::setprecision(3); + + test<2, 3, double>(); + test<2, 3, float>(); + test<3, 3, double>(); +} diff --git a/tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=1.output b/tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..57a015933a --- /dev/null +++ b/tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=1.output @@ -0,0 +1,10 @@ + +DEAL::Testing FE_Q<2>(3) +DEAL::Error in 1x merged operation: 0.00 +DEAL::Error in 2x merged operation: 0.00 +DEAL::Testing FE_Q<2>(3) +DEAL::Error in 1x merged operation: 0.00 +DEAL::Error in 2x merged operation: 0.00 +DEAL::Testing FE_Q<3>(3) +DEAL::Error in 1x merged operation: 0.00 +DEAL::Error in 2x merged operation: 0.00 diff --git a/tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=3.output b/tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..57a015933a --- /dev/null +++ b/tests/matrix_free/pre_and_post_loops_06.with_mpi=true.with_p4est=true.mpirun=3.output @@ -0,0 +1,10 @@ + +DEAL::Testing FE_Q<2>(3) +DEAL::Error in 1x merged operation: 0.00 +DEAL::Error in 2x merged operation: 0.00 +DEAL::Testing FE_Q<2>(3) +DEAL::Error in 1x merged operation: 0.00 +DEAL::Error in 2x merged operation: 0.00 +DEAL::Testing FE_Q<3>(3) +DEAL::Error in 1x merged operation: 0.00 +DEAL::Error in 2x merged operation: 0.00 -- 2.39.5