From 47356575379fdca81d126d573de2814505797d95 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 10 May 2011 15:27:04 +0000 Subject: [PATCH] Fix a bug where we get column indices wrong in TrilinosWrappers::SparseMatrix::print. git-svn-id: https://svn.dealii.org/trunk@23685 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 5 + deal.II/source/lac/trilinos_sparse_matrix.cc | 2 +- tests/mpi/habib_01.cc | 87 +++++++++++++++ tests/mpi/habib_01/ncpu_10/cmp/generic | 12 ++ tests/mpi/habib_01/ncpu_3/cmp/generic | 5 + tests/mpi/habib_01/ncpu_4/cmp/generic | 6 + tests/mpi/trilinos_sparse_matrix_print_01.cc | 104 ++++++++++++++++++ .../ncpu_1/cmp/generic | 2 + .../ncpu_2/cmp/generic | 2 + 9 files changed, 224 insertions(+), 1 deletion(-) create mode 100644 tests/mpi/habib_01.cc create mode 100644 tests/mpi/habib_01/ncpu_10/cmp/generic create mode 100644 tests/mpi/habib_01/ncpu_3/cmp/generic create mode 100644 tests/mpi/habib_01/ncpu_4/cmp/generic create mode 100644 tests/mpi/trilinos_sparse_matrix_print_01.cc create mode 100644 tests/mpi/trilinos_sparse_matrix_print_01/ncpu_1/cmp/generic create mode 100644 tests/mpi/trilinos_sparse_matrix_print_01/ncpu_2/cmp/generic diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 5e7fa9daac..cdfa497921 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -100,6 +100,11 @@ should be fixed now.

Specific improvements

    +
  1. Fixed: The TrilinosWrappers::SparseMatrix::print() function +didn't get column indices right. This is now fixed. +
    +(Habib Talavatifard, Wolfgang Bangerth 2011/05/10) +
  2. Fixed: The TrilinosWrappers::SparseMatrix::operator() and TrilinosWrappers::SparseMatrix::el() functions sometimes produced wrong results for rectangular matrices. This is now fixed. diff --git a/deal.II/source/lac/trilinos_sparse_matrix.cc b/deal.II/source/lac/trilinos_sparse_matrix.cc index 7449c8accf..274c4c8a07 100644 --- a/deal.II/source/lac/trilinos_sparse_matrix.cc +++ b/deal.II/source/lac/trilinos_sparse_matrix.cc @@ -1330,7 +1330,7 @@ namespace TrilinosWrappers { matrix->ExtractMyRowView (i, num_entries, values, indices); for (int j=0; jGRID(i) << "," << indices[matrix->GCID(j)] << ") " + out << "(" << matrix->GRID(i) << "," << matrix->GCID(indices[j]) << ") " << values[j] << std::endl; } } diff --git a/tests/mpi/habib_01.cc b/tests/mpi/habib_01.cc new file mode 100644 index 0000000000..1219d575e9 --- /dev/null +++ b/tests/mpi/habib_01.cc @@ -0,0 +1,87 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- + + +// + +#include "../tests.h" +#include +#include + +#include +//#include + + +void test_mpi() +{ + Assert( Utilities::System::job_supports_mpi(), ExcInternalError()); + + + unsigned int myid = Utilities::System::get_this_mpi_process (MPI_COMM_WORLD); + const unsigned int numprocs = Utilities::System::get_n_mpi_processes (MPI_COMM_WORLD); + + if (myid==0) + deallog << "Running on " << numprocs << " CPU(s)." << std::endl; + + for (unsigned int i=1;i +#include +#include +#include +#include +#include +#include + + +void test () +{ + const unsigned int n_procs = Utilities::System::get_n_mpi_processes(MPI_COMM_WORLD); + const unsigned int my_id = Utilities::System::get_this_mpi_process(MPI_COMM_WORLD); + + const unsigned int n_rows = 3; + const unsigned int n_cols = 4; + + IndexSet row_partitioning (n_rows); + IndexSet col_partitioning (n_cols); + + if (n_procs == 1) + { + row_partitioning.add_range(0, n_rows); + col_partitioning.add_range(0, n_cols); + } + else if (n_procs == 2) + { + // row_partitioning should be { [0, 2), [2, n_rows) } + // col_partitioning should be { [0, 2), [2, n_cols) } + // col_relevant_set should be { [0, 3), [1, n_cols) } + if (my_id == 0) + { + row_partitioning.add_range(0, 2); + col_partitioning.add_range(0, 2); + } + else if (my_id == 1) + { + row_partitioning.add_range(2, n_rows); + col_partitioning.add_range(2, n_cols); + } + } + else + Assert (false, ExcNotImplemented()); + + TrilinosWrappers::SparsityPattern sp (row_partitioning, + col_partitioning, MPI_COMM_WORLD); + if ((n_procs == 1) || (my_id == 1)) + sp.add(2,3); + sp.compress(); + + TrilinosWrappers::SparseMatrix A; + A.clear (); + A.reinit (sp); + if ((n_procs == 1) || (my_id == 1)) + A.add(2,3, 2.0); + A.compress(); + + if ((n_procs == 1) || (my_id == 1)) + A.print(deallog.get_file_stream()); +} + + + +int main (int argc, char **argv) +{ + Utilities::System::MPI_InitFinalize mpi_initialization(argc, argv); + + const unsigned int n_procs = Utilities::System::get_n_mpi_processes(MPI_COMM_WORLD); + unsigned int myid = Utilities::System::get_this_mpi_process (MPI_COMM_WORLD); + deallog.push(Utilities::int_to_string(myid)); + + // let processor 1 speak if we run + // in parallel + if ((n_procs == 1) || (myid == 1)) + { + std::ofstream logfile(output_file_for_mpi("trilinos_sparse_matrix_print_01").c_str()); + deallog.attach(logfile); + deallog << std::setprecision(4); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test(); + } + else + test(); + +} diff --git a/tests/mpi/trilinos_sparse_matrix_print_01/ncpu_1/cmp/generic b/tests/mpi/trilinos_sparse_matrix_print_01/ncpu_1/cmp/generic new file mode 100644 index 0000000000..dfa3dbd283 --- /dev/null +++ b/tests/mpi/trilinos_sparse_matrix_print_01/ncpu_1/cmp/generic @@ -0,0 +1,2 @@ + +(2,3) 2.00000 diff --git a/tests/mpi/trilinos_sparse_matrix_print_01/ncpu_2/cmp/generic b/tests/mpi/trilinos_sparse_matrix_print_01/ncpu_2/cmp/generic new file mode 100644 index 0000000000..dfa3dbd283 --- /dev/null +++ b/tests/mpi/trilinos_sparse_matrix_print_01/ncpu_2/cmp/generic @@ -0,0 +1,2 @@ + +(2,3) 2.00000 -- 2.39.5