]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clone the _04 test and also check const_iterators.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 9 Mar 2015 19:15:11 +0000 (14:15 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 9 Mar 2015 19:15:11 +0000 (14:15 -0500)
tests/mpi/trilinos_sparse_matrix_05.cc [new file with mode: 0644]
tests/mpi/trilinos_sparse_matrix_05.with_trilinos=true.mpirun=2.output [new file with mode: 0644]

diff --git a/tests/mpi/trilinos_sparse_matrix_05.cc b/tests/mpi/trilinos_sparse_matrix_05.cc
new file mode 100644 (file)
index 0000000..3ff863f
--- /dev/null
@@ -0,0 +1,146 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2004 - 2015 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 TrilinosWrappers::SparseMatrix::iterator semantics. make sure
+// that rows not stored locally look like they're empty
+//
+// this test is like _04 but uses the const_iterator type which goes
+// through a different code path in one place
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/index_set.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/trilinos_vector.h>
+#include <deal.II/lac/parallel_vector.h>
+#include <fstream>
+#include <iostream>
+#include <vector>
+
+
+void test ()
+{
+  const unsigned int n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  const unsigned int my_id = Utilities::MPI::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 (my_id == 0)
+    {
+      sp.add (0, 0);
+      sp.add (0, 2);
+    }
+  if ((n_procs == 1) || (my_id == 1))
+    sp.add(2,3);
+  sp.compress();
+
+  TrilinosWrappers::SparseMatrix A;
+  A.reinit (sp);
+  if (my_id==0)
+    {
+      A.set(0, 0, 0.1);
+      A.set(0, 2, 0.2);
+    }
+  if ((n_procs == 1) || (my_id == 1))
+    {
+      A.set(0, 0, 0.1);
+      A.set(0, 2, 0.2);
+      A.set(2,3, 0.3);
+    }
+
+  A.compress(VectorOperation::insert);
+
+  const TrilinosWrappers::SparseMatrix &B = A;
+  
+  // now access elements by iterator. ensure that we can iterate over
+  // all rows but that iterators into rows not stored locally just
+  // look empty
+  for (TrilinosWrappers::SparseMatrix::const_iterator p=B.begin();
+       p != B.end();
+       ++p)
+    if (my_id == 0)
+      {
+       deallog << "Looking at entry (" << p->row() << ','
+               << p->column() << ") with value "
+               << p->value()
+               << std::endl;
+       AssertThrow (p->row() == 0, ExcInternalError());
+      }
+    else
+      {
+       AssertThrow (p->row() == 2, ExcInternalError());
+      }  
+  
+  if (my_id == 0)
+    deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int);
+
+  const unsigned int n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile("output");
+      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_05.with_trilinos=true.mpirun=2.output b/tests/mpi/trilinos_sparse_matrix_05.with_trilinos=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..bd2018a
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL:0::Looking at entry (0,0) with value 0.1000
+DEAL:0::Looking at entry (0,2) with value 0.2000
+DEAL:0::OK

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.