]> https://gitweb.dealii.org/ - dealii.git/commitdiff
ReadWriteVector and TrilinosWrappers::MPI::Vector print in parallel
authorTimo Heister <timo.heister@gmail.com>
Sun, 3 Sep 2017 19:09:41 +0000 (15:09 -0400)
committerTimo Heister <timo.heister@gmail.com>
Sun, 3 Sep 2017 19:19:07 +0000 (15:19 -0400)
Print index and value for ::print() of ReadWriteVector and
TrilinosWrappers::MPI::Vector.

include/deal.II/lac/read_write_vector.templates.h
source/lac/trilinos_vector.cc
tests/lac/readwritevector_print.cc [new file with mode: 0644]
tests/lac/readwritevector_print.output [new file with mode: 0644]
tests/trilinos/print_01.cc [new file with mode: 0644]
tests/trilinos/print_01.mpirun=2.output [new file with mode: 0644]

index 7938f76efa50ebfaf97fad168376b6fc18e06174..722ca3903094e131ea99a99bffe0a2ca8ee97fdc 100644 (file)
@@ -25,6 +25,8 @@
 #include <deal.II/lac/vector_operations_internal.h>
 #include <deal.II/lac/la_parallel_vector.h>
 
+#include <boost/io/ios_state.hpp>
+
 #ifdef DEAL_II_WITH_PETSC
 #  include <deal.II/lac/petsc_parallel_vector.h>
 #endif
@@ -444,8 +446,7 @@ namespace LinearAlgebra
                                   const bool         scientific) const
   {
     AssertThrow (out, ExcIO());
-    std::ios::fmtflags old_flags = out.flags();
-    unsigned int old_precision = out.precision (precision);
+    boost::io::ios_flags_saver restore_flags(out);
 
     out.precision (precision);
     if (scientific)
@@ -456,15 +457,12 @@ namespace LinearAlgebra
     out << "IndexSet: ";
     stored_elements.print(out);
     out << std::endl;
-    for (unsigned int i=0; i<this->n_elements(); ++i)
-      out << val[i] << std::endl;
+    unsigned int i=0;
+    for (const auto &idx : this->stored_elements)
+      out << "[" << idx << "]: " << val[i++] << '\n';
     out << std::flush;
 
-
     AssertThrow (out, ExcIO());
-    // reset output format
-    out.flags (old_flags);
-    out.precision(old_precision);
   }
 
 
index daecaa7e31596a1fae65d40bc27a5a64162cba99..1e731a0143d04ab95b17a60b72a35305c96fa0b0 100644 (file)
@@ -761,33 +761,39 @@ namespace TrilinosWrappers
       AssertThrow (out, ExcIO());
       boost::io::ios_flags_saver restore_flags(out);
 
-      // get a representation of the
-      // vector and loop over all
-      // the elements TODO: up to
-      // now only local data printed
-      // out! Find a way to neatly
-      // output distributed data...
-      TrilinosScalar *val;
-      int leading_dimension;
-      int ierr = vector->ExtractView (&val, &leading_dimension);
 
-      AssertThrow (ierr == 0, ExcTrilinosError(ierr));
       out.precision (precision);
       if (scientific)
         out.setf (std::ios::scientific, std::ios::floatfield);
       else
         out.setf (std::ios::fixed, std::ios::floatfield);
 
-      if (across)
-        for (size_type i=0; i<size(); ++i)
-          out << static_cast<double>(val[i]) << ' ';
+      if (size() != local_size())
+        {
+          auto global_id = [&] (const size_type index)
+          {
+            return vector->Map().GID(static_cast<TrilinosWrappers::types::int_type>(index));
+          };
+          out << "size:" << size() << " local_size:" << local_size() << " :" << std::endl;
+          for (size_type i=0; i<local_size(); ++i)
+            out << "[" << global_id(i) << "]: " << (*(vector))[0][i] << std::endl;
+        }
       else
-        for (size_type i=0; i<size(); ++i)
-          out << static_cast<double>(val[i]) << std::endl;
-      out << std::endl;
+        {
+          TrilinosScalar *val;
+          int leading_dimension;
+          int ierr = vector->ExtractView (&val, &leading_dimension);
+
+          AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+          if (across)
+            for (size_type i=0; i<size(); ++i)
+              out << static_cast<double>(val[i]) << ' ';
+          else
+            for (size_type i=0; i<size(); ++i)
+              out << static_cast<double>(val[i]) << std::endl;
+          out << std::endl;
+        }
 
-      // restore the representation
-      // of the vector
       AssertThrow (out, ExcIO());
     }
 
diff --git a/tests/lac/readwritevector_print.cc b/tests/lac/readwritevector_print.cc
new file mode 100644 (file)
index 0000000..2226204
--- /dev/null
@@ -0,0 +1,44 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 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.
+//
+// ---------------------------------------------------------------------
+
+// Check printing
+
+#include "../tests.h"
+#include <deal.II/base/index_set.h>
+#include <deal.II/lac/read_write_vector.h>
+
+
+void test()
+{
+  IndexSet is(10);
+  is.add_range(1,3);
+  is.add_index(7);
+  LinearAlgebra::ReadWriteVector<double> vec(is);
+  deallog << "size: " << vec.n_elements() <<std::endl;
+
+  vec = 0.;
+  for (unsigned int i=0; i<vec.n_elements(); ++i)
+    vec.local_element(i) += i;
+
+  vec.print(deallog.get_file_stream());
+}
+
+int main()
+{
+  initlog();
+  test();
+
+  return 0;
+}
diff --git a/tests/lac/readwritevector_print.output b/tests/lac/readwritevector_print.output
new file mode 100644 (file)
index 0000000..b863c17
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL::size: 3
+IndexSet: {[1,2], 7}
+
+[1]: 0.000e+00
+[2]: 1.000e+00
+[7]: 2.000e+00
diff --git a/tests/trilinos/print_01.cc b/tests/trilinos/print_01.cc
new file mode 100644 (file)
index 0000000..608dbee
--- /dev/null
@@ -0,0 +1,68 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 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 TrilinosWrappers::MPI::Vector::print in parallel
+
+#include "../tests.h"
+#include <deal.II/base/index_set.h>
+#include <deal.II/lac/read_write_vector.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_vector.h>
+
+#include <vector>
+
+void test()
+{
+  IndexSet is(8);
+  unsigned int rank = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  if (rank==0)
+    is.add_range(0,4);
+  if (rank==1)
+    is.add_range(4,8);
+  is.compress();
+
+  deallog << "owned IndexSet: ";
+  is.print(deallog);
+
+  IndexSet is_ghosted(8);
+  is_ghosted.add_range(2,6);
+  deallog << "ghosted IndexSet: ";
+  is_ghosted.print(deallog);
+
+  TrilinosWrappers::MPI::Vector tril_vector(is);
+  TrilinosWrappers::MPI::Vector tril_vector_ghosted;
+  tril_vector_ghosted.reinit(is, is_ghosted, MPI_COMM_WORLD);
+  for (unsigned int i=0; i<8; ++i)
+    tril_vector[i] = i;
+  tril_vector.compress(VectorOperation::insert);
+
+  deallog << "trilinos vec:" << std::endl;
+  tril_vector.print(deallog.get_file_stream());
+
+  tril_vector_ghosted = tril_vector;
+  deallog << "trilinos vec ghosted:" << std::endl;
+  tril_vector_ghosted.print(deallog.get_file_stream());
+
+  deallog << "OK" <<std::endl;
+}
+
+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/print_01.mpirun=2.output b/tests/trilinos/print_01.mpirun=2.output
new file mode 100644 (file)
index 0000000..c363023
--- /dev/null
@@ -0,0 +1,37 @@
+
+DEAL:0::owned IndexSet: {[0,3]}
+DEAL:0::ghosted IndexSet: {[2,5]}
+DEAL:0::trilinos vec:
+size:8 local_size:4 :
+[0]: 0.000e+00
+[1]: 1.000e+00
+[2]: 2.000e+00
+[3]: 3.000e+00
+DEAL:0::trilinos vec ghosted:
+size:8 local_size:6 :
+[0]: 0.000e+00
+[1]: 1.000e+00
+[2]: 2.000e+00
+[3]: 3.000e+00
+[4]: 4.000e+00
+[5]: 5.000e+00
+DEAL:0::OK
+
+DEAL:1::owned IndexSet: {[4,7]}
+DEAL:1::ghosted IndexSet: {[2,5]}
+DEAL:1::trilinos vec:
+size:8 local_size:4 :
+[4]: 4.000e+00
+[5]: 5.000e+00
+[6]: 6.000e+00
+[7]: 7.000e+00
+DEAL:1::trilinos vec ghosted:
+size:8 local_size:6 :
+[2]: 2.000e+00
+[3]: 3.000e+00
+[4]: 4.000e+00
+[5]: 5.000e+00
+[6]: 6.000e+00
+[7]: 7.000e+00
+DEAL:1::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.