const char *zero_string = " ",
const double denominator = 1.) 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
+ * computers using a different
+ * operating system of number
+ * format.
+ *
+ * The purpose of this function
+ * is that you can swap out
+ * matrices and sparsity pattern
+ * if you are short of memory,
+ * want to communicate between
+ * different programs, or allow
+ * objects to be persistent
+ * across different runs of the
+ * program.
+ */
+ void block_write (ostream &out) const;
+
+ /**
+ * Read data that has previously
+ * been written by
+ * @p{block_write} en block from
+ * a file. This is done using the
+ * inverse operations to the
+ * above function, so it is
+ * reasonably fast because the
+ * bitstream is not interpreted
+ * except for a few numbers up
+ * front.
+ *
+ * The object is resized on this
+ * operation, and all previous
+ * contents are lost. Note,
+ * however, that no checks are
+ * performed whether new data and
+ * the underlying
+ * @ref{SparsityPattern} object
+ * fit together. It is your
+ * responsibility to make sure
+ * that the sparsity pattern and
+ * the data to be read match.
+ *
+ * A primitive form of error
+ * checking is performed which
+ * will recognize the bluntest
+ * attempts to interpret some
+ * data as a vector stored
+ * bitwise to a file, but not
+ * more.
+ */
+ void block_read (std::istream &in);
+
/**
* Determine an estimate for the
* memory consumption (in bytes)
number *val;
/**
- * Allocated size of @p{val}. This
- * can be larger than the
- * actually used part if the size
- * of the matrix was reduced
- * somewhen in the past by
- * associating a sparsity pattern
- * with a smaller size to this
- * object somewhen, using the
+ * Allocated size of
+ * @p{val}. This can be larger
+ * than the actually used part if
+ * the size of the matrix was
+ * reduced somewhen in the past
+ * by associating a sparsity
+ * pattern with a smaller size to
+ * this object, using the
* @p{reinit} function.
*/
unsigned int max_len;
+template <typename number>
+void
+SparseMatrix<number>::block_write (ostream &out) const
+{
+ AssertThrow (out, ExcIO());
+
+ // first the simple objects,
+ // bracketed in [...]
+ out << '[' << max_len << "][";
+ // then write out real data
+ out.write (reinterpret_cast<const char*>(&val[0]),
+ reinterpret_cast<const char*>(&val[max_len])
+ - reinterpret_cast<const char*>(&val[0]));
+ out << ']';
+
+ AssertThrow (out, ExcIO());
+};
+
+
+
+template <typename number>
+void
+SparseMatrix<number>::block_read (istream &in)
+{
+ AssertThrow (in, ExcIO());
+
+ char c;
+
+ // first read in simple data
+ in >> c;
+ AssertThrow (c == '[', ExcIO());
+ in >> max_len;
+
+ in >> c;
+ AssertThrow (c == ']', ExcIO());
+ in >> c;
+ AssertThrow (c == '[', ExcIO());
+
+ // reallocate space
+ delete[] val;
+ val = new number[max_len];
+
+ // then read data
+ in.read (reinterpret_cast<char*>(&val[0]),
+ reinterpret_cast<char*>(&val[max_len])
+ - reinterpret_cast<char*>(&val[0]));
+ in >> c;
+ AssertThrow (c == ']', ExcIO());
+};
+
+
+
template <typename number>
unsigned int
SparseMatrix<number>::memory_consumption () const
if (std::fabs(A_res[i] - E_res[i]) > 1.e-14)
deallog << "SparseMatrix and SparseMatrixEZ differ!!!"
<< std::endl;
+
+ // dump A into a file, and re-read
+ // it, then check equality
+ std::ofstream tmp_write ("sparse_matrices.tmp");
+ A.block_write (tmp_write);
+ tmp_write.close ();
+
+ std::ifstream tmp_read ("sparse_matrices.tmp");
+ SparseMatrix<double> A_tmp;
+ A_tmp.reinit (A.get_sparsity_pattern());
+ A_tmp.block_read (tmp_read);
+ tmp_read.close ();
+
+ for (unsigned int i=0; i<A.n_nonzero_elements(); ++i)
+ Assert (std::fabs(A.global_entry(i) - A_tmp.global_entry(i)) <=
+ std::fabs(1e-14*A.global_entry(i)),
+ ExcInternalError());
}