1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 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 #ifndef LIBMESH_PETSC_SHELL_MATRIX_H
19 #define LIBMESH_PETSC_SHELL_MATRIX_H
20
21
22 #include "libmesh/libmesh_config.h"
23 #ifdef LIBMESH_HAVE_PETSC
24
25
26 // Local includes
27 #include "libmesh/libmesh_common.h"
28 #include "libmesh/reference_counted_object.h"
29 #include "libmesh/libmesh.h"
30 #include "libmesh/shell_matrix.h"
31 #include "libmesh/petsc_macro.h"
32 #include "libmesh/petsc_solver_exception.h"
33 #include "libmesh/petsc_vector.h"
34 #include "libmesh/libmesh_common.h"
35
36 // Petsc include files.
37 #ifdef I
38 # define LIBMESH_SAW_I
39 #endif
40 #include <petscmat.h>
41 #ifndef LIBMESH_SAW_I
42 # undef I // Avoid complex.h contamination
43 #endif
44
45 namespace libMesh
46 {
47
48 /**
49 * This class allows to use a PETSc shell matrix.
50 * All overridden virtual functions are documented in
51 * shell_matrix.h.
52 *
53 * \author Fande Kong (fdkong.jd@gmail.com)
54 * \date 2019
55 */
56 template <typename T>
57 class PetscShellMatrix : public ShellMatrix<T>
58 {
59 public:
60
61 PetscShellMatrix (const Parallel::Communicator & comm_in);
62
63 /**
64 * Destructor.
65 */
66 virtual ~PetscShellMatrix ();
67
68 virtual numeric_index_type m () const override;
69
70 virtual numeric_index_type n () const override;
71
72 virtual numeric_index_type local_m () const;
73
74 virtual numeric_index_type local_n () const;
75
76 virtual void vector_mult (NumericVector<T> & dest,
77 const NumericVector<T> & arg) const override;
78
79 virtual void vector_mult_add (NumericVector<T> & dest,
80 const NumericVector<T> & arg) const override;
81
82 virtual void get_diagonal (NumericVector<T> & dest) const override;
83
84 virtual void clear () override;
85
86 virtual void init () override;
87
88 /**
89 * \returns \p true if the matrix has been initialized,
90 * \p false otherwise.
91 */
92 virtual bool initialized() const;
93
94 /**
95 * Returns a pointer to the underlying PETSc Mat object. Must call
96 * init() before this.
97 */
98 Mat mat();
99
100 protected:
101
102 /**
103 * Petsc Shell Matrix
104 */
105 Mat _mat;
106
107 bool _is_initialized;
108 };
109
110
111 //-----------------------------------------------------------------------
112 // PetscShellMatrix inline members
113 template <typename T>
114 inline
PetscShellMatrix(const Parallel::Communicator & comm_in)115 PetscShellMatrix<T>::PetscShellMatrix (const Parallel::Communicator & comm_in):
116 ShellMatrix<T>(comm_in),
117 _mat(nullptr),
118 _is_initialized(false)
119 {}
120
121
122
123 template <typename T>
124 inline
~PetscShellMatrix()125 PetscShellMatrix<T>::~PetscShellMatrix ()
126 {
127 this->clear();
128 }
129
130
131
132 template <typename T>
133 inline
m()134 numeric_index_type PetscShellMatrix<T>::m () const
135 {
136 PetscErrorCode ierr;
137 PetscInt m;
138
139 ierr = MatGetSize(_mat,&m,nullptr);
140 LIBMESH_CHKERR(ierr);
141
142 return m;
143 }
144
145
146
147 template <typename T>
148 inline
n()149 numeric_index_type PetscShellMatrix<T>::n () const
150 {
151 PetscErrorCode ierr;
152 PetscInt n;
153
154 ierr = MatGetSize(_mat,nullptr,&n);
155 LIBMESH_CHKERR(ierr);
156
157 return n;
158 }
159
160
161 template <typename T>
162 inline
local_m()163 numeric_index_type PetscShellMatrix<T>::local_m () const
164 {
165 PetscErrorCode ierr;
166 PetscInt m;
167
168 ierr = MatGetLocalSize(_mat,&m,nullptr);
169 LIBMESH_CHKERR(ierr);
170
171 return m;
172 }
173
174
175
176 template <typename T>
177 inline
local_n()178 numeric_index_type PetscShellMatrix<T>::local_n () const
179 {
180 PetscErrorCode ierr;
181 PetscInt n;
182
183 ierr = MatGetLocalSize(_mat,nullptr,&n);
184 LIBMESH_CHKERR(ierr);
185
186 return n;
187 }
188
189
190 template <typename T>
191 inline
get_diagonal(NumericVector<T> & dest)192 void PetscShellMatrix<T>::get_diagonal (NumericVector<T> & dest) const
193 {
194 // Make sure the NumericVector passed in is really a PetscVector
195 PetscVector<T> & petsc_dest = cast_ref<PetscVector<T> &>(dest);
196
197 PetscErrorCode ierr =
198 MatGetDiagonal(_mat,petsc_dest.vec()); LIBMESH_CHKERR(ierr);
199 }
200
201
202
203 } // namespace libMesh
204
205 #endif // LIBMESH_HAVE_PETSC
206 #endif // LIBMESH_SPARSE_SHELL_MATRIX_H
207