sparse_matrix.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_SPARSE_MATRIX_H
21 #define LIBMESH_SPARSE_MATRIX_H
22 
23 
24 // Local includes
25 #include "libmesh/libmesh.h"
26 #include "libmesh/libmesh_common.h"
27 #include "libmesh/auto_ptr.h" // deprecated
28 #include "libmesh/id_types.h"
31 
32 // C++ includes
33 #include <cstddef>
34 #include <iomanip>
35 #include <vector>
36 #include <memory>
37 
38 namespace libMesh
39 {
40 
41 // forward declarations
42 template <typename T> class SparseMatrix;
43 template <typename T> class DenseMatrix;
44 class DofMap;
45 namespace SparsityPattern { class Graph; }
46 template <typename T> class NumericVector;
47 
48 // This template helper function must be declared before it
49 // can be defined below.
50 template <typename T>
51 std::ostream & operator << (std::ostream & os, const SparseMatrix<T> & m);
52 
53 
63 template <typename T>
64 class SparseMatrix : public ReferenceCountedObject<SparseMatrix<T>>,
65  public ParallelObject
66 {
67 public:
79  explicit
80  SparseMatrix (const Parallel::Communicator & comm
81  LIBMESH_CAN_DEFAULT_TO_COMMWORLD);
82 
87  virtual ~SparseMatrix ();
88 
93  static std::unique_ptr<SparseMatrix<T>>
94  build(const Parallel::Communicator & comm,
95  const SolverPackage solver_package = libMesh::default_solver_package());
96 
101  virtual bool initialized() const { return _is_initialized; }
102 
106  void attach_dof_map (const DofMap & dof_map)
107  { _dof_map = &dof_map; }
108 
118  virtual bool need_full_sparsity_pattern() const
119  { return false; }
120 
127 
140  virtual void init (const numeric_index_type m,
141  const numeric_index_type n,
142  const numeric_index_type m_l,
143  const numeric_index_type n_l,
144  const numeric_index_type nnz=30,
145  const numeric_index_type noz=10,
146  const numeric_index_type blocksize=1) = 0;
147 
151  virtual void init () = 0;
152 
156  virtual void clear () = 0;
157 
161  virtual void zero () = 0;
162 
166  virtual void zero_rows (std::vector<numeric_index_type> & rows, T diag_value = 0.0);
167 
172  virtual void close () = 0;
173 
179  virtual void flush () { close(); }
180 
184  virtual numeric_index_type m () const = 0;
185 
189  virtual numeric_index_type n () const = 0;
190 
195  virtual numeric_index_type row_start () const = 0;
196 
201  virtual numeric_index_type row_stop () const = 0;
202 
208  virtual void set (const numeric_index_type i,
209  const numeric_index_type j,
210  const T value) = 0;
211 
217  virtual void add (const numeric_index_type i,
218  const numeric_index_type j,
219  const T value) = 0;
220 
225  virtual void add_matrix (const DenseMatrix<T> & dm,
226  const std::vector<numeric_index_type> & rows,
227  const std::vector<numeric_index_type> & cols) = 0;
228 
233  virtual void add_matrix (const DenseMatrix<T> & dm,
234  const std::vector<numeric_index_type> & dof_indices) = 0;
235 
242  virtual void add_block_matrix (const DenseMatrix<T> & dm,
243  const std::vector<numeric_index_type> & brows,
244  const std::vector<numeric_index_type> & bcols);
245 
250  virtual void add_block_matrix (const DenseMatrix<T> & dm,
251  const std::vector<numeric_index_type> & dof_indices)
252  { this->add_block_matrix (dm, dof_indices, dof_indices); }
253 
257  virtual void add (const T a, SparseMatrix<T> & X) = 0;
258 
265  virtual T operator () (const numeric_index_type i,
266  const numeric_index_type j) const = 0;
267 
276  virtual Real l1_norm () const = 0;
277 
288  virtual Real linfty_norm () const = 0;
289 
293  virtual bool closed() const = 0;
294 
299  void print(std::ostream & os=libMesh::out, const bool sparse=false) const;
300 
328  friend std::ostream & operator << <>(std::ostream & os, const SparseMatrix<T> & m);
329 
334  virtual void print_personal(std::ostream & os=libMesh::out) const = 0;
335 
341  virtual void print_matlab(const std::string & /*name*/ = "") const
342  {
343  libmesh_not_implemented();
344  }
345 
351  virtual void create_submatrix(SparseMatrix<T> & submatrix,
352  const std::vector<numeric_index_type> & rows,
353  const std::vector<numeric_index_type> & cols) const
354  {
355  this->_get_submatrix(submatrix,
356  rows,
357  cols,
358  false); // false means DO NOT REUSE submatrix
359  }
360 
367  virtual void reinit_submatrix(SparseMatrix<T> & submatrix,
368  const std::vector<numeric_index_type> & rows,
369  const std::vector<numeric_index_type> & cols) const
370  {
371  this->_get_submatrix(submatrix,
372  rows,
373  cols,
374  true); // true means REUSE submatrix
375  }
376 
381  void vector_mult (NumericVector<T> & dest,
382  const NumericVector<T> & arg) const;
383 
388  void vector_mult_add (NumericVector<T> & dest,
389  const NumericVector<T> & arg) const;
390 
394  virtual void get_diagonal (NumericVector<T> & dest) const = 0;
395 
400  virtual void get_transpose (SparseMatrix<T> & dest) const = 0;
401 
402 protected:
403 
411  virtual void _get_submatrix(SparseMatrix<T> & /*submatrix*/,
412  const std::vector<numeric_index_type> & /*rows*/,
413  const std::vector<numeric_index_type> & /*cols*/,
414  const bool /*reuse_submatrix*/) const
415  {
416  libmesh_not_implemented();
417  }
418 
422  DofMap const * _dof_map;
423 
428 };
429 
430 
431 
432 //-----------------------------------------------------------------------
433 // SparseMatrix inline members
434 
435 // For SGI MIPSpro this implementation must occur after
436 // the full specialization of the print() member.
437 //
438 // It's generally easier to define these friend functions in the header
439 // file.
440 template <typename T>
441 std::ostream & operator << (std::ostream & os, const SparseMatrix<T> & m)
442 {
443  m.print(os);
444  return os;
445 }
446 
447 
448 } // namespace libMesh
449 
450 
451 #endif // LIBMESH_SPARSE_MATRIX_H
bool closed()
Definition: libmesh.C:280
virtual void print_matlab(const std::string &="") const
const Number zero
Definition: libmesh.h:178
virtual bool initialized() const
SolverPackage default_solver_package()
Definition: libmesh.C:994
Manages the degrees of freedom (DOFs) in a simulation.
Definition: dof_map.h:168
virtual void flush()
virtual void create_submatrix(SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const
dof_id_type numeric_index_type
Definition: id_types.h:92
void init(triangulateio &t)
virtual bool need_full_sparsity_pattern() const
virtual void reinit_submatrix(SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const
void attach_dof_map(const DofMap &dof_map)
DofMap const * _dof_map
virtual void update_sparsity_pattern(const SparsityPattern::Graph &)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const bool value
Definition: xdr_io.C:108
virtual void _get_submatrix(SparseMatrix< T > &, const std::vector< numeric_index_type > &, const std::vector< numeric_index_type > &, const bool) const
A matrix object used for finite element assembly and numerics.
Definition: dense_matrix.h:54
OStreamProxy out(std::cout)
virtual void add_block_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices)