#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
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)
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);
}
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());
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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;
+}
--- /dev/null
+
+DEAL::size: 3
+IndexSet: {[1,2], 7}
+
+[1]: 0.000e+00
+[2]: 1.000e+00
+[7]: 2.000e+00
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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();
+}
--- /dev/null
+
+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
+