sparse_matrix.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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"
28 #include "libmesh/id_types.h"
31 
32 // C++ includes
33 #include <cstddef>
34 #include <iomanip>
35 #include <vector>
36 
37 namespace libMesh
38 {
39 
40 // forward declarations
41 template <typename T> class SparseMatrix;
42 template <typename T> class DenseMatrix;
43 class DofMap;
44 namespace SparsityPattern { class Graph; }
45 template <typename T> class NumericVector;
46 
47 // This template helper function must be declared before it
48 // can be defined below.
49 template <typename T>
50 std::ostream & operator << (std::ostream & os, const SparseMatrix<T> & m);
51 
52 
64 template <typename T>
65 class SparseMatrix : public ReferenceCountedObject<SparseMatrix<T> >,
66  public ParallelObject
67 {
68 public:
84  explicit
85  SparseMatrix (const Parallel::Communicator & comm
86  LIBMESH_CAN_DEFAULT_TO_COMMWORLD);
87 
93  virtual ~SparseMatrix ();
94 
99  static UniquePtr<SparseMatrix<T> >
100  build(const Parallel::Communicator & comm,
101  const SolverPackage solver_package = libMesh::default_solver_package());
102 
107  virtual bool initialized() const { return _is_initialized; }
108 
112  void attach_dof_map (const DofMap & dof_map)
113  { _dof_map = &dof_map; }
114 
122  virtual bool need_full_sparsity_pattern() const
123  { return false; }
124 
131 
142  virtual void init (const numeric_index_type m,
143  const numeric_index_type n,
144  const numeric_index_type m_l,
145  const numeric_index_type n_l,
146  const numeric_index_type nnz=30,
147  const numeric_index_type noz=10,
148  const numeric_index_type blocksize=1) = 0;
149 
153  virtual void init () = 0;
154 
161  virtual void clear () = 0;
162 
166  virtual void zero () = 0;
167 
171  virtual void zero_rows (std::vector<numeric_index_type> & rows, T diag_value = 0.0);
172 
178  virtual void close () const = 0;
179 
184  virtual numeric_index_type m () const = 0;
185 
190  virtual numeric_index_type n () const = 0;
191 
196  virtual numeric_index_type row_start () const = 0;
197 
202  virtual numeric_index_type row_stop () const = 0;
203 
210  virtual void set (const numeric_index_type i,
211  const numeric_index_type j,
212  const T value) = 0;
213 
222  virtual void add (const numeric_index_type i,
223  const numeric_index_type j,
224  const T value) = 0;
225 
232  virtual void add_matrix (const DenseMatrix<T> & dm,
233  const std::vector<numeric_index_type> & rows,
234  const std::vector<numeric_index_type> & cols) = 0;
235 
240  virtual void add_matrix (const DenseMatrix<T> & dm,
241  const std::vector<numeric_index_type> & dof_indices) = 0;
242 
250  virtual void add_block_matrix (const DenseMatrix<T> & dm,
251  const std::vector<numeric_index_type> & brows,
252  const std::vector<numeric_index_type> & bcols);
253 
258  virtual void add_block_matrix (const DenseMatrix<T> & dm,
259  const std::vector<numeric_index_type> & dof_indices)
260  { this->add_block_matrix (dm, dof_indices, dof_indices); }
261 
267  virtual void add (const T, SparseMatrix<T> &) = 0;
268 
276  virtual T operator () (const numeric_index_type i,
277  const numeric_index_type j) const = 0;
278 
289  virtual Real l1_norm () const = 0;
290 
302  virtual Real linfty_norm () const = 0;
303 
308  virtual bool closed() const = 0;
309 
315  void print(std::ostream & os=libMesh::out, const bool sparse=false) const;
316 
344  friend std::ostream & operator << <>(std::ostream & os, const SparseMatrix<T> & m);
345 
350  virtual void print_personal(std::ostream & os=libMesh::out) const = 0;
351 
358  virtual void print_matlab(const std::string & /*name*/ = "") const
359  {
360  libmesh_not_implemented();
361  }
362 
368  virtual void create_submatrix(SparseMatrix<T> & submatrix,
369  const std::vector<numeric_index_type> & rows,
370  const std::vector<numeric_index_type> & cols) const
371  {
372  this->_get_submatrix(submatrix,
373  rows,
374  cols,
375  false); // false means DO NOT REUSE submatrix
376  }
377 
384  virtual void reinit_submatrix(SparseMatrix<T> & submatrix,
385  const std::vector<numeric_index_type> & rows,
386  const std::vector<numeric_index_type> & cols) const
387  {
388  this->_get_submatrix(submatrix,
389  rows,
390  cols,
391  true); // true means REUSE submatrix
392  }
393 
398  void vector_mult (NumericVector<T> & dest,
399  const NumericVector<T> & arg) const;
400 
404  void vector_mult_add (NumericVector<T> & dest,
405  const NumericVector<T> & arg) const;
406 
410  virtual void get_diagonal (NumericVector<T> & dest) const = 0;
411 
416  virtual void get_transpose (SparseMatrix<T> & dest) const = 0;
417 
418 protected:
419 
426  const std::vector<numeric_index_type> & ,
427  const std::vector<numeric_index_type> & ,
428  const bool) const
429  {
430  libmesh_not_implemented();
431  }
432 
436  DofMap const * _dof_map;
437 
443 };
444 
445 
446 
447 //-----------------------------------------------------------------------
448 // SparseMatrix inline members
449 
450 // For SGI MIPSpro this implementation must occur after
451 // the full specialization of the print() member.
452 //
453 // It's generally easier to define these friend functions in the header
454 // file.
455 template <typename T>
456 std::ostream & operator << (std::ostream & os, const SparseMatrix<T> & m)
457 {
458  m.print(os);
459  return os;
460 }
461 
462 
463 } // namespace libMesh
464 
465 
466 #endif // LIBMESH_SPARSE_MATRIX_H
bool closed()
Definition: libmesh.C:279
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:988
Manages the degrees of freedom (DOFs) in a simulation.
Definition: dof_map.h:167
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
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:53
OStreamProxy out(std::cout)
virtual void add_block_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices)