OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
sparsemtrx.h
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #ifndef sparsemtrx_h
36 #define sparsemtrx_h
37 
38 #include "oofemcfg.h"
39 #include "floatarray.h"
40 #include "floatmatrix.h"
41 #include "intarray.h"
42 #include "error.h"
43 #include "sparsemtrxtype.h"
44 
45 namespace oofem {
46 class EngngModel;
47 class TimeStep;
48 class UnknownNumberingScheme;
49 
60 class OOFEM_EXPORT SparseMtrx
61 {
62 public:
63  typedef long SparseMtrxVersionType;
64 
65 protected:
67  int nRows;
69  int nColumns;
70 
80  SparseMtrxVersionType version;
81 
82 public:
87  SparseMtrx(int n, int m) : nRows(n), nColumns(m), version(0) { }
89  SparseMtrx() : nRows(0), nColumns(0), version(0) { }
91  virtual ~SparseMtrx() { }
92 
94  SparseMtrxVersionType giveVersion() { return this->version; }
95 
102  void checkBounds(int i, int j) const {
103  if ( i <= 0 ) {
104  OOFEM_ERROR("matrix error on rows : %d <= 0", i);
105  } else if ( j <= 0 ) {
106  OOFEM_ERROR("matrix error on columns : %d <= 0", j);
107  } else if ( i > nRows ) {
108  OOFEM_ERROR("matrix error on rows : %d => %d", i, nRows);
109  } else if ( j > nColumns ) {
110  OOFEM_ERROR("matrix error on columns : %d => %d", j, nColumns);
111  }
112  }
114  int giveNumberOfRows() const { return nRows; }
116  int giveNumberOfColumns() const { return nColumns; }
118  bool isSquare() const { return nRows == nColumns; }
120  bool isNotEmpty() const { return nRows > 0 && nColumns > 0; }
121 
127  virtual SparseMtrx *GiveCopy() const {
128  OOFEM_ERROR("Not implemented");
129  return NULL;
130  }
131 
137  virtual void times(const FloatArray &x, FloatArray &answer) const { OOFEM_ERROR("Not implemented"); }
143  virtual void timesT(const FloatArray &x, FloatArray &answer) const { OOFEM_ERROR("Not implemented"); }
149  virtual void times(const FloatMatrix &B, FloatMatrix &answer) const { OOFEM_ERROR("Not implemented"); }
155  virtual void timesT(const FloatMatrix &B, FloatMatrix &answer) const { OOFEM_ERROR("Not implemented"); }
160  virtual void times(double x) { OOFEM_ERROR("Not implemented"); }
166  virtual void add(double x, SparseMtrx &m) { OOFEM_ERROR("Not implemented"); }
172  virtual void addDiagonal(double x, FloatArray &m) {
173  for ( int i = 1; i <= m.giveSize(); ++i ) {
174  this->at(i, i) += x * m.at(i);
175  }
176  }
185  virtual int buildInternalStructure(EngngModel *eModel, int n, int m, const IntArray &I, const IntArray &J) { OOFEM_ERROR("Not implemented"); }
198  virtual int buildInternalStructure(EngngModel *eModel, int di, const UnknownNumberingScheme &s) = 0;
208  virtual int buildInternalStructure(EngngModel *eModel, int di, const UnknownNumberingScheme &r_s,
209  const UnknownNumberingScheme &c_s) {
210  OOFEM_ERROR("Not implemented");
211  return 0;
212  }
221  virtual int assemble(const IntArray &loc, const FloatMatrix &mat) = 0;
232  virtual int assemble(const IntArray &rloc, const IntArray &cloc, const FloatMatrix &mat) = 0;
233 
235  virtual int assembleBegin() { return 1; }
237  virtual int assembleEnd() { return 1; }
238 
240  virtual bool canBeFactorized() const = 0;
245  virtual SparseMtrx *factorized() { return NULL; }
253  virtual FloatArray *backSubstitutionWith(FloatArray &y) const { return NULL; }
255  virtual void zero() = 0;
256 
258  virtual double computeNorm() const {
259  OOFEM_ERROR("Not implemented");
260  return 0.0;
261  }
262 
263  virtual SparseMtrx *giveSubMatrix(const IntArray &rows, const IntArray &cols)
264  { OOFEM_ERROR("Not implemented"); return NULL; }
265 
267  virtual double &at(int i, int j) = 0;
269  virtual double at(int i, int j) const = 0;
271  virtual bool isAllocatedAt(int i, int j) const { return false; }
273  virtual void toFloatMatrix(FloatMatrix &answer) const { OOFEM_ERROR("Not implemented"); }
275  virtual void printStatistics() const { OOFEM_LOG_INFO("Not implemented file %s, line %d\n", __FILE__, __LINE__); }
277  virtual void printYourself() const { OOFEM_LOG_INFO("Not implemented file %s, line %d\n", __FILE__, __LINE__); }
279  virtual void writeToFile(const char *fname) const { OOFEM_LOG_INFO("Not implemented file %s, line %d\n", __FILE__, __LINE__); }
281  virtual SparseMtrxType giveType() const = 0;
283  virtual bool isAsymmetric() const = 0;
284 
285  virtual const char *giveClassName() const = 0;
287  std :: string errorInfo(const char *func) const { return std :: string(giveClassName()) + func; }
289 
290  FloatArray operator *( const FloatArray & x ) const
292  {
293  FloatArray answer;
294  this->times(x, answer);
295  return answer;
296  }
299  {
300  FloatArray answer;
301  this->timesT(x, answer);
302  return answer;
303  }
305 };
306 } // end namespace oofem
307 #endif // sparsemtrx_h
bool isNotEmpty() const
Tests for empty matrix.
Definition: sparsemtrx.h:120
int nColumns
Number of columns.
Definition: sparsemtrx.h:69
virtual void timesT(const FloatArray &x, FloatArray &answer) const
Evaluates .
Definition: sparsemtrx.h:143
Base class for all matrices stored in sparse format.
Definition: sparsemtrx.h:60
virtual FloatArray * backSubstitutionWith(FloatArray &y) const
Computes the solution of linear system where A is receiver.
Definition: sparsemtrx.h:253
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual double computeNorm() const
Returns the norm of receiver.
Definition: sparsemtrx.h:258
SparseMtrx()
Constructor.
Definition: sparsemtrx.h:89
Class implementing an array of integers.
Definition: intarray.h:61
virtual void add(double x, SparseMtrx &m)
Adds x * m.
Definition: sparsemtrx.h:166
virtual SparseMtrx * giveSubMatrix(const IntArray &rows, const IntArray &cols)
Definition: sparsemtrx.h:263
virtual void printStatistics() const
Prints the receiver statistics (one-line) to stdout.
Definition: sparsemtrx.h:275
virtual void addDiagonal(double x, FloatArray &m)
Adds x * m (treats m as a diagonal matrix, stored as an array)
Definition: sparsemtrx.h:172
int nRows
Number of rows.
Definition: sparsemtrx.h:67
int giveNumberOfRows() const
Returns number of rows of receiver.
Definition: sparsemtrx.h:114
virtual void toFloatMatrix(FloatMatrix &answer) const
Converts receiving sparse matrix to a dense float matrix.
Definition: sparsemtrx.h:273
#define OOFEM_LOG_INFO(...)
Definition: logger.h:127
virtual int buildInternalStructure(EngngModel *eModel, int n, int m, const IntArray &I, const IntArray &J)
Builds internal structure of receiver based on I and J.
Definition: sparsemtrx.h:185
virtual int assembleBegin()
Starts assembling the elements.
Definition: sparsemtrx.h:235
virtual void printYourself() const
Prints receiver to stdout. Works only for relatively small matrices.
Definition: sparsemtrx.h:277
void checkBounds(int i, int j) const
Checks size of receiver towards requested bounds.
Definition: sparsemtrx.h:102
#define OOFEM_ERROR(...)
Definition: error.h:61
long SparseMtrxVersionType
Definition: sparsemtrx.h:63
SparseMtrxType
Enumerative type used to identify the sparse matrix type.
Abstract base class allowing to control the way, how equations are assigned to individual DOFs...
virtual void times(const FloatMatrix &B, FloatMatrix &answer) const
Evaluates .
Definition: sparsemtrx.h:149
virtual int buildInternalStructure(EngngModel *eModel, int di, const UnknownNumberingScheme &r_s, const UnknownNumberingScheme &c_s)
Build internal structure of receiver.
Definition: sparsemtrx.h:208
SparseMtrxVersionType version
Allows to track if receiver changes.
Definition: sparsemtrx.h:80
virtual SparseMtrx * GiveCopy() const
Returns a newly allocated copy of receiver.
Definition: sparsemtrx.h:127
Class representing vector of real numbers.
Definition: floatarray.h:82
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
FloatArray trans_mult(const FloatArray &x) const
IML compatibility, .
Definition: sparsemtrx.h:298
virtual SparseMtrx * factorized()
Returns the receiver factorized.
Definition: sparsemtrx.h:245
bool isSquare() const
Returns nonzero if receiver is square matrix.
Definition: sparsemtrx.h:118
virtual void timesT(const FloatMatrix &B, FloatMatrix &answer) const
Evaluates .
Definition: sparsemtrx.h:155
FloatArray operator*(const double &a, const FloatArray &x)
Definition: floatarray.C:940
virtual int assembleEnd()
Returns when assemble is completed.
Definition: sparsemtrx.h:237
std::string errorInfo(const char *func) const
Error printing helper.
Definition: sparsemtrx.h:287
int giveNumberOfColumns() const
Returns number of columns of receiver.
Definition: sparsemtrx.h:116
SparseMtrx(int n, int m)
Constructor, creates (n,m) sparse matrix.
Definition: sparsemtrx.h:87
virtual void times(double x)
Multiplies receiver by scalar value.
Definition: sparsemtrx.h:160
Abstract base class representing the "problem" under consideration.
Definition: engngm.h:181
virtual bool isAllocatedAt(int i, int j) const
Checks whether memory is allocated at position (i,j).
Definition: sparsemtrx.h:271
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
virtual void times(const FloatArray &x, FloatArray &answer) const
Evaluates .
Definition: sparsemtrx.h:137
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual ~SparseMtrx()
Destructor.
Definition: sparsemtrx.h:91
virtual void writeToFile(const char *fname) const
Helpful for debugging, writes the matrix to given file.
Definition: sparsemtrx.h:279
SparseMtrxVersionType giveVersion()
Return receiver version.
Definition: sparsemtrx.h:94

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:31 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011