]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
SparseMatrix::block_read/write
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Sep 2002 23:18:33 +0000 (23:18 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Sep 2002 23:18:33 +0000 (23:18 +0000)
git-svn-id: https://svn.dealii.org/trunk@6379 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/2002/c-3-4.html
deal.II/lac/include/lac/sparse_matrix.h
deal.II/lac/include/lac/sparse_matrix.templates.h
tests/lac/sparse_matrices.cc

index bbea7bf7574f77a631f2ca8ce19a4ee5bd529385..3284318ab0d397f84a3282b9949659f079a72bcf 100644 (file)
@@ -207,10 +207,11 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 
 <ol>
   <li> <p>
-       New: <code class="class">SparsityPattern</code> now has functions <code
+       New: Classes <code class="class">SparsityPattern</code> and <code
+       class="class">SparseMatrix</code> now have functions <code 
        class="member">block_write/block_read</code>, allowing to dump the data
-       of this object in to a file in binary format, and later to re-read it
-       without much need for much parsing.
+       of these objects into a file in binary format, and later to re-read it
+       without much need for parsing.
        <br>
        (WB 2002/09/09)
        </p>
index 211a3517e4ffc2ab6f1d028b78621f74d18ab607..38d3b8967b66034a473a592c026d23642525315c 100644 (file)
@@ -774,6 +774,62 @@ class SparseMatrix : public Subscriptor
                          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)
@@ -849,14 +905,14 @@ class SparseMatrix : public Subscriptor
     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;
index eb72506f775911e6b306c94c7b4f3dc6bdf85f66..fcb345952c55bddc6c7ec98056a06fb49ca8c032 100644 (file)
@@ -1172,6 +1172,58 @@ void SparseMatrix<number>::print_formatted (std::ostream &out,
 
 
 
+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
index 9cfa6f85c86148046954e300fc4212709d36c2a7..94c7f47f831f7d2edb6763af8656a6cc0d43be20 100644 (file)
@@ -155,4 +155,21 @@ int main()
     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());
 }

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.