]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix a bug where we get column indices wrong in TrilinosWrappers::SparseMatrix::print.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 10 May 2011 15:27:04 +0000 (15:27 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 10 May 2011 15:27:04 +0000 (15:27 +0000)
git-svn-id: https://svn.dealii.org/trunk@23685 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/source/lac/trilinos_sparse_matrix.cc
tests/mpi/habib_01.cc [new file with mode: 0644]
tests/mpi/habib_01/ncpu_10/cmp/generic [new file with mode: 0644]
tests/mpi/habib_01/ncpu_3/cmp/generic [new file with mode: 0644]
tests/mpi/habib_01/ncpu_4/cmp/generic [new file with mode: 0644]
tests/mpi/trilinos_sparse_matrix_print_01.cc [new file with mode: 0644]
tests/mpi/trilinos_sparse_matrix_print_01/ncpu_1/cmp/generic [new file with mode: 0644]
tests/mpi/trilinos_sparse_matrix_print_01/ncpu_2/cmp/generic [new file with mode: 0644]

index 5e7fa9daacb31cb2352e3ef437be26e763f8f5ee..cdfa4979216be4acf2ecc21836fd1fb52a1f8c20 100644 (file)
@@ -100,6 +100,11 @@ should be fixed now.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Fixed: The TrilinosWrappers::SparseMatrix::print() function
+didn't get column indices right. This is now fixed.
+<br>
+(Habib Talavatifard, Wolfgang Bangerth 2011/05/10)
+
 <li> Fixed: The TrilinosWrappers::SparseMatrix::operator() and
 TrilinosWrappers::SparseMatrix::el() functions sometimes produced
 wrong results for rectangular matrices. This is now fixed.
index 7449c8accfb87614393f6409f293ee828137601b..274c4c8a077415051edeb457d8c17422d9aefc06 100644 (file)
@@ -1330,7 +1330,7 @@ namespace TrilinosWrappers
          {
            matrix->ExtractMyRowView (i, num_entries, values, indices);
            for (int j=0; j<num_entries; ++j)
-             out << "(" << matrix->GRID(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 (file)
index 0000000..1219d57
--- /dev/null
@@ -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 <base/logstream.h>
+#include <base/utilities.h>
+
+#include <fstream>
+//#include <mpi.h>
+
+
+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<numprocs;++i)
+    {
+      MPI_Barrier(MPI_COMM_WORLD);
+//      system("sleep 1");
+
+      if (myid==0)
+       {
+         unsigned int buf=numbers::invalid_unsigned_int;
+         MPI_Status status;
+         MPI_Recv(&buf, 1, MPI_UNSIGNED, i, 1, MPI_COMM_WORLD, &status);
+         deallog << "got message '" << buf << "' from CPU " << i+1 << "!" << std::endl;
+         Assert(buf == i, ExcInternalError());
+
+       }
+      else if (myid==i )
+       {
+         MPI_Send(&myid, 1, MPI_INT, 0, 1, MPI_COMM_WORLD);
+       }
+
+
+    }
+  if (myid==0)
+    deallog << "done" << std::endl;
+
+}
+
+
+int main(int argc, char *argv[])
+{
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+  Utilities::System::MPI_InitFinalize mpi (argc, argv);
+#else
+  (void)argc;
+  (void)argv;
+  compile_time_error;
+
+#endif
+
+  if (Utilities::System::get_this_mpi_process (MPI_COMM_WORLD) == 0)
+    {
+      std::ofstream logfile(output_file_for_mpi("simple_mpi_01").c_str());
+      deallog.attach(logfile);
+      deallog.depth_console(0);
+      deallog.threshold_double(1.e-10);
+
+      deallog.push("mpi");
+      test_mpi();
+      deallog.pop();
+    }
+  else
+    test_mpi();
+}
diff --git a/tests/mpi/habib_01/ncpu_10/cmp/generic b/tests/mpi/habib_01/ncpu_10/cmp/generic
new file mode 100644 (file)
index 0000000..30b5ea4
--- /dev/null
@@ -0,0 +1,12 @@
+
+DEAL:mpi::Running on 10 CPU(s).
+DEAL:mpi::got message '1' from CPU 2!
+DEAL:mpi::got message '2' from CPU 3!
+DEAL:mpi::got message '3' from CPU 4!
+DEAL:mpi::got message '4' from CPU 5!
+DEAL:mpi::got message '5' from CPU 6!
+DEAL:mpi::got message '6' from CPU 7!
+DEAL:mpi::got message '7' from CPU 8!
+DEAL:mpi::got message '8' from CPU 9!
+DEAL:mpi::got message '9' from CPU 10!
+DEAL:mpi::done
diff --git a/tests/mpi/habib_01/ncpu_3/cmp/generic b/tests/mpi/habib_01/ncpu_3/cmp/generic
new file mode 100644 (file)
index 0000000..588f4ae
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL:mpi::Running on 3 CPU(s).
+DEAL:mpi::got message '1' from CPU 2!
+DEAL:mpi::got message '2' from CPU 3!
+DEAL:mpi::done
diff --git a/tests/mpi/habib_01/ncpu_4/cmp/generic b/tests/mpi/habib_01/ncpu_4/cmp/generic
new file mode 100644 (file)
index 0000000..8dcacbe
--- /dev/null
@@ -0,0 +1,6 @@
+
+DEAL:mpi::Running on 4 CPU(s).
+DEAL:mpi::got message '1' from CPU 2!
+DEAL:mpi::got message '2' from CPU 3!
+DEAL:mpi::got message '3' from CPU 4!
+DEAL:mpi::done
diff --git a/tests/mpi/trilinos_sparse_matrix_print_01.cc b/tests/mpi/trilinos_sparse_matrix_print_01.cc
new file mode 100644 (file)
index 0000000..41048fe
--- /dev/null
@@ -0,0 +1,104 @@
+//----------------------------  trilinos_vector_equality_4.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004, 2005, 2008, 2010, 2011 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.
+//
+//----------------------------  trilinos_vector_equality_4.cc  ---------------------------
+
+
+// TrilinosWrappers::SparseMatrix::print got column indices wrong
+
+#include "../tests.h"
+#include <base/utilities.h>
+#include <base/index_set.h>
+#include <lac/trilinos_sparse_matrix.h>
+#include <lac/trilinos_sparsity_pattern.h>
+#include <fstream>
+#include <iostream>
+#include <vector>
+
+
+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 (file)
index 0000000..dfa3dbd
--- /dev/null
@@ -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 (file)
index 0000000..dfa3dbd
--- /dev/null
@@ -0,0 +1,2 @@
+
+(2,3) 2.00000

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.