From 9381371badad66b427af4ba6607ba3818d11d29f Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Fri, 22 Dec 2017 09:14:09 -0500 Subject: [PATCH] implement missing TrilinosWrappers::SparsityPattern::copy_from() also test the different copy_from implementations. --- doc/news/changes/minor/20171222Heister | 3 + source/lac/trilinos_sparsity_pattern.cc | 14 +++++ tests/trilinos/sparsity_pattern_05.cc | 73 +++++++++++++++++++++++ tests/trilinos/sparsity_pattern_05.output | 10 ++++ 4 files changed, 100 insertions(+) create mode 100644 doc/news/changes/minor/20171222Heister create mode 100644 tests/trilinos/sparsity_pattern_05.cc create mode 100644 tests/trilinos/sparsity_pattern_05.output diff --git a/doc/news/changes/minor/20171222Heister b/doc/news/changes/minor/20171222Heister new file mode 100644 index 0000000000..88f1a2d206 --- /dev/null +++ b/doc/news/changes/minor/20171222Heister @@ -0,0 +1,3 @@ +Fixed: implement missing instantiation of TrilinosWrappers::SparsityPattern::copy_from(). +
+(Timo Heister, 2017/12/22) diff --git a/source/lac/trilinos_sparsity_pattern.cc b/source/lac/trilinos_sparsity_pattern.cc index dbf1802db7..050876e1e7 100644 --- a/source/lac/trilinos_sparsity_pattern.cc +++ b/source/lac/trilinos_sparsity_pattern.cc @@ -659,6 +659,20 @@ namespace TrilinosWrappers + void + SparsityPattern::copy_from (const SparsityPattern &sp) + { + column_space_map.reset (new Epetra_Map (*sp.column_space_map)); + graph.reset (new Epetra_FECrsGraph(*sp.graph)); + + if (sp.nonlocal_graph.get()!=nullptr) + nonlocal_graph.reset(new Epetra_CrsGraph(*sp.nonlocal_graph)); + else + nonlocal_graph.reset(); + } + + + template void SparsityPattern::copy_from (const SparsityPatternType &sp) diff --git a/tests/trilinos/sparsity_pattern_05.cc b/tests/trilinos/sparsity_pattern_05.cc new file mode 100644 index 0000000000..4976ed7294 --- /dev/null +++ b/tests/trilinos/sparsity_pattern_05.cc @@ -0,0 +1,73 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Test copy_from(T) + +#include "../tests.h" +#include +#include +#include +#include + +int main (int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads()); + + IndexSet partitioning(3); + + partitioning.add_range(0, 3); + + // Add element (2,1) to the matrix + TrilinosWrappers::SparsityPattern A(partitioning); + A.add(2,1); + A.compress(); + + // Check copy_from(TrilinosWrappers::SparsityPattern): + TrilinosWrappers::SparsityPattern B; + B.copy_from(A); + for (unsigned int i=0; i<3; ++i) + for (unsigned int j=0; j<3; ++j) + { + if ((i == 2) && (j == 1)) + { + AssertThrow(B.exists(i,j) == true, ExcInternalError()); + } + else + { + AssertThrow(B.exists(i,j) == false, ExcInternalError()); + } + } + deallog << "OK" << std::endl; + + // copy_from(DynamicSparsityPattern) + DynamicSparsityPattern dsp(4,4); + dsp.add(2,3); + B.copy_from(dsp); + B.print(deallog.get_file_stream()); + + deallog << "OK" << std::endl; + + SparsityPattern sp(4,4); + sp.add(1,2); + sp.add(3,3); + sp.compress(); + B.copy_from(sp); + B.print(deallog.get_file_stream()); + + deallog << "OK" << std::endl; +} diff --git a/tests/trilinos/sparsity_pattern_05.output b/tests/trilinos/sparsity_pattern_05.output new file mode 100644 index 0000000000..9d9e86c73f --- /dev/null +++ b/tests/trilinos/sparsity_pattern_05.output @@ -0,0 +1,10 @@ + +DEAL::OK +(2,3) +DEAL::OK +(0,0) +(1,1) +(1,2) +(2,2) +(3,3) +DEAL::OK -- 2.39.5