From 9d63315e9e3b8bad58b3b3ecc8ad50d41d4e35c9 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 18 Sep 2019 22:46:02 +0200 Subject: [PATCH] Test case --- tests/matrix_free/parallel_multigrid_mf_02.cc | 277 ++++++++++++++++++ ...4est=true.with_lapack=true.mpirun=3.output | 33 +++ 2 files changed, 310 insertions(+) create mode 100644 tests/matrix_free/parallel_multigrid_mf_02.cc create mode 100644 tests/matrix_free/parallel_multigrid_mf_02.with_mpi=true.with_p4est=true.with_lapack=true.mpirun=3.output diff --git a/tests/matrix_free/parallel_multigrid_mf_02.cc b/tests/matrix_free/parallel_multigrid_mf_02.cc new file mode 100644 index 0000000000..857765176c --- /dev/null +++ b/tests/matrix_free/parallel_multigrid_mf_02.cc @@ -0,0 +1,277 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2014 - 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. +// +// --------------------------------------------------------------------- + + +// same test as parallel_multigrid_mf but using the Laplace operator from +// MatrixFreeOperators::Laplace and handing the vector partitioning of those +// operators to the MG transfer + +#include + +#include + +#include + +#include +#include + +#include +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include +#include + +#include + +#include "../tests.h" + + + +template +class MGCoarseIterative + : public MGCoarseGridBase> +{ +public: + MGCoarseIterative() + {} + + void + initialize(const MatrixType &matrix) + { + coarse_matrix = &matrix; + } + + virtual void + operator()(const unsigned int level, + LinearAlgebra::distributed::Vector & dst, + const LinearAlgebra::distributed::Vector &src) const + { + ReductionControl solver_control(1e4, 1e-50, 1e-10); + SolverCG> solver_coarse( + solver_control); + solver_coarse.solve(*coarse_matrix, dst, src, PreconditionIdentity()); + } + + const MatrixType *coarse_matrix; +}; + + + +template +void +do_test(const DoFHandler &dof) +{ + deallog << "Testing " << dof.get_fe().get_name(); + deallog << std::endl; + deallog << "Number of degrees of freedom: " << dof.n_dofs() << std::endl; + + MappingQ mapping(fe_degree + 1); + MatrixFreeOperators::LaplaceOperator< + dim, + fe_degree, + n_q_points_1d, + 1, + LinearAlgebra::distributed::Vector> + fine_matrix; + + // Dirichlet BC + Functions::ZeroFunction zero_function; + std::map *> dirichlet_boundary; + dirichlet_boundary[0] = &zero_function; + + // fine-level constraints + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dof, locally_relevant_dofs); + AffineConstraints constraints; + constraints.reinit(locally_relevant_dofs); + VectorTools::interpolate_boundary_values(dof, + dirichlet_boundary, + constraints); + constraints.close(); + + std::shared_ptr> fine_level_data( + new MatrixFree()); + + typename MatrixFree::AdditionalData fine_level_additional_data; + fine_level_data->reinit(mapping, + dof, + constraints, + QGauss<1>(n_q_points_1d), + fine_level_additional_data); + + fine_matrix.initialize(fine_level_data); + + LinearAlgebra::distributed::Vector in, sol; + fine_matrix.initialize_dof_vector(in); + fine_matrix.initialize_dof_vector(sol); + + // set constant rhs vector + in = 1.; + + // set up multigrid in analogy to step-37 + MGConstrainedDoFs mg_constrained_dofs; + mg_constrained_dofs.initialize(dof, dirichlet_boundary); + + typedef MatrixFreeOperators::LaplaceOperator< + dim, + fe_degree, + n_q_points_1d, + 1, + LinearAlgebra::distributed::Vector> + LevelMatrixType; + + MGLevelObject mg_matrices; + MGLevelObject>> mg_level_data; + mg_matrices.resize(0, dof.get_triangulation().n_global_levels() - 1); + mg_level_data.resize(0, dof.get_triangulation().n_global_levels() - 1); + for (unsigned int level = 0; + level < dof.get_triangulation().n_global_levels(); + ++level) + { + typename MatrixFree::AdditionalData mg_additional_data; + mg_additional_data.tasks_parallel_scheme = + MatrixFree::AdditionalData::none; + mg_additional_data.tasks_block_size = 3; + mg_additional_data.level_mg_handler = level; + + AffineConstraints level_constraints; + IndexSet relevant_dofs; + DoFTools::extract_locally_relevant_level_dofs(dof, level, relevant_dofs); + level_constraints.reinit(relevant_dofs); + level_constraints.add_lines( + mg_constrained_dofs.get_boundary_indices(level)); + level_constraints.close(); + + mg_level_data[level].reset(new MatrixFree()); + mg_level_data[level]->reinit(mapping, + dof, + level_constraints, + QGauss<1>(n_q_points_1d), + mg_additional_data); + mg_matrices[level].initialize(mg_level_data[level], + mg_constrained_dofs, + level); + mg_matrices[level].compute_diagonal(); + } + + MGTransferMatrixFree mg_transfer; + std::vector> partitioners( + dof.get_triangulation().n_global_levels()); + for (unsigned int level = 0; level < partitioners.size(); ++level) + partitioners[level] = + mg_level_data[level]->get_dof_info().vector_partitioner; + mg_transfer.build(dof, partitioners); + + MGCoarseIterative mg_coarse; + mg_coarse.initialize(mg_matrices[0]); + + typedef PreconditionChebyshev> + SMOOTHER; + MGSmootherPrecondition> + mg_smoother; + + MGLevelObject smoother_data; + smoother_data.resize(0, dof.get_triangulation().n_global_levels() - 1); + for (unsigned int level = 0; + level < dof.get_triangulation().n_global_levels(); + ++level) + { + smoother_data[level].smoothing_range = 15.; + smoother_data[level].degree = 5; + smoother_data[level].eig_cg_n_iterations = 15; + smoother_data[level].preconditioner = + mg_matrices[level].get_matrix_diagonal_inverse(); + } + + // temporarily disable deallog for the setup of the preconditioner that + // involves a CG solver for eigenvalue estimation + mg_smoother.initialize(mg_matrices, smoother_data); + + mg::Matrix> mg_matrix(mg_matrices); + + Multigrid> mg( + dof, mg_matrix, mg_coarse, mg_transfer, mg_smoother, mg_smoother); + PreconditionMG, + MGTransferMatrixFree> + preconditioner(dof, mg, mg_transfer); + + { + // avoid output from inner (coarse-level) solver + deallog.depth_file(2); + ReductionControl control(30, 1e-20, 1e-7); + SolverCG> solver(control); + solver.solve(fine_matrix, sol, in, preconditioner); + } +} + + + +template +void +test() +{ + for (unsigned int i = 5; i < 7; ++i) + { + parallel::distributed::Triangulation tria( + MPI_COMM_WORLD, + Triangulation::limit_level_difference_at_vertices, + parallel::distributed::Triangulation< + dim>::construct_multigrid_hierarchy); + GridGenerator::hyper_cube(tria); + tria.refine_global(i - dim); + + FE_Q fe(fe_degree); + DoFHandler dof(tria); + dof.distribute_dofs(fe); + dof.distribute_mg_dofs(); + + do_test(dof); + } +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_init(argc, argv, 1); + + mpi_initlog(); + + { + test<2, 1, double>(); + test<2, 2, float>(); + + test<3, 1, double>(); + test<3, 2, float>(); + } +} diff --git a/tests/matrix_free/parallel_multigrid_mf_02.with_mpi=true.with_p4est=true.with_lapack=true.mpirun=3.output b/tests/matrix_free/parallel_multigrid_mf_02.with_mpi=true.with_p4est=true.with_lapack=true.mpirun=3.output new file mode 100644 index 0000000000..604fbf3e47 --- /dev/null +++ b/tests/matrix_free/parallel_multigrid_mf_02.with_mpi=true.with_p4est=true.with_lapack=true.mpirun=3.output @@ -0,0 +1,33 @@ + +DEAL::Testing FE_Q<2>(1) +DEAL::Number of degrees of freedom: 81 +DEAL:cg::Starting value 9.00000 +DEAL:cg::Convergence step 4 value 2.07070e-08 +DEAL::Testing FE_Q<2>(1) +DEAL::Number of degrees of freedom: 289 +DEAL:cg::Starting value 17.0000 +DEAL:cg::Convergence step 4 value 2.77789e-08 +DEAL::Testing FE_Q<2>(2) +DEAL::Number of degrees of freedom: 289 +DEAL:cg::Starting value 17.0000 +DEAL:cg::Convergence step 4 value 3.93059e-08 +DEAL::Testing FE_Q<2>(2) +DEAL::Number of degrees of freedom: 1089 +DEAL:cg::Starting value 33.0000 +DEAL:cg::Convergence step 4 value 8.11904e-08 +DEAL::Testing FE_Q<3>(1) +DEAL::Number of degrees of freedom: 125 +DEAL:cg::Starting value 11.1803 +DEAL:cg::Convergence step 4 value 3.45976e-07 +DEAL::Testing FE_Q<3>(1) +DEAL::Number of degrees of freedom: 729 +DEAL:cg::Starting value 27.0000 +DEAL:cg::Convergence step 4 value 1.07425e-07 +DEAL::Testing FE_Q<3>(2) +DEAL::Number of degrees of freedom: 729 +DEAL:cg::Starting value 27.0000 +DEAL:cg::Convergence step 4 value 2.46943e-06 +DEAL::Testing FE_Q<3>(2) +DEAL::Number of degrees of freedom: 4913 +DEAL:cg::Starting value 70.0928 +DEAL:cg::Convergence step 4 value 4.57526e-06 -- 2.39.5