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 
62 template <typename T>
63 class SparseMatrix : public ReferenceCountedObject<SparseMatrix<T> >,
64  public ParallelObject
65 {
66 public:
78  explicit
79  SparseMatrix (const Parallel::Communicator & comm
80  LIBMESH_CAN_DEFAULT_TO_COMMWORLD);
81 
86  virtual ~SparseMatrix ();
87 
92  static UniquePtr<SparseMatrix<T> >
93  build(const Parallel::Communicator & comm,
94  const SolverPackage solver_package = libMesh::default_solver_package());
95 
100  virtual bool initialized() const { return _is_initialized; }
101 
105  void attach_dof_map (const DofMap & dof_map)
106  { _dof_map = &dof_map; }
107 
117  virtual bool need_full_sparsity_pattern() const
118  { return false; }
119 
126 
139  virtual void init (const numeric_index_type m,
140  const numeric_index_type n,
141  const numeric_index_type m_l,
142  const numeric_index_type n_l,
143  const numeric_index_type nnz=30,
144  const numeric_index_type noz=10,
145  const numeric_index_type blocksize=1) = 0;
146 
150  virtual void init () = 0;
151 
155  virtual void clear () = 0;
156 
160  virtual void zero () = 0;
161 
165  virtual void zero_rows (std::vector<numeric_index_type> & rows, T diag_value = 0.0);
166 
171  virtual void close () = 0;
172 
176  virtual numeric_index_type m () const = 0;
177 
181  virtual numeric_index_type n () const = 0;
182 
187  virtual numeric_index_type row_start () const = 0;
188 
193  virtual numeric_index_type row_stop () const = 0;
194 
200  virtual void set (const numeric_index_type i,
201  const numeric_index_type j,
202  const T value) = 0;
203 
209  virtual void add (const numeric_index_type i,
210  const numeric_index_type j,
211  const T value) = 0;
212 
217  virtual void add_matrix (const DenseMatrix<T> & dm,
218  const std::vector<numeric_index_type> & rows,
219  const std::vector<numeric_index_type> & cols) = 0;
220 
225  virtual void add_matrix (const DenseMatrix<T> & dm,
226  const std::vector<numeric_index_type> & dof_indices) = 0;
227 
234  virtual void add_block_matrix (const DenseMatrix<T> & dm,
235  const std::vector<numeric_index_type> & brows,
236  const std::vector<numeric_index_type> & bcols);
237 
242  virtual void add_block_matrix (const DenseMatrix<T> & dm,
243  const std::vector<numeric_index_type> & dof_indices)
244  { this->add_block_matrix (dm, dof_indices, dof_indices); }
245 
249  virtual void add (const T a, SparseMatrix<T> & X) = 0;
250 
257  virtual T operator () (const numeric_index_type i,
258  const numeric_index_type j) const = 0;
259 
268  virtual Real l1_norm () const = 0;
269 
280  virtual Real linfty_norm () const = 0;
281 
285  virtual bool closed() const = 0;
286 
291  void print(std::ostream & os=libMesh::out, const bool sparse=false) const;
292 
320  friend std::ostream & operator << <>(std::ostream & os, const SparseMatrix<T> & m);
321 
326  virtual void print_personal(std::ostream & os=libMesh::out) const = 0;
327 
333  virtual void print_matlab(const std::string & /*name*/ = "") const
334  {
335  libmesh_not_implemented();
336  }
337 
343  virtual void create_submatrix(SparseMatrix<T> & submatrix,
344  const std::vector<numeric_index_type> & rows,
345  const std::vector<numeric_index_type> & cols) const
346  {
347  this->_get_submatrix(submatrix,
348  rows,
349  cols,
350  false); // false means DO NOT REUSE submatrix
351  }
352 
359  virtual void reinit_submatrix(SparseMatrix<T> & submatrix,
360  const std::vector<numeric_index_type> & rows,
361  const std::vector<numeric_index_type> & cols) const
362  {
363  this->_get_submatrix(submatrix,
364  rows,
365  cols,
366  true); // true means REUSE submatrix
367  }
368 
373  void vector_mult (NumericVector<T> & dest,
374  const NumericVector<T> & arg) const;
375 
380  void vector_mult_add (NumericVector<T> & dest,
381  const NumericVector<T> & arg) const;
382 
386  virtual void get_diagonal (NumericVector<T> & dest) const = 0;
387 
392  virtual void get_transpose (SparseMatrix<T> & dest) const = 0;
393 
394 protected:
395 
403  virtual void _get_submatrix(SparseMatrix<T> & /*submatrix*/,
404  const std::vector<numeric_index_type> & /*rows*/,
405  const std::vector<numeric_index_type> & /*cols*/,
406  const bool /*reuse_submatrix*/) const
407  {
408  libmesh_not_implemented();
409  }
410 
414  DofMap const * _dof_map;
415 
420 };
421 
422 
423 
424 //-----------------------------------------------------------------------
425 // SparseMatrix inline members
426 
427 // For SGI MIPSpro this implementation must occur after
428 // the full specialization of the print() member.
429 //
430 // It's generally easier to define these friend functions in the header
431 // file.
432 template <typename T>
433 std::ostream & operator << (std::ostream & os, const SparseMatrix<T> & m)
434 {
435  m.print(os);
436  return os;
437 }
438 
439 
440 } // namespace libMesh
441 
442 
443 #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:54
OStreamProxy out(std::cout)
virtual void add_block_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices)