From e867bae7031e13d3999f702402a941d8cb00de94 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Sat, 30 Jul 2016 20:24:12 +0200 Subject: [PATCH] Fix print functions of Trilinos sparsity pattern --- doc/news/changes.h | 21 +- source/lac/trilinos_sparsity_pattern.cc | 17 +- .../trilinos/trilinos_sparsity_pattern_05.cc | 71 ++++ ...pattern_05.with_p4est=true.mpirun=3.output | 343 ++++++++++++++++++ 4 files changed, 437 insertions(+), 15 deletions(-) create mode 100644 tests/trilinos/trilinos_sparsity_pattern_05.cc create mode 100644 tests/trilinos/trilinos_sparsity_pattern_05.with_p4est=true.mpirun=3.output diff --git a/doc/news/changes.h b/doc/news/changes.h index 761121a9cd..d917083a7d 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -395,14 +395,11 @@ inconvenience this causes.

Specific improvements

    -
  1. New: There is now a new DoFTools::make_flux_sparsity_pattern() - which takes a constraint matrix and flux and internal dof masks, in - parallel. This is useful in the case where some components of a - finite element are continuous and some discontinuous, allowing - constraints to be imposed on the continuous part while also building - building the flux terms needed for the discontinuous part. +
  2. Fixed: The TrilinosWrappers::SparsityPattern::print() and + TrilinosWrappers::SparsityPattern::print_gnuplot() methods did not produce + correct output on distributed computations. This is now fixed.
    - (Sam Cox, 2016/07/25) + (Martin Kronbichler, 2016/07/30)
  3. Fixed: Level indices for geometric multigrid queried through @@ -412,6 +409,16 @@ inconvenience this causes. (Martin Kronbichler, 2016/07/27)
  4. +
  5. New: There is now a new DoFTools::make_flux_sparsity_pattern() + which takes a constraint matrix and flux and internal dof masks, in + parallel. This is useful in the case where some components of a + finite element are continuous and some discontinuous, allowing + constraints to be imposed on the continuous part while also building + building the flux terms needed for the discontinuous part. +
    + (Sam Cox, 2016/07/25) +
  6. +
  7. Improved: Allow for including dofs for individual components on boundary in DoFTools::make_vertex_patches().
    diff --git a/source/lac/trilinos_sparsity_pattern.cc b/source/lac/trilinos_sparsity_pattern.cc index 040c5382bd..bdc7f9c735 100644 --- a/source/lac/trilinos_sparsity_pattern.cc +++ b/source/lac/trilinos_sparsity_pattern.cc @@ -67,9 +67,9 @@ namespace TrilinosWrappers return graph.NumGlobalEntries(); } - int global_row_index(const Epetra_CrsGraph &graph, int i) + int local_to_global_index(const Epetra_BlockMap &map, int i) { - return graph.GRID(i); + return map.GID(i); } #else long long int n_global_elements (const Epetra_BlockMap &map) @@ -102,9 +102,9 @@ namespace TrilinosWrappers return graph.NumGlobalEntries64(); } - long long int global_row_index(const Epetra_CrsGraph &graph, int i) + long long int local_to_global_index(const Epetra_BlockMap &map, int i) { - return graph.GRID64(i); + return map.GID64(i); } #endif } @@ -1097,7 +1097,8 @@ namespace TrilinosWrappers { graph->ExtractMyRowView (i, num_entries, indices); for (int j=0; jRowMap(), i) + << "," << local_to_global_index(graph->ColMap(), indices[j]) << ") " << std::endl; } } @@ -1117,14 +1118,14 @@ namespace TrilinosWrappers int num_entries; graph->ExtractMyRowView (row, num_entries, indices); - for (unsigned int j=0; j<(unsigned int)num_entries; ++j) + for (int j=0; j(j))] - << " " << -static_cast(row) << std::endl; + out << static_cast(local_to_global_index(graph->ColMap(), indices[j])) + << " " << -static_cast(local_to_global_index(graph->RowMap(), row)) << std::endl; } AssertThrow (out, ExcIO()); diff --git a/tests/trilinos/trilinos_sparsity_pattern_05.cc b/tests/trilinos/trilinos_sparsity_pattern_05.cc new file mode 100644 index 0000000000..7152036720 --- /dev/null +++ b/tests/trilinos/trilinos_sparsity_pattern_05.cc @@ -0,0 +1,71 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 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 print functions of Trilinos sparsity pattern + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +void test () +{ + const int dim = 2; + //setup system + dealii::parallel::distributed::Triangulation triangulation(MPI_COMM_WORLD); + + GridGenerator::hyper_cube (triangulation); + + triangulation.refine_global(2); + + const FE_Q fe_system(1); + DoFHandler dh (triangulation); + dh.distribute_dofs (fe_system); + + IndexSet relevant_partitioning (dh.n_dofs()); + DoFTools::extract_locally_relevant_dofs (dh, relevant_partitioning); + + //generate empty constraints + ConstraintMatrix constraints; + + //generate sparsity pattern + TrilinosWrappers::SparsityPattern sp (dh.locally_owned_dofs(), + dh.locally_owned_dofs(), + relevant_partitioning, + MPI_COMM_WORLD); + + DoFTools::make_sparsity_pattern (dh, sp, + constraints,true, + Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)); + sp.compress(); + + //output + sp.print_gnuplot(deallog.get_file_stream()); + sp.print(deallog.get_file_stream()); +} + +int main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads()); + MPILogInitAll log; + + test(); +} diff --git a/tests/trilinos/trilinos_sparsity_pattern_05.with_p4est=true.mpirun=3.output b/tests/trilinos/trilinos_sparsity_pattern_05.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..53d639715c --- /dev/null +++ b/tests/trilinos/trilinos_sparsity_pattern_05.with_p4est=true.mpirun=3.output @@ -0,0 +1,343 @@ + +0 0 +1 0 +2 0 +3 0 +0 -1 +1 -1 +2 -1 +3 -1 +4 -1 +5 -1 +0 -2 +1 -2 +2 -2 +3 -2 +6 -2 +7 -2 +0 -3 +1 -3 +2 -3 +3 -3 +4 -3 +5 -3 +6 -3 +7 -3 +8 -3 +1 -4 +3 -4 +4 -4 +5 -4 +9 -4 +10 -4 +1 -5 +3 -5 +4 -5 +5 -5 +7 -5 +8 -5 +9 -5 +10 -5 +13 -5 +2 -6 +3 -6 +6 -6 +7 -6 +15 -6 +16 -6 +2 -7 +3 -7 +5 -7 +6 -7 +7 -7 +8 -7 +15 -7 +16 -7 +17 -7 +3 -8 +5 -8 +7 -8 +8 -8 +10 -8 +13 -8 +16 -8 +17 -8 +21 -8 +(0,0) +(0,1) +(0,2) +(0,3) +(1,0) +(1,1) +(1,2) +(1,3) +(1,4) +(1,5) +(2,0) +(2,1) +(2,2) +(2,3) +(2,6) +(2,7) +(3,0) +(3,1) +(3,2) +(3,3) +(3,4) +(3,5) +(3,6) +(3,7) +(3,8) +(4,1) +(4,3) +(4,4) +(4,5) +(4,9) +(4,10) +(5,1) +(5,3) +(5,4) +(5,5) +(5,7) +(5,8) +(5,9) +(5,10) +(5,13) +(6,2) +(6,3) +(6,6) +(6,7) +(6,15) +(6,16) +(7,2) +(7,3) +(7,5) +(7,6) +(7,7) +(7,8) +(7,15) +(7,16) +(7,17) +(8,3) +(8,5) +(8,7) +(8,8) +(8,10) +(8,13) +(8,16) +(8,17) +(8,21) + +9 -9 +10 -9 +11 -9 +12 -9 +4 -9 +5 -9 +9 -10 +10 -10 +11 -10 +12 -10 +13 -10 +14 -10 +4 -10 +5 -10 +8 -10 +9 -11 +10 -11 +11 -11 +12 -11 +9 -12 +10 -12 +11 -12 +12 -12 +13 -12 +14 -12 +10 -13 +12 -13 +13 -13 +14 -13 +17 -13 +5 -13 +8 -13 +21 -13 +22 -13 +10 -14 +12 -14 +13 -14 +14 -14 +21 -14 +22 -14 +15 -15 +16 -15 +18 -15 +19 -15 +6 -15 +7 -15 +15 -16 +16 -16 +17 -16 +18 -16 +19 -16 +20 -16 +8 -16 +6 -16 +7 -16 +13 -17 +16 -17 +17 -17 +19 -17 +20 -17 +8 -17 +7 -17 +21 -17 +23 -17 +15 -18 +16 -18 +18 -18 +19 -18 +15 -19 +16 -19 +17 -19 +18 -19 +19 -19 +20 -19 +16 -20 +17 -20 +19 -20 +20 -20 +21 -20 +23 -20 +(9,9) +(9,10) +(9,11) +(9,12) +(9,4) +(9,5) +(10,9) +(10,10) +(10,11) +(10,12) +(10,13) +(10,14) +(10,4) +(10,5) +(10,8) +(11,9) +(11,10) +(11,11) +(11,12) +(12,9) +(12,10) +(12,11) +(12,12) +(12,13) +(12,14) +(13,10) +(13,12) +(13,13) +(13,14) +(13,17) +(13,5) +(13,8) +(13,21) +(13,22) +(14,10) +(14,12) +(14,13) +(14,14) +(14,21) +(14,22) +(15,15) +(15,16) +(15,18) +(15,19) +(15,6) +(15,7) +(16,15) +(16,16) +(16,17) +(16,18) +(16,19) +(16,20) +(16,8) +(16,6) +(16,7) +(17,13) +(17,16) +(17,17) +(17,19) +(17,20) +(17,8) +(17,7) +(17,21) +(17,23) +(18,15) +(18,16) +(18,18) +(18,19) +(19,15) +(19,16) +(19,17) +(19,18) +(19,19) +(19,20) +(20,16) +(20,17) +(20,19) +(20,20) +(20,21) +(20,23) + + +21 -21 +22 -21 +23 -21 +24 -21 +8 -21 +13 -21 +14 -21 +17 -21 +20 -21 +21 -22 +22 -22 +23 -22 +24 -22 +13 -22 +14 -22 +21 -23 +22 -23 +23 -23 +24 -23 +17 -23 +20 -23 +21 -24 +22 -24 +23 -24 +24 -24 +(21,21) +(21,22) +(21,23) +(21,24) +(21,8) +(21,13) +(21,14) +(21,17) +(21,20) +(22,21) +(22,22) +(22,23) +(22,24) +(22,13) +(22,14) +(23,21) +(23,22) +(23,23) +(23,24) +(23,17) +(23,20) +(24,21) +(24,22) +(24,23) +(24,24) + -- 2.39.5