dune-istl 3.0-git
io.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_ISTL_IO_HH
4#define DUNE_ISTL_IO_HH
5
6#include <cmath>
7#include <complex>
8#include <limits>
9#include <ios>
10#include <iomanip>
11#include <fstream>
12#include <string>
13
14#include "matrixutils.hh"
15#include "istlexception.hh"
16#include <dune/common/fvector.hh>
17#include <dune/common/fmatrix.hh>
18#include <dune/common/dynmatrix.hh>
19#include <dune/common/diagonalmatrix.hh>
20#include <dune/common/unused.hh>
21
22#include <dune/istl/matrix.hh>
24#include "bcrsmatrix.hh"
25
26
27namespace Dune {
28
41 //
42 // pretty printing of vectors
43 //
44
52 template<class V>
53 void recursive_printvector (std::ostream& s, const V& v, std::string rowtext,
54 int& counter, int columns, int width,
55 int precision)
56 {
57 for (typename V::ConstIterator i=v.begin(); i!=v.end(); ++i)
59 }
60
68 template<class K, int n>
69 void recursive_printvector (std::ostream& s, const FieldVector<K,n>& v,
70 std::string rowtext, int& counter, int columns,
71 int width, int precision)
72 {
74 // we now can print n numbers
75 for (int i=0; i<n; i++)
76 {
77 if (counter%columns==0)
78 {
79 s << rowtext; // start a new row
80 s << " "; // space in front of each entry
81 s.width(4); // set width for counter
82 s << counter; // number of first entry in a line
83 }
84 s << " "; // space in front of each entry
85 s.width(width); // set width for each entry anew
86 s << v[i]; // yeah, the number !
87 counter++; // increment the counter
88 if (counter%columns==0)
89 s << std::endl; // start a new line
90 }
91 }
92
93
101 template<class V>
102 void printvector (std::ostream& s, const V& v, std::string title,
103 std::string rowtext, int columns=1, int width=10,
104 int precision=2)
105 {
106 // count the numbers printed to make columns
107 int counter=0;
108
109 // remember old flags
110 std::ios_base::fmtflags oldflags = s.flags();
111
112 // set the output format
113 s.setf(std::ios_base::scientific, std::ios_base::floatfield);
114 int oldprec = s.precision();
115 s.precision(precision);
116
117 // print title
118 s << title << " [blocks=" << v.N() << ",dimension=" << v.dim() << "]"
119 << std::endl;
120
121 // print data from all blocks
123
124 // check if new line is required
125 if (counter%columns!=0)
126 s << std::endl;
127
128 // reset the output format
129 s.flags(oldflags);
130 s.precision(oldprec);
131 }
132
133
135 //
136 // pretty printing of matrices
137 //
138
146 inline void fill_row (std::ostream& s, int m, int width, int precision)
147 {
149 for (int j=0; j<m; j++)
150 {
151 s << " "; // space in front of each entry
152 s.width(width); // set width for each entry anew
153 s << "."; // yeah, the number !
154 }
155 }
156
164 template<class M>
165 void print_row (std::ostream& s, const M& A, typename M::size_type I,
166 typename M::size_type J, typename M::size_type therow,
167 int width, int precision)
168 {
169 typename M::size_type i0=I;
170 for (typename M::size_type i=0; i<A.N(); i++)
171 {
173 {
174 // the row is in this block row !
175 typename M::size_type j0=J;
176 for (typename M::size_type j=0; j<A.M(); j++)
177 {
178 // find this block
179 typename M::ConstColIterator it = A[i].find(j);
180
181 // print row or filler
182 if (it!=A[i].end())
184 else
186
187 // advance columns
189 }
190 }
191 // advance rows
193 }
194 }
195
203 template<class K, int n, int m>
204 void print_row (std::ostream& s, const FieldMatrix<K,n,m>& A,
208 int precision)
209 {
212
213 typedef typename FieldMatrix<K,n,m>::size_type size_type;
214
215 for (size_type i=0; i<n; i++)
216 if (I+i==therow)
217 for (int j=0; j<m; j++)
218 {
219 s << " "; // space in front of each entry
220 s.width(width); // set width for each entry anew
221 s << A[i][j]; // yeah, the number !
222 }
223 }
224
232 template<class K>
233 void print_row (std::ostream& s, const FieldMatrix<K,1,1>& A,
237 int width, int precision)
238 {
241
242 if (I==therow)
243 {
244 s << " "; // space in front of each entry
245 s.width(width); // set width for each entry anew
246 s << static_cast<K>(A); // yeah, the number !
247 }
248 }
249
258 template<class M>
259 void printmatrix (std::ostream& s, const M& A, std::string title,
260 std::string rowtext, int width=10, int precision=2)
261 {
262
263 // remember old flags
264 std::ios_base::fmtflags oldflags = s.flags();
265
266 // set the output format
267 s.setf(std::ios_base::scientific, std::ios_base::floatfield);
268 int oldprec = s.precision();
269 s.precision(precision);
270
271 // print title
272 s << title
273 << " [n=" << A.N()
274 << ",m=" << A.M()
275 << ",rowdim=" << MatrixDimension<M>::rowdim(A)
276 << ",coldim=" << MatrixDimension<M>::coldim(A)
277 << "]" << std::endl;
278
279 // print all rows
280 for (typename M::size_type i=0; i<MatrixDimension<M>::rowdim(A); i++)
281 {
282 s << rowtext; // start a new row
283 s << " "; // space in front of each entry
284 s.width(4); // set width for counter
285 s << i; // number of first entry in a line
286 print_row(s,A,0,0,i,width,precision); // generic print
287 s << std::endl; // start a new line
288 }
289
290 // reset the output format
291 s.flags(oldflags);
292 s.precision(oldprec);
293 }
294
316 template<class B, int n, int m, class A>
317 void printSparseMatrix(std::ostream& s,
319 std::string title, std::string rowtext,
320 int width=3, int precision=2)
321 {
323 // remember old flags
324 std::ios_base::fmtflags oldflags = s.flags();
325 // set the output format
326 s.setf(std::ios_base::scientific, std::ios_base::floatfield);
327 int oldprec = s.precision();
328 s.precision(precision);
329 // print title
330 s << title
331 << " [n=" << mat.N()
332 << ",m=" << mat.M()
333 << ",rowdim=" << MatrixDimension<Matrix>::rowdim(mat)
334 << ",coldim=" << MatrixDimension<Matrix>::coldim(mat)
335 << "]" << std::endl;
336
337 typedef typename Matrix::ConstRowIterator Row;
338
339 for(Row row=mat.begin(); row != mat.end(); ++row) {
340 int skipcols=0;
341 bool reachedEnd=false;
342
343 while(!reachedEnd) {
344 for(int innerrow=0; innerrow<n; ++innerrow) {
345 int count=0;
346 typedef typename Matrix::ConstColIterator Col;
347 Col col=row->begin();
348 for(; col != row->end(); ++col,++count) {
349 if(count<skipcols)
350 continue;
351 if(count>=skipcols+width)
352 break;
353 if(innerrow==0) {
354 if(count==skipcols) {
355 s << rowtext; // start a new row
356 s << " "; // space in front of each entry
357 s.width(4); // set width for counter
358 s << row.index()<<": "; // number of first entry in a line
359 }
360 s.width(4);
361 s<<col.index()<<": |";
362 } else {
363 if(count==skipcols) {
364 for(typename std::string::size_type i=0; i < rowtext.length(); i++)
365 s<<" ";
366 s<<" ";
367 }
368 s<<" |";
369 }
370 for(int innercol=0; innercol < m; ++innercol) {
371 s.width(9);
372 s<<(*col)[innerrow][innercol]<<" ";
373 }
374
375 s<<"|";
376 }
377 if(innerrow==n-1 && col==row->end())
378 reachedEnd = true;
379 else
380 s << std::endl;
381 }
382 skipcols += width;
383 s << std::endl;
384 }
385 s << std::endl;
386 }
387
388 // reset the output format
389 s.flags(oldflags);
390 s.precision(oldprec);
391 }
392
393 namespace
394 {
395 template<typename T>
396 struct MatlabPODWriter
397 {
398 static std::ostream& write(const T& t, std::ostream& s)
399 {
400 s << t;
401 return s;
402 }
403 };
404 template<typename T>
405 struct MatlabPODWriter<std::complex<T> >
406 {
407 static std::ostream& write(const std::complex<T>& t, std::ostream& s)
408 {
409 s << t.real() << " " << t.imag();
410 return s;
411 }
412 };
413 } // anonymous namespace
414
424 template <class FieldType, int dim>
426 {
427 for (int i=0; i<dim; i++)
428 {
429 //+1 for Matlab numbering
430 s << rowOffset + i + 1 << " " << colOffset + i + 1 << " ";
431 MatlabPODWriter<FieldType>::write(matrix.scalar(), s)<< std::endl;
432 }
433 }
434
444 template <class FieldType, int dim>
446 {
447 for (int i=0; i<dim; i++)
448 {
449 //+1 for Matlab numbering
450 s << rowOffset + i + 1 << " " << colOffset + i + 1 << " ";
451 MatlabPODWriter<FieldType>::write(matrix.diagonal(i), s)<< std::endl;
452 }
453 }
454
464 template <class FieldType, int rows, int cols>
466 ( const FieldMatrix<FieldType,rows,cols>& matrix, int rowOffset,
467 int colOffset, std::ostream& s)
468 {
469 for (int i=0; i<rows; i++)
470 for (int j=0; j<cols; j++) {
471 //+1 for Matlab numbering
472 s << rowOffset + i + 1 << " " << colOffset + j + 1 << " ";
473 MatlabPODWriter<FieldType>::write(matrix[i][j], s)<< std::endl;
474 }
475 }
476
486 template <class FieldType>
488 int colOffset, std::ostream& s)
489 {
490 for (int i=0; i<matrix.N(); i++)
491 for (int j=0; j<matrix.M(); j++) {
492 //+1 for Matlab numbering
493 s << rowOffset + i + 1 << " " << colOffset + j + 1 << " ";
494 MatlabPODWriter<FieldType>::write(matrix[i][j], s)<< std::endl;
495 }
496 }
497
505 template <class MatrixType>
506 void writeMatrixToMatlabHelper(const MatrixType& matrix,
508 std::ostream& s)
509 {
510 // Precompute the accumulated sizes of the columns
511 std::vector<typename MatrixType::size_type> colOffset(matrix.M());
512 if (colOffset.size() > 0)
513 colOffset[0] = 0;
514
515 for (typename MatrixType::size_type i=0; i<matrix.M()-1; i++)
516 colOffset[i+1] = colOffset[i] +
518
519 typename MatrixType::size_type rowOffset = 0;
520
521 // Loop over all matrix rows
522 for (typename MatrixType::size_type rowIdx=0; rowIdx<matrix.N(); rowIdx++)
523 {
524
525 const typename MatrixType::row_type& row = matrix[rowIdx];
526
527 typename MatrixType::row_type::ConstIterator cIt = row.begin();
528 typename MatrixType::row_type::ConstIterator cEndIt = row.end();
529
530 // Loop over all columns in this row
531 for (; cIt!=cEndIt; ++cIt)
535 s);
536
538 }
539
540 }
541
561 template <class MatrixType>
562 void writeMatrixToMatlab(const MatrixType& matrix,
563 const std::string& filename, int outputPrecision = 18)
564 {
565 std::ofstream outStream(filename.c_str());
566 int oldPrecision = outStream.precision();
567 outStream.precision(outputPrecision);
568
570 outStream.precision(oldPrecision);
571 }
572
575} // namespace Dune
576
577#endif
Implementation of the BCRSMatrix class.
A dynamic dense block matrix class.
Some handy generic functions for ISTL matrices.
This file implements a quadratic matrix of fixed size which is a multiple of the identity.
Col col
Definition matrixmatrix.hh:347
Matrix & mat
Definition matrixmatrix.hh:343
void writeMatrixToMatlab(const MatrixType &matrix, const std::string &filename, int outputPrecision=18)
Writes sparse matrix in a Matlab-readable format.
Definition io.hh:562
void recursive_printvector(std::ostream &s, const V &v, std::string rowtext, int &counter, int columns, int width, int precision)
Recursively print all the blocks.
Definition io.hh:53
void printmatrix(std::ostream &s, const M &A, std::string title, std::string rowtext, int width=10, int precision=2)
Print a generic block matrix.
Definition io.hh:259
void printvector(std::ostream &s, const V &v, std::string title, std::string rowtext, int columns=1, int width=10, int precision=2)
Print an ISTL vector.
Definition io.hh:102
void print_row(std::ostream &s, const M &A, typename M::size_type I, typename M::size_type J, typename M::size_type therow, int width, int precision)
Print one row of a matrix.
Definition io.hh:165
void printSparseMatrix(std::ostream &s, const BCRSMatrix< FieldMatrix< B, n, m >, A > &mat, std::string title, std::string rowtext, int width=3, int precision=2)
Prints a BCRSMatrix with fixed sized blocks.
Definition io.hh:317
void fill_row(std::ostream &s, int m, int width, int precision)
Print a row of zeros for a non-existing block.
Definition io.hh:146
void writeMatrixToMatlabHelper(const ScaledIdentityMatrix< FieldType, dim > &matrix, int rowOffset, int colOffset, std::ostream &s)
Helper method for the writeMatrixToMatlab routine.
Definition io.hh:425
Definition basearray.hh:19
Iterator implementation class
Definition basearray.hh:84
Definition matrixutils.hh:231
Statistics about compression achieved in implicit mode.
Definition bcrsmatrix.hh:81
A sparse block matrix with compressed row storage.
Definition bcrsmatrix.hh:412
ConstIterator class for sequential access.
Definition matrix.hh:398
A generic dynamic dense matrix.
Definition matrix.hh:555
RowIterator end()
Get iterator to one beyond last row.
Definition matrix.hh:615
RowIterator begin()
Get iterator to first row.
Definition matrix.hh:609
size_type M() const
Return the number of columns.
Definition matrix.hh:695
size_type N() const
Return the number of rows.
Definition matrix.hh:690
Definition matrixutils.hh:25