petsc_matrix.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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_PETSC_MATRIX_H
21 #define LIBMESH_PETSC_MATRIX_H
22 
23 #include "libmesh/libmesh_common.h"
24 
25 #ifdef LIBMESH_HAVE_PETSC
26 
27 // Local includes
28 #include "libmesh/sparse_matrix.h"
29 #include "libmesh/petsc_macro.h"
30 #include "libmesh/parallel.h"
31 
32 // C++ includes
33 #include <algorithm>
34 
35 // Macro to identify and debug functions which should be called in
36 // parallel on parallel matrices but which may be called in serial on
37 // serial matrices. This macro will only be valid inside non-static
38 // PetscMatrix methods
39 #undef semiparallel_only
40 #ifndef NDEBUG
41  #include <cstring>
42 
43  #define semiparallel_only() do { if (this->initialized()) { const char *mytype; \
44  MatGetType(_mat,&mytype); \
45  if (!strcmp(mytype, MATSEQAIJ)) \
46  parallel_object_only(); } } while (0)
47 #else
48  #define semiparallel_only()
49 #endif
50 
51 
55 EXTERN_C_FOR_PETSC_BEGIN
56 # include <petscmat.h>
57 EXTERN_C_FOR_PETSC_END
58 
59 
60 
61 namespace libMesh
62 {
63 
64 // Forward Declarations
65 template <typename T> class DenseMatrix;
66 
67 
76 template <typename T>
77 class PetscMatrix : public SparseMatrix<T>
78 {
79 public:
95  explicit
96  PetscMatrix (const Parallel::Communicator &comm
97  LIBMESH_CAN_DEFAULT_TO_COMMWORLD);
98 
106  explicit
107  PetscMatrix (Mat m,
108  const Parallel::Communicator &comm
109  LIBMESH_CAN_DEFAULT_TO_COMMWORLD);
110 
116  ~PetscMatrix ();
117 
128  void init (const numeric_index_type m,
129  const numeric_index_type n,
130  const numeric_index_type m_l,
131  const numeric_index_type n_l,
132  const numeric_index_type nnz=30,
133  const numeric_index_type noz=10,
134  const numeric_index_type blocksize=1);
135 
146  void init (const numeric_index_type m,
147  const numeric_index_type n,
148  const numeric_index_type m_l,
149  const numeric_index_type n_l,
150  const std::vector<numeric_index_type>& n_nz,
151  const std::vector<numeric_index_type>& n_oz,
152  const numeric_index_type blocksize=1);
153 
157  void init ();
158 
165  void clear ();
166 
171  void zero ();
172 
176  void zero_rows (std::vector<numeric_index_type> & rows, T diag_value = 0.0);
177 
183  void close () const;
184 
189  numeric_index_type m () const;
190 
195  numeric_index_type n () const;
196 
201  numeric_index_type row_start () const;
202 
207  numeric_index_type row_stop () const;
208 
215  void set (const numeric_index_type i,
216  const numeric_index_type j,
217  const T value);
218 
227  void add (const numeric_index_type i,
228  const numeric_index_type j,
229  const T value);
230 
238  void add_matrix (const DenseMatrix<T> &dm,
239  const std::vector<numeric_index_type> &rows,
240  const std::vector<numeric_index_type> &cols);
241 
246  void add_matrix (const DenseMatrix<T> &dm,
247  const std::vector<numeric_index_type> &dof_indices);
248 
256  virtual void add_block_matrix (const DenseMatrix<T> &dm,
257  const std::vector<numeric_index_type> &brows,
258  const std::vector<numeric_index_type> &bcols);
259 
264  virtual void add_block_matrix (const DenseMatrix<T> &dm,
265  const std::vector<numeric_index_type> &dof_indices)
266  { this->add_block_matrix (dm, dof_indices, dof_indices); }
267 
279  void add (const T a, SparseMatrix<T> &X);
280 
289  const numeric_index_type j) const;
290 
302  Real l1_norm () const;
303 
316  Real linfty_norm () const;
317 
322  bool closed() const;
323 
331  void print_personal(std::ostream& os=libMesh::out) const;
332 
339  void print_matlab(const std::string name="NULL") const;
340 
344  virtual void get_diagonal (NumericVector<T>& dest) const;
345 
350  virtual void get_transpose (SparseMatrix<T>& dest) const;
351 
355  void swap (PetscMatrix<T> &);
356 
362  Mat mat () { libmesh_assert (_mat); return _mat; }
363 
364 protected:
365 
375  virtual void _get_submatrix(SparseMatrix<T>& submatrix,
376  const std::vector<numeric_index_type>& rows,
377  const std::vector<numeric_index_type>& cols,
378  const bool reuse_submatrix) const;
379 
380 private:
381 
385  Mat _mat;
386 
392 };
393 
394 } // namespace libMesh
395 
396 #endif // #ifdef LIBMESH_HAVE_PETSC
397 #endif // LIBMESH_PETSC_MATRIX_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:06 UTC

Hosted By:
SourceForge.net Logo