<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.
{
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;
}
}
--- /dev/null
+//---------------------------------------------------------------------------
+// $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();
+}
--- /dev/null
+
+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
--- /dev/null
+
+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
--- /dev/null
+
+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
--- /dev/null
+//---------------------------- 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();
+
+}
--- /dev/null
+
+(2,3) 2.00000
--- /dev/null
+
+(2,3) 2.00000