From 7ef352a9a076dee0e88f2677b8d0f47db3cfa60f Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Fri, 4 Mar 2022 20:10:57 -0700 Subject: [PATCH] Tests for n_nonzero_elements of BlockSparseMatrix from third party libraries. --- ...block_matrices.cc => block_matrices_01.cc} | 0 ...trices.output => block_matrices_01.output} | 0 tests/petsc/block_matrices_01.cc | 84 +++++++++++++++++++ tests/petsc/block_matrices_01.mpirun=1.output | 3 + tests/trilinos/block_matrices_01.cc | 72 ++++++++++++++++ .../block_matrices_01.mpirun=1.output | 3 + 6 files changed, 162 insertions(+) rename tests/lac/{block_matrices.cc => block_matrices_01.cc} (100%) rename tests/lac/{block_matrices.output => block_matrices_01.output} (100%) create mode 100644 tests/petsc/block_matrices_01.cc create mode 100644 tests/petsc/block_matrices_01.mpirun=1.output create mode 100644 tests/trilinos/block_matrices_01.cc create mode 100644 tests/trilinos/block_matrices_01.mpirun=1.output diff --git a/tests/lac/block_matrices.cc b/tests/lac/block_matrices_01.cc similarity index 100% rename from tests/lac/block_matrices.cc rename to tests/lac/block_matrices_01.cc diff --git a/tests/lac/block_matrices.output b/tests/lac/block_matrices_01.output similarity index 100% rename from tests/lac/block_matrices.output rename to tests/lac/block_matrices_01.output diff --git a/tests/petsc/block_matrices_01.cc b/tests/petsc/block_matrices_01.cc new file mode 100644 index 0000000000..51dc0875dd --- /dev/null +++ b/tests/petsc/block_matrices_01.cc @@ -0,0 +1,84 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 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. +// +// --------------------------------------------------------------------- + + +// Verify n_nonzero_elements for PETSc BlockSparseMatrix. + +#include + +#include +#include + +#include "../tests.h" + + + +void +test() +{ + // create block sparsity pattern + BlockDynamicSparsityPattern bdsp(2, 2); + bdsp.block(0, 0).reinit(2, 2); + bdsp.block(0, 1).reinit(2, 3); + bdsp.block(1, 0).reinit(3, 2); + bdsp.block(1, 1).reinit(3, 3); + bdsp.collect_sizes(); + + // add triangular sparsity pattern to each block + // block (0, 0) + for (unsigned int row = 0; row < 2; ++row) + for (unsigned int col = 0; col <= row; ++col) + bdsp.add(row, col); + // block (1, 1) + for (unsigned int row = 2; row < 5; ++row) + for (unsigned int col = 2; col <= row; ++col) + bdsp.add(row, col); + + bdsp.compress(); + deallog << "nonzeros BlockSparsityPattern: " << bdsp.n_nonzero_elements() + << std::endl; + + // create block sparse matrix + std::vector rows(2); + rows[0].set_size(2); + rows[0].add_range(0, 2); + rows[1].set_size(3); + rows[1].add_range(0, 3); + + std::vector cols(2); + cols[0].set_size(2); + cols[0].add_range(0, 2); + cols[1].set_size(3); + cols[1].add_range(0, 3); + + PETScWrappers::MPI::BlockSparseMatrix pbsm; + pbsm.reinit(rows, cols, bdsp, MPI_COMM_WORLD); + deallog << "nonzeros BlockSparseMatrix: " << pbsm.n_nonzero_elements() + << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + + test(); + + return 0; +} diff --git a/tests/petsc/block_matrices_01.mpirun=1.output b/tests/petsc/block_matrices_01.mpirun=1.output new file mode 100644 index 0000000000..a9c5e15408 --- /dev/null +++ b/tests/petsc/block_matrices_01.mpirun=1.output @@ -0,0 +1,3 @@ + +DEAL::nonzeros BlockSparsityPattern: 9 +DEAL::nonzeros BlockSparseMatrix: 9 diff --git a/tests/trilinos/block_matrices_01.cc b/tests/trilinos/block_matrices_01.cc new file mode 100644 index 0000000000..118d39bc94 --- /dev/null +++ b/tests/trilinos/block_matrices_01.cc @@ -0,0 +1,72 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 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. +// +// --------------------------------------------------------------------- + + +// Verify n_nonzero_elements for Trilinos BlockSparseMatrix. + +#include + +#include +#include + +#include "../tests.h" + + + +void +test() +{ + // create block sparsity pattern + BlockDynamicSparsityPattern bdsp(2, 2); + bdsp.block(0, 0).reinit(2, 2); + bdsp.block(0, 1).reinit(2, 3); + bdsp.block(1, 0).reinit(3, 2); + bdsp.block(1, 1).reinit(3, 3); + bdsp.collect_sizes(); + + // add triangular sparsity pattern to each block + // block (0, 0) + for (unsigned int row = 0; row < 2; ++row) + for (unsigned int col = 0; col <= row; ++col) + bdsp.add(row, col); + // block (1, 1) + for (unsigned int row = 2; row < 5; ++row) + for (unsigned int col = 2; col <= row; ++col) + bdsp.add(row, col); + + bdsp.compress(); + deallog << "nonzeros BlockSparsityPattern: " << bdsp.n_nonzero_elements() + << std::endl; + + // create block sparse matrix + TrilinosWrappers::BlockSparseMatrix tbsm; + tbsm.reinit(bdsp); + deallog << "nonzeros BlockSparseMatrix: " << tbsm.n_nonzero_elements() + << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + + test(); + + return 0; +} diff --git a/tests/trilinos/block_matrices_01.mpirun=1.output b/tests/trilinos/block_matrices_01.mpirun=1.output new file mode 100644 index 0000000000..a9c5e15408 --- /dev/null +++ b/tests/trilinos/block_matrices_01.mpirun=1.output @@ -0,0 +1,3 @@ + +DEAL::nonzeros BlockSparsityPattern: 9 +DEAL::nonzeros BlockSparseMatrix: 9 -- 2.39.5