From 198c2b336589fea59676e9086a68ad29f874ebcb Mon Sep 17 00:00:00 2001 From: Qingyuan Shi Date: Sun, 27 Apr 2025 18:56:54 +0800 Subject: [PATCH] Fix: explicitly instantiate BlockSparseMatrix::reinit(vector, BlockDynamicSparsityPattern, MPI_Comm, bool) --- doc/news/changes/minor/20250427Shi | 6 ++ .../trilinos_tpetra_block_sparse_matrix.cc | 8 ++ tests/trilinos_tpetra/block_matrices_01.cc | 91 +++++++++++++++++++ .../trilinos_tpetra/block_matrices_01.output | 4 + 4 files changed, 109 insertions(+) create mode 100644 doc/news/changes/minor/20250427Shi create mode 100644 tests/trilinos_tpetra/block_matrices_01.cc create mode 100644 tests/trilinos_tpetra/block_matrices_01.output diff --git a/doc/news/changes/minor/20250427Shi b/doc/news/changes/minor/20250427Shi new file mode 100644 index 0000000000..d881958a9a --- /dev/null +++ b/doc/news/changes/minor/20250427Shi @@ -0,0 +1,6 @@ +Fix: Add explicit instantiation for TpetraWrappers::BlockSparseMatrix::reinit +(std::vector, BlockDynamicSparsityPattern, MPI_Comm, bool) to fix +missing symbol error when initializing Tpetra block matrices from +distributed sparsity patterns. +
+(Qingyuan Shi, 2025/04/27) diff --git a/source/lac/trilinos_tpetra_block_sparse_matrix.cc b/source/lac/trilinos_tpetra_block_sparse_matrix.cc index a51c197ce9..87c37f1925 100644 --- a/source/lac/trilinos_tpetra_block_sparse_matrix.cc +++ b/source/lac/trilinos_tpetra_block_sparse_matrix.cc @@ -32,6 +32,14 @@ namespace LinearAlgebra BlockSparseMatrix::reinit( const ::dealii::BlockDynamicSparsityPattern &); + template void + BlockSparseMatrix::reinit< + ::dealii::BlockDynamicSparsityPattern>( + const std::vector &, + const ::dealii::BlockDynamicSparsityPattern &, + MPI_Comm, + bool); + template void BlockSparseMatrix::vmult( TpetraWrappers::Vector &, diff --git a/tests/trilinos_tpetra/block_matrices_01.cc b/tests/trilinos_tpetra/block_matrices_01.cc new file mode 100644 index 0000000000..183d094107 --- /dev/null +++ b/tests/trilinos_tpetra/block_matrices_01.cc @@ -0,0 +1,91 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2022 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + +// Verify n_nonzero_elements and reinit() for Tpetra 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; + + + + // create partitioning + std::vector partitioning(2); + partitioning[0] = dealii::IndexSet(2); // 2 dofs + partitioning[1] = dealii::IndexSet(3); // 3 dofs + + partitioning[0].add_range(0, 2); + partitioning[1].add_range(0, 3); + + for (auto &set : partitioning) + set.compress(); + + // create block sparse matrix using another reinit() function + TrilinosWrappers::BlockSparseMatrix tbsm2; + tbsm2.reinit(partitioning, bdsp, MPI_COMM_SELF, false); + + 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_tpetra/block_matrices_01.output b/tests/trilinos_tpetra/block_matrices_01.output new file mode 100644 index 0000000000..bdea2497ce --- /dev/null +++ b/tests/trilinos_tpetra/block_matrices_01.output @@ -0,0 +1,4 @@ + +DEAL::nonzeros BlockSparsityPattern: 9 +DEAL::nonzeros BlockSparseMatrix: 9 +DEAL::nonzeros BlockSparseMatrix: 9 -- 2.39.5