From: Bruno Turcksin <bruno.turcksin@gmail.com>
Date: Sun, 2 Jun 2019 22:58:10 +0000 (-0400)
Subject: Add function to print a sparse matrix in a format readable by numpy
X-Git-Tag: v9.2.0-rc1~1420^2~4
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=39aca32692b99ad44857f4cd904905fac98f7764;p=dealii.git

Add function to print a sparse matrix in a format readable by numpy
---

diff --git a/include/deal.II/lac/sparse_matrix.h b/include/deal.II/lac/sparse_matrix.h
index eeddd5de0e..b7e6a9b34f 100644
--- a/include/deal.II/lac/sparse_matrix.h
+++ b/include/deal.II/lac/sparse_matrix.h
@@ -1583,6 +1583,17 @@ public:
   void
   print_pattern(std::ostream &out, const double threshold = 0.) const;
 
+  /**
+   * Print the vector to the output stream @p out in a format that can be
+   * read by numpy::readtxt(). To load the matrix in python just do
+   * <code>
+   *  [data, row, column] = numpy.loadtxt('my_matrix.txt')
+   *  sparse_matrix = scipy.sparse.csr_matrix((data, (row, column)))
+   * </code>
+   */
+  void
+  print_numpy_array(std::ostream &out, const unsigned int precision = 9) const;
+
   /**
    * Write the data of this object en bloc to a file. This is done in a binary
    * mode, so the output is neither readable by humans nor (probably) by other
diff --git a/include/deal.II/lac/sparse_matrix.templates.h b/include/deal.II/lac/sparse_matrix.templates.h
index a1b62de95b..a9b7bd40d2 100644
--- a/include/deal.II/lac/sparse_matrix.templates.h
+++ b/include/deal.II/lac/sparse_matrix.templates.h
@@ -1965,6 +1965,54 @@ SparseMatrix<number>::print_pattern(std::ostream &out,
 
 
 
+template <typename number>
+void
+SparseMatrix<number>::print_numpy_array(std::ostream &     out,
+                                        const unsigned int precision) const
+{
+  AssertThrow(out, ExcIO());
+  boost::io::ios_flags_saver restore_flags(out);
+
+  out.precision(precision);
+
+  Assert(cols != nullptr, ExcNotInitialized());
+  Assert(val != nullptr, ExcNotInitialized());
+
+  std::vector<number> rows;
+  std::vector<number> columns;
+  std::vector<number> values;
+  rows.reserve(n_nonzero_elements());
+  columns.reserve(n_nonzero_elements());
+  values.reserve(n_nonzero_elements());
+
+  for (size_type i = 0; i < cols->rows; ++i)
+    {
+      for (size_type j = cols->rowstart[i]; j < cols->rowstart[i + 1]; ++j)
+        {
+          rows.push_back(i);
+          columns.push_back(cols->colnums[j]);
+          values.push_back(j);
+        }
+    }
+
+  for (auto d : values)
+    out << d << ' ';
+  out << '\n';
+
+  for (auto r : rows)
+    out << r << ' ';
+  out << '\n';
+
+  for (auto c : columns)
+    out << c << ' ';
+  out << '\n';
+  out << std::flush;
+
+  AssertThrow(out, ExcIO());
+}
+
+
+
 template <typename number>
 void
 SparseMatrix<number>::block_write(std::ostream &out) const